v0.13.1
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 /**
107 * \brief get measure of element
108 * @return area of face
109 */
110 inline double getMeasure();
111
112 /** \brief get triangle normal
113 */
114 inline VectorDouble &getNormal();
115
116 /** \brief get triangle tangent 1
117 */
118 inline VectorDouble &getTangent1();
119
120 /** \brief get triangle tangent 2
121 */
122 inline VectorDouble &getTangent2();
123
124 /** \brief get normal as tensor
125 */
126 inline auto getFTensor1Normal();
127
128 /** \brief get tangentOne as tensor
129 */
130 inline auto getFTensor1Tangent1();
131
132 /** \brief get tangentTwo as tensor
133 */
134 inline auto getFTensor1Tangent2();
135
136 /** \brief get element number of nodes
137 */
138 inline int getNumNodes();
139
140 /** \brief get element connectivity
141 */
142 inline const EntityHandle *getConn();
143
144 /** \brief get triangle coordinates
145 */
146 inline VectorDouble &getCoords();
147
148 /**
149 * \brief get get coords at gauss points
150
151 \code
152 FTensor::Index<'i',3> i;
153 FTensor::Tensor1<double,3> t_center;
154 auto t_coords = getFTensor1Coords();
155 t_center(i) = 0;
156 for(int nn = 0;nn!=3;nn++) {
157 t_center(i) += t_coords(i);
158 ++t_coords;
159 }
160 t_center(i) /= 3;
161 \endcode
162
163 */
164 inline auto getFTensor1Coords();
165
166 /** \brief if higher order geometry return normals at Gauss pts.
167
168 Note: returned matrix has size 0 in rows and columns if no HO approximation
169 of geometry is available.
170
171 */
173
174 /** \brief if higher order geometry return normals at Gauss pts.
175 *
176 * \param gg gauss point number
177 */
178 inline ublas::matrix_row<MatrixDouble> getNormalsAtGaussPts(const int gg);
179
180 /** \brief if higher order geometry return tangent vector to triangle at
181 Gauss pts.
182
183 Note: returned matrix has size 0 in rows and columns if no HO approximation
184 of geometry is avaliable.
185
186 */
188
189 /** \brief if higher order geometry return tangent vector to triangle at
190 Gauss pts.
191
192 Note: returned matrix has size 0 in rows and columns if no HO approximation
193 of geometry is avaliable.
194
195 */
197
198 /** \brief get normal at integration points
199
200 Example:
201 \code
202 double nrm2;
203 FTensor::Index<'i',3> i;
204 auto t_normal = getFTensor1NormalsAtGaussPts();
205 for(int gg = gg!=data.getN().size1();gg++) {
206 nrm2 = sqrt(t_normal(i)*t_normal(i));
207 ++t_normal;
208 }
209 \endcode
210
211 */
212 inline auto getFTensor1NormalsAtGaussPts();
213
214 /** \brief get tangent 1 at integration points
215
216 */
217 inline auto getFTensor1Tangent1AtGaussPts();
218
219 /** \brief get tangent 2 at integration points
220
221 */
222 inline auto getFTensor1Tangent2AtGaussPts();
223
224 /** \brief return pointer to Generic Triangle Finite Element object
225 */
227
228 /**
229 *
230 * User call this function to loop over elements on the side of face. This
231 * function calls MoFEM::VolumeElementForcesAndSourcesCoreOnSide with is
232 * operator to do calculations.
233 *
234 * @param fe_name Name of the element
235 * @param method Finite element object
236 * @return error code
237 */
239 loopSideVolumes(const string fe_name,
241
242private:
244};
245
246double FaceElementForcesAndSourcesCore::UserDataOperator::getArea() {
247 return static_cast<FaceElementForcesAndSourcesCore *>(ptrFE)->aRea;
248}
249
250double FaceElementForcesAndSourcesCore::UserDataOperator::getMeasure() {
251 return getArea();
252}
253
254VectorDouble &FaceElementForcesAndSourcesCore::UserDataOperator::getNormal() {
255 return static_cast<FaceElementForcesAndSourcesCore *>(ptrFE)->nOrmal;
256}
257
258VectorDouble &FaceElementForcesAndSourcesCore::UserDataOperator::getTangent1() {
259 return static_cast<FaceElementForcesAndSourcesCore *>(ptrFE)->tangentOne;
260}
261
262VectorDouble &FaceElementForcesAndSourcesCore::UserDataOperator::getTangent2() {
263 return static_cast<FaceElementForcesAndSourcesCore *>(ptrFE)->tangentTwo;
264}
265
266auto FaceElementForcesAndSourcesCore::UserDataOperator::getFTensor1Normal() {
267 double *ptr = &*getNormal().data().begin();
268 return FTensor::Tensor1<double *, 3>(ptr, &ptr[1], &ptr[2]);
269}
270
271auto FaceElementForcesAndSourcesCore::UserDataOperator::getFTensor1Tangent1() {
272 double *ptr = &*getTangent1().data().begin();
273 return FTensor::Tensor1<double *, 3>(ptr, &ptr[1], &ptr[2]);
274}
275
276auto FaceElementForcesAndSourcesCore::UserDataOperator::getFTensor1Tangent2() {
277 double *ptr = &*getTangent2().data().begin();
278 return FTensor::Tensor1<double *, 3>(ptr, &ptr[1], &ptr[2]);
279}
280
281int FaceElementForcesAndSourcesCore::UserDataOperator::getNumNodes() {
282 return static_cast<FaceElementForcesAndSourcesCore *>(ptrFE)->num_nodes;
283}
284
285const EntityHandle *
286FaceElementForcesAndSourcesCore::UserDataOperator::getConn() {
287 return static_cast<FaceElementForcesAndSourcesCore *>(ptrFE)->conn;
288}
289
290VectorDouble &FaceElementForcesAndSourcesCore::UserDataOperator::getCoords() {
291 return static_cast<FaceElementForcesAndSourcesCore *>(ptrFE)->coords;
292}
293
294auto FaceElementForcesAndSourcesCore::UserDataOperator::getFTensor1Coords() {
295 double *ptr = &*getCoords().data().begin();
297 &ptr[2]);
298}
299
301FaceElementForcesAndSourcesCore::UserDataOperator::getNormalsAtGaussPts() {
302 return static_cast<FaceElementForcesAndSourcesCore *>(ptrFE)
304}
305
306ublas::matrix_row<MatrixDouble>
307FaceElementForcesAndSourcesCore::UserDataOperator::getNormalsAtGaussPts(
308 const int gg) {
309 return ublas::matrix_row<MatrixDouble>(
311 gg);
312}
313
315FaceElementForcesAndSourcesCore::UserDataOperator::getTangent1AtGaussPts() {
316 return static_cast<FaceElementForcesAndSourcesCore *>(ptrFE)
318}
319
321FaceElementForcesAndSourcesCore::UserDataOperator::getTangent2AtGaussPts() {
322 return static_cast<FaceElementForcesAndSourcesCore *>(ptrFE)
324}
325
326auto FaceElementForcesAndSourcesCore::UserDataOperator::
327 getFTensor1NormalsAtGaussPts() {
328 double *ptr = &*getNormalsAtGaussPts().data().begin();
330 &ptr[2]);
331}
332
333auto FaceElementForcesAndSourcesCore::UserDataOperator::
334 getFTensor1Tangent1AtGaussPts() {
335 double *ptr = &*getTangent1AtGaussPts().data().begin();
337 &ptr[2]);
338}
339
340auto FaceElementForcesAndSourcesCore::UserDataOperator::
341 getFTensor1Tangent2AtGaussPts() {
342 double *ptr = &*getTangent2AtGaussPts().data().begin();
344 &ptr[2]);
345}
346
348FaceElementForcesAndSourcesCore::UserDataOperator::getFaceFE() {
349 return static_cast<FaceElementForcesAndSourcesCore *>(ptrFE);
350}
351
352/** \deprecated use FaceElementForcesAndSourcesCore
353 */
356
357} // namespace MoFEM
358
359#endif //__FACEELEMENTFORCESANDSOURCESCORE_HPP__
360
361/**
362 * \defgroup mofem_forces_and_sources_tri_element Face Element
363 * \brief Implementation of face element
364 *
365 * \ingroup mofem_forces_and_sources
366 **/
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.