v0.14.0
Loading...
Searching...
No Matches
FaceElementForcesAndSourcesCore.hpp
Go to the documentation of this file.
1/** \file FaceElementForcesAndSourcesCore.hpp
2 \brief Implementation of face element.
3
4*/
5
6
7
8#ifndef __FACEELEMENTFORCESANDSOURCESCORE_HPP__
9#define __FACEELEMENTFORCESANDSOURCESCORE_HPP__
10
11using namespace boost::numeric;
12
13namespace MoFEM {
14
15struct VolumeElementForcesAndSourcesCoreOnSide;
16
17/** \brief Face finite element
18 \ingroup mofem_forces_and_sources_tri_element
19
20 User is implementing own operator at Gauss point level, by own object
21 derived from FaceElementForcesAndSourcesCoreL::UserDataOperator. Arbitrary
22 number of operator added pushing objects to OpPtrVector
23
24 */
26
27 /**
28 * @deprecated not used anumore, will be removed in next versions
29 *
30 */
32
33 /** \brief default operator for TRI element
34 * \ingroup mofem_forces_and_sources_tri_element
35 */
36 struct UserDataOperator;
37
39
41
42protected:
43 /**
44 * \brief Calculate element area and normal of the face at integration points
45 *
46 * @return Error code
47 */
49
50 /**
51 * \brief Calculate element area and normal of the face
52 *
53 * Note that at that point is assumed that geometry is exclusively defined
54 * by corner nodes.
55 *
56 * @return Error code
57 */
59
60 /**
61 * \brief Set integration points
62 * @return Error code
63 */
65
66 /**
67 * \brief Determine approximation space and order of base functions
68 * @return Error code
69 */
71
72 /**
73 * \brief Calculate coordinate at integration points
74 * @return Error code
75 */
77
78 double &aRea;
83
87
88 friend class UserDataOperator;
90};
91
92/** \brief default operator for TRI element
93 * \ingroup mofem_forces_and_sources_tri_element
94 */
97
98 using ForcesAndSourcesCore::UserDataOperator::UserDataOperator;
99
100 /**
101 * \brief get area of face
102 * @return area of face
103 */
104 inline double getArea();
105
106 /** \brief get triangle normal
107 */
108 inline VectorDouble &getNormal();
109
110 /** \brief get triangle tangent 1
111 */
112 inline VectorDouble &getTangent1();
113
114 /** \brief get triangle tangent 2
115 */
116 inline VectorDouble &getTangent2();
117
118 /** \brief get normal as tensor
119 */
120 inline auto getFTensor1Normal();
121
122 /** \brief get tangentOne as tensor
123 */
124 inline auto getFTensor1Tangent1();
125
126 /** \brief get tangentTwo as tensor
127 */
128 inline auto getFTensor1Tangent2();
129
130 /** \brief get element number of nodes
131 */
132 inline int getNumNodes();
133
134 /** \brief get element connectivity
135 */
136 inline const EntityHandle *getConn();
137
138 /** \brief get triangle coordinates
139 */
140 inline VectorDouble &getCoords();
141
142 /**
143 * \brief get get coords at gauss points
144
145 \code
146 FTensor::Index<'i',3> i;
147 FTensor::Tensor1<double,3> t_center;
148 auto t_coords = getFTensor1Coords();
149 t_center(i) = 0;
150 for(int nn = 0;nn!=3;nn++) {
151 t_center(i) += t_coords(i);
152 ++t_coords;
153 }
154 t_center(i) /= 3;
155 \endcode
156
157 */
158 inline auto getFTensor1Coords();
159
160 /** \brief if higher order geometry return normals at Gauss pts.
161
162 Note: returned matrix has size 0 in rows and columns if no HO approximation
163 of geometry is available.
164
165 */
167
168 /** \brief if higher order geometry return normals at Gauss pts.
169 *
170 * \param gg gauss point number
171 */
172 inline ublas::matrix_row<MatrixDouble> getNormalsAtGaussPts(const int gg);
173
174 /** \brief if higher order geometry return tangent vector to triangle at
175 Gauss pts.
176
177 Note: returned matrix has size 0 in rows and columns if no HO approximation
178 of geometry is avaliable.
179
180 */
182
183 /** \brief if higher order geometry return tangent vector to triangle at
184 Gauss pts.
185
186 Note: returned matrix has size 0 in rows and columns if no HO approximation
187 of geometry is avaliable.
188
189 */
191
192 /** \brief get normal at integration points
193
194 Example:
195 \code
196 double nrm2;
197 FTensor::Index<'i',3> i;
198 auto t_normal = getFTensor1NormalsAtGaussPts();
199 for(int gg = gg!=data.getN().size1();gg++) {
200 nrm2 = sqrt(t_normal(i)*t_normal(i));
201 ++t_normal;
202 }
203 \endcode
204
205 */
206 inline auto getFTensor1NormalsAtGaussPts();
207
208 /** \brief get tangent 1 at integration points
209
210 */
211 inline auto getFTensor1Tangent1AtGaussPts();
212
213 /** \brief get tangent 2 at integration points
214
215 */
216 inline auto getFTensor1Tangent2AtGaussPts();
217
218 /** \brief return pointer to Generic Triangle Finite Element object
219 */
221
222 /**
223 *
224 * User call this function to loop over elements on the side of face. This
225 * function calls MoFEM::VolumeElementForcesAndSourcesCoreOnSide with is
226 * operator to do calculations.
227 *
228 * @param fe_name Name of the element
229 * @param method Finite element object
230 * @return error code
231 */
233 loopSideVolumes(const string fe_name,
235
236private:
238};
239
240double FaceElementForcesAndSourcesCore::UserDataOperator::getArea() {
241 return static_cast<FaceElementForcesAndSourcesCore *>(ptrFE)->aRea;
242}
243
244VectorDouble &FaceElementForcesAndSourcesCore::UserDataOperator::getNormal() {
245 return static_cast<FaceElementForcesAndSourcesCore *>(ptrFE)->nOrmal;
246}
247
248VectorDouble &FaceElementForcesAndSourcesCore::UserDataOperator::getTangent1() {
249 return static_cast<FaceElementForcesAndSourcesCore *>(ptrFE)->tangentOne;
250}
251
252VectorDouble &FaceElementForcesAndSourcesCore::UserDataOperator::getTangent2() {
253 return static_cast<FaceElementForcesAndSourcesCore *>(ptrFE)->tangentTwo;
254}
255
256auto FaceElementForcesAndSourcesCore::UserDataOperator::getFTensor1Normal() {
257 double *ptr = &*getNormal().data().begin();
258 return FTensor::Tensor1<double *, 3>(ptr, &ptr[1], &ptr[2]);
259}
260
261auto FaceElementForcesAndSourcesCore::UserDataOperator::getFTensor1Tangent1() {
262 double *ptr = &*getTangent1().data().begin();
263 return FTensor::Tensor1<double *, 3>(ptr, &ptr[1], &ptr[2]);
264}
265
266auto FaceElementForcesAndSourcesCore::UserDataOperator::getFTensor1Tangent2() {
267 double *ptr = &*getTangent2().data().begin();
268 return FTensor::Tensor1<double *, 3>(ptr, &ptr[1], &ptr[2]);
269}
270
271int FaceElementForcesAndSourcesCore::UserDataOperator::getNumNodes() {
272 return static_cast<FaceElementForcesAndSourcesCore *>(ptrFE)->num_nodes;
273}
274
275const EntityHandle *
276FaceElementForcesAndSourcesCore::UserDataOperator::getConn() {
277 return static_cast<FaceElementForcesAndSourcesCore *>(ptrFE)->conn;
278}
279
280VectorDouble &FaceElementForcesAndSourcesCore::UserDataOperator::getCoords() {
281 return static_cast<FaceElementForcesAndSourcesCore *>(ptrFE)->coords;
282}
283
284auto FaceElementForcesAndSourcesCore::UserDataOperator::getFTensor1Coords() {
285 double *ptr = &*getCoords().data().begin();
287 &ptr[2]);
288}
289
291FaceElementForcesAndSourcesCore::UserDataOperator::getNormalsAtGaussPts() {
292 return static_cast<FaceElementForcesAndSourcesCore *>(ptrFE)
294}
295
296ublas::matrix_row<MatrixDouble>
297FaceElementForcesAndSourcesCore::UserDataOperator::getNormalsAtGaussPts(
298 const int gg) {
299 return ublas::matrix_row<MatrixDouble>(
301 gg);
302}
303
305FaceElementForcesAndSourcesCore::UserDataOperator::getTangent1AtGaussPts() {
306 return static_cast<FaceElementForcesAndSourcesCore *>(ptrFE)
308}
309
311FaceElementForcesAndSourcesCore::UserDataOperator::getTangent2AtGaussPts() {
312 return static_cast<FaceElementForcesAndSourcesCore *>(ptrFE)
314}
315
316auto FaceElementForcesAndSourcesCore::UserDataOperator::
317 getFTensor1NormalsAtGaussPts() {
318 double *ptr = &*getNormalsAtGaussPts().data().begin();
320 &ptr[2]);
321}
322
323auto FaceElementForcesAndSourcesCore::UserDataOperator::
324 getFTensor1Tangent1AtGaussPts() {
325 double *ptr = &*getTangent1AtGaussPts().data().begin();
327 &ptr[2]);
328}
329
330auto FaceElementForcesAndSourcesCore::UserDataOperator::
331 getFTensor1Tangent2AtGaussPts() {
332 double *ptr = &*getTangent2AtGaussPts().data().begin();
334 &ptr[2]);
335}
336
338FaceElementForcesAndSourcesCore::UserDataOperator::getFaceFE() {
339 return static_cast<FaceElementForcesAndSourcesCore *>(ptrFE);
340}
341
342/** \deprecated use FaceElementForcesAndSourcesCore
343 */
346
347} // namespace MoFEM
348
349#endif //__FACEELEMENTFORCESANDSOURCESCORE_HPP__
350
351/**
352 * \defgroup mofem_forces_and_sources_tri_element Face Element
353 * \brief Implementation of face element
354 *
355 * \ingroup mofem_forces_and_sources
356 **/
T data[Tensor_Dim]
#define DEPRECATED
Definition: definitions.h:17
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:10
DEPRECATED typedef FaceElementForcesAndSourcesCore FaceElementForcesAndSourcesCoreBase
Deprecated interface functions.
MatrixDouble & getNormalsAtGaussPts()
if higher order geometry return normals at Gauss pts.
MoFEMErrorCode loopSideVolumes(const string fe_name, VolumeElementForcesAndSourcesCoreOnSide &fe_method)
FaceElementForcesAndSourcesCore * getFaceFE()
return pointer to Generic Triangle Finite Element object
MatrixDouble & getTangent2AtGaussPts()
if higher order geometry return tangent vector to triangle at Gauss pts.
MatrixDouble & getTangent1AtGaussPts()
if higher order geometry return tangent vector to triangle at Gauss pts.
virtual MoFEMErrorCode calculateAreaAndNormal()
Calculate element area and normal of the face.
MoFEMErrorCode operator()()
function is run for every finite element
virtual MoFEMErrorCode getSpaceBaseAndOrderOnElement()
Determine approximation space and order of base functions.
virtual MoFEMErrorCode calculateAreaAndNormalAtIntegrationPts()
Calculate element area and normal of the face at integration points.
virtual MoFEMErrorCode setIntegrationPts()
Set integration points.
virtual MoFEMErrorCode calculateCoordinatesAtGaussPts()
Calculate coordinate at integration points.
structure to get information form mofem into EntitiesFieldData
Base volume element used to integrate on skeleton.