v0.13.1
FlatPrismElementForcesAndSourcesCore.hpp
Go to the documentation of this file.
1/** \file ForcesAndSourcesCore.hpp
2 \brief Implementation of elements on entities.
3
4 Those element are inherited by user to implement specific implementation of
5 particular problem.
6
7*/
8
9
10
11using namespace boost::numeric;
12
13#ifndef __FLATPRISMELEMENTFORCESANDSURCESCORE_HPP__
14#define __FLATPRISMELEMENTFORCESANDSURCESCORE_HPP__
15
16namespace MoFEM {
17
18/** \brief FlatPrism finite element
19 \ingroup mofem_forces_and_sources_prism_element
20
21 User is implementing own operator at Gauss points level, by own object
22 derived from FlatPrismElementForcesAndSourcesCoreL::UserDataOperator. Arbitrary
23 number of operator added pushing objects to rowOpPtrVector and
24 rowColOpPtrVector.
25
26 */
28
30
31 /** \brief default operator for Flat Prism element
32 * \ingroup mofem_forces_and_sources_prism_element
33 */
35
37
38 /** \brief get face aRea
39 \param dd if dd == 0 it is for face F3 if dd == 1 is for face F4
40 */
41 inline double getArea(const int dd);
42
43 inline double getAreaF3();
44
45 inline double getAreaF4();
46
47 /** \brief get triangle normal
48
49 Normal has 6 elements, first 3 are for face F3 another three for face F4
50
51 */
52 inline VectorDouble &getNormal();
53
55
57
58 /** \brief get triangle coordinates
59
60 Vector has 6 elements, i.e. coordinates on face F3 and F4
61
62 */
63 inline VectorDouble &getCoords();
64
65 /** \brief get coordinates at Gauss pts.
66
67 Matrix has size (nb integration points)x(coordinates on F3 and F4 = 6),
68 i.e. coordinates on face F3 and F4
69
70 */
72
73 /** \brief coordinate at Gauss points on face 3 (if hierarchical
74 * approximation of element geometry)
75 */
77
78 /** \brief coordinate at Gauss points on face 4 (if hierarchical
79 * approximation of element geometry)
80 */
82
83 /** \brief if higher order geometry return normals at face F3 at Gauss pts.
84 *
85 * Face 3 is top face in canonical triangle numeration, see \cite
86 * tautges2010canonical
87 *
88 */
90
91 /** \brief if higher order geometry return normals at face F4 at Gauss pts.
92 *
93 * Face 4 is top face in canonical triangle numeration, see \cite
94 * tautges2010canonical
95 *
96 */
98
99 /** \brief if higher order geometry return normals at Gauss pts.
100 *
101 * Face 3 is top face in canonical triangle numeration, see \cite
102 * tautges2010canonical
103 *
104 * \param gg gauss point number
105 */
106 inline ublas::matrix_row<MatrixDouble> getNormalsAtGaussPtsF3(const int gg);
107
108 /** \brief if higher order geometry return normals at Gauss pts.
109 *
110 * Face 3 is top face in canonical triangle numeration, see \cite
111 * tautges2010canonical
112 *
113 * \param gg gauss point number
114 */
115 inline ublas::matrix_row<MatrixDouble> getNormalsAtGaussPtsF4(const int gg);
116
117 /** \brief if higher order geometry return tangent vector to triangle at
118 * Gauss pts.
119 */
121
122 /** \brief if higher order geometry return tangent vector to triangle at
123 * Gauss pts.
124 */
126
127 /** \brief if higher order geometry return tangent vector to triangle at
128 * Gauss pts.
129 */
131
132 /** \brief if higher order geometry return tangent vector to triangle at
133 * Gauss pts.
134 */
136
137 /** \brief return pointer to triangle finite element object
138 */
141
142 protected:
144
145 };
146
148
149protected:
150 double aRea[2];
154
156
166
167 friend class UserDataOperator;
168};
169
170/// \brief USe FlatPrismElementForcesAndSourcesCore
173
174inline double
176 return static_cast<FlatPrismElementForcesAndSourcesCore *>(ptrFE)->aRea[0];
177}
178
179inline double
181 return static_cast<FlatPrismElementForcesAndSourcesCore *>(ptrFE)->aRea[0];
182}
183inline double
185 return static_cast<FlatPrismElementForcesAndSourcesCore *>(ptrFE)->aRea[1];
186}
187
188inline VectorDouble &
190 return static_cast<FlatPrismElementForcesAndSourcesCore *>(ptrFE)->normal;
191}
192
193inline VectorAdaptor
195 double *data =
196 &(static_cast<FlatPrismElementForcesAndSourcesCore *>(ptrFE)->normal[0]);
197 return VectorAdaptor(3, ublas::shallow_array_adaptor<double>(3, data));
198}
199
200inline VectorAdaptor
202 double *data =
203 &(static_cast<FlatPrismElementForcesAndSourcesCore *>(ptrFE)->normal[3]);
204 return VectorAdaptor(3, ublas::shallow_array_adaptor<double>(3, data));
205}
206
207inline VectorDouble &
209 return static_cast<FlatPrismElementForcesAndSourcesCore *>(ptrFE)->coords;
210}
211
212inline MatrixDouble &
214 return static_cast<FlatPrismElementForcesAndSourcesCore *>(ptrFE)
216}
217
220 return static_cast<FlatPrismElementForcesAndSourcesCore *>(ptrFE)
222}
223
226 return static_cast<FlatPrismElementForcesAndSourcesCore *>(ptrFE)
228}
229
232 return static_cast<FlatPrismElementForcesAndSourcesCore *>(ptrFE)
234}
235
238 return static_cast<FlatPrismElementForcesAndSourcesCore *>(ptrFE)
240}
241
242inline ublas::matrix_row<MatrixDouble>
244 const int gg) {
245 return ublas::matrix_row<MatrixDouble>(
246 static_cast<FlatPrismElementForcesAndSourcesCore *>(ptrFE)
248 gg);
249}
250
251inline ublas::matrix_row<MatrixDouble>
253 const int gg) {
254 return ublas::matrix_row<MatrixDouble>(
255 static_cast<FlatPrismElementForcesAndSourcesCore *>(ptrFE)
257 gg);
258}
259
262 return static_cast<FlatPrismElementForcesAndSourcesCore *>(ptrFE)
264}
265
268 return static_cast<FlatPrismElementForcesAndSourcesCore *>(ptrFE)
270}
271
274 return static_cast<FlatPrismElementForcesAndSourcesCore *>(ptrFE)
276}
277
280 return static_cast<FlatPrismElementForcesAndSourcesCore *>(ptrFE)
282}
283
287 return static_cast<FlatPrismElementForcesAndSourcesCore *>(ptrFE);
288}
289
290} // namespace MoFEM
291
292#endif //__FLATPRISMELEMENTFORCESANDSURCESCORE_HPP__
293
294/**
295 * \defgroup mofem_forces_and_sources_prism_element Prism Element
296 * \ingroup mofem_forces_and_sources
297 **/
#define DEPRECATED
Definition: definitions.h:17
const Tensor2_symmetric_Expr< const ddTensor0< T, Dim, i, j >, typename promote< T, double >::V, Dim, i, j > dd(const Tensor0< T * > &a, const Index< i, Dim > index1, const Index< j, Dim > index2, const Tensor1< int, Dim > &d_ijk, const Tensor1< double, Dim > &d_xyz)
Definition: ddTensor0.hpp:33
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
VectorShallowArrayAdaptor< double > VectorAdaptor
Definition: Types.hpp:115
implementation of Data Operators for Forces and Sources
Definition: MoFEM.hpp:24
DEPRECATED typedef FlatPrismElementForcesAndSourcesCore FlatPrismElementForcesAndSurcesCore
USe FlatPrismElementForcesAndSourcesCore.
Deprecated interface functions.
MatrixDouble & getTangent2AtGaussPtF3()
if higher order geometry return tangent vector to triangle at Gauss pts.
MatrixDouble & getTangent2AtGaussPtF4()
if higher order geometry return tangent vector to triangle at Gauss pts.
MatrixDouble & getHOCoordsAtGaussPtsF3()
coordinate at Gauss points on face 3 (if hierarchical approximation of element geometry)
MatrixDouble & getTangent1AtGaussPtF3()
if higher order geometry return tangent vector to triangle at Gauss pts.
MatrixDouble & getHOCoordsAtGaussPtsF4()
coordinate at Gauss points on face 4 (if hierarchical approximation of element geometry)
MatrixDouble & getNormalsAtGaussPtsF4()
if higher order geometry return normals at face F4 at Gauss pts.
const FlatPrismElementForcesAndSourcesCore * getFlatPrismElementForcesAndSourcesCore()
return pointer to triangle finite element object
MatrixDouble & getTangent1AtGaussPtF4()
if higher order geometry return tangent vector to triangle at Gauss pts.
MatrixDouble & getNormalsAtGaussPtsF3()
if higher order geometry return normals at face F3 at Gauss pts.
MoFEMErrorCode operator()()
function is run for every finite element
UserDataOperator(const FieldSpace space, const char type=OPSPACE, const bool symm=true)
structure to get information form mofem into EntitiesFieldData
calculate normals at Gauss points of triangle element