v0.13.1
FatPrismElementForcesAndSourcesCore.hpp
Go to the documentation of this file.
1/** \file ForcesAndSourcesCore.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/* This file is part of MoFEM.
11 * MoFEM is free software: you can redistribute it and/or modify it under
12 * the terms of the GNU Lesser General Public License as published by the
13 * Free Software Foundation, either version 3 of the License, or (at your
14 * option) any later version.
15 *
16 * MoFEM is distributed in the hope that it will be useful, but WITHOUT
17 * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18 * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
19 * License for more details.
20 *
21 * You should have received a copy of the GNU Lesser General Public
22 * License along with MoFEM. If not, see <http://www.gnu.org/licenses/>. */
23
24#ifndef __FATPRISMELEMENTFORCESANDSURCESCORE_HPP__
25#define __FATPRISMELEMENTFORCESANDSURCESCORE_HPP__
26
27using namespace boost::numeric;
28
29namespace MoFEM {
30
31/** \brief FatPrism finite element
32 \ingroup mofem_forces_and_sources_prism_element
33
34 User is implementing own operator at Gauss points level, by own object
35 derived from FatPrismElementForcesAndSourcesCoreL::UserDataOperator. Arbitrary
36 number of operator added pushing objects to rowOpPtrVector and
37 rowColOpPtrVector.
38
39 \todo Need to implement operators that will make this element work as Volume
40 element
41
42 */
45
47
48 virtual int getRuleTrianglesOnly(int order) { return 2 * order; };
49 virtual int getRuleThroughThickness(int order) { return 2 * order; };
50
51 virtual MoFEMErrorCode setGaussPtsTrianglesOnly(int order_triangles_only) {
53 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
55 }
56
57 virtual MoFEMErrorCode setGaussPtsThroughThickness(int order_thickness) {
59 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
61 }
62
63 /** \brief default operator for Flat Prism element
64 * \ingroup mofem_forces_and_sources_prism_element
65 */
68
70
71 /** \brief get face aRea
72 \param dd if dd == 0 it is for face F3 if dd == 1 is for face F4
73 */
74 inline double getArea(const int dd);
75
76 inline double getAreaF3();
77 inline double getAreaF4();
78
79 /** \brief get triangle normal
80
81 Normal has 6 elements, first 3 are for face F3 another three for face F4
82
83 */
84 inline VectorDouble &getNormal();
85
87
89
90 /** \brief get Gauss pts. in the prism
91 */
92 inline MatrixDouble &getGaussPts();
93
94 /** \brief get Gauss pts. on triangles
95 */
97
98 /** \brief get Gauss pts. through thickness
99 */
101
102 /** \brief get coordinates at Gauss pts.
103
104 Matrix has size (nb integration points)x(coordinates on F3 and F4 = 6),
105 i.e. coordinates on face F3 and F4
106
107 */
109
110 /** \brief get coordinates at Gauss pts.
111
112 Matrix has size (nb integration points)x(coordinates on F3 and F4 = 6),
113 i.e. coordinates on face F3 and F4
114
115 */
117
118 /** \brief coordinate at Gauss points on face 3 (if hierarchical
119 * approximation of element geometry)
120 */
122
123 /** \brief coordinate at Gauss points on face 4 (if hierarchical
124 * approximation of element geometry)
125 */
127
128 /** \brief if higher order geometry return normals at face F3 at Gauss pts.
129 *
130 * Face 3 is top face in canonical triangle numeration, see \cite
131 * tautges2010canonical
132 *
133 */
135
136 /** \brief if higher order geometry return normals at face F4 at Gauss pts.
137 *
138 * Face 4 is top face in canonical triangle numeration, see \cite
139 * tautges2010canonical
140 *
141 */
143
144 /** \brief if higher order geometry return normals at Gauss pts.
145 *
146 * Face 3 is top face in canonical triangle numeration, see \cite
147 * tautges2010canonical
148 *
149 * \param gg gauss point number
150 */
151 inline ublas::matrix_row<MatrixDouble> getNormalsAtGaussPtsF3(const int gg);
152
153 /** \brief if higher order geometry return normals at Gauss pts.
154 *
155 * Face 3 is top face in canonical triangle numeration, see \cite
156 * tautges2010canonical
157 *
158 * \param gg gauss point number
159 */
160 inline ublas::matrix_row<MatrixDouble> getNormalsAtGaussPtsF4(const int gg);
161
162 /** \brief if higher order geometry return tangent vector to triangle at
163 * Gauss pts.
164 */
166
167 /** \brief if higher order geometry return tangent vector to triangle at
168 * Gauss pts.
169 */
171
172 /** \brief if higher order geometry return tangent vector to triangle at
173 * Gauss pts.
174 */
176
177 /** \brief if higher order geometry return tangent vector to triangle at
178 * Gauss pts.
179 */
181
183
185
186 /** \brief return pointer to fat prism finite element
187 */
189
190 protected:
192
193 };
194
196
197protected:
198 double aRea[2];
200
204
207
217
218 friend class UserDataOperator;
219};
220
221inline double
223 return static_cast<FatPrismElementForcesAndSourcesCore *>(ptrFE)->aRea[0];
224}
225
226inline double
228 return static_cast<FatPrismElementForcesAndSourcesCore *>(ptrFE)->aRea[0];
229}
230inline double
232 return static_cast<FatPrismElementForcesAndSourcesCore *>(ptrFE)->aRea[1];
233}
234
235inline VectorDouble &
237 return static_cast<FatPrismElementForcesAndSourcesCore *>(ptrFE)->normal;
238}
239
240inline VectorAdaptor
242 double *data =
243 &(static_cast<FatPrismElementForcesAndSourcesCore *>(ptrFE)->normal[0]);
244 return VectorAdaptor(3, ublas::shallow_array_adaptor<double>(3, data));
245}
246
247inline VectorAdaptor
249 double *data =
250 &(static_cast<FatPrismElementForcesAndSourcesCore *>(ptrFE)->normal[3]);
251 return VectorAdaptor(3, ublas::shallow_array_adaptor<double>(3, data));
252}
253
254inline MatrixDouble &
256 return static_cast<FatPrismElementForcesAndSourcesCore *>(ptrFE)->gaussPts;
257}
258
261 return static_cast<FatPrismElementForcesAndSourcesCore *>(ptrFE)
263}
264
267 return static_cast<FatPrismElementForcesAndSourcesCore *>(ptrFE)
269}
270
271inline MatrixDouble &
273 return static_cast<FatPrismElementForcesAndSourcesCore *>(ptrFE)
275}
276
279 return static_cast<FatPrismElementForcesAndSourcesCore *>(ptrFE)
281}
282
285 return static_cast<FatPrismElementForcesAndSourcesCore *>(ptrFE)
287}
288
291 return static_cast<FatPrismElementForcesAndSourcesCore *>(ptrFE)
293}
294
295inline MatrixDouble &
297 return static_cast<FatPrismElementForcesAndSourcesCore *>(ptrFE)
299}
300
301inline MatrixDouble &
303 return static_cast<FatPrismElementForcesAndSourcesCore *>(ptrFE)
305}
306
307inline ublas::matrix_row<MatrixDouble>
309 const int gg) {
310 return ublas::matrix_row<MatrixDouble>(
311 static_cast<FatPrismElementForcesAndSourcesCore *>(ptrFE)
313 gg);
314}
315
316inline ublas::matrix_row<MatrixDouble>
318 const int gg) {
319 return ublas::matrix_row<MatrixDouble>(
320 static_cast<FatPrismElementForcesAndSourcesCore *>(ptrFE)
322 gg);
323}
324
327 return static_cast<FatPrismElementForcesAndSourcesCore *>(ptrFE)
329}
330
333 return static_cast<FatPrismElementForcesAndSourcesCore *>(ptrFE)
335}
336
339 return static_cast<FatPrismElementForcesAndSourcesCore *>(ptrFE)
341}
342
345 return static_cast<FatPrismElementForcesAndSourcesCore *>(ptrFE)
347}
348
351 return static_cast<FatPrismElementForcesAndSourcesCore *>(ptrFE)
353}
354
357 return static_cast<FatPrismElementForcesAndSourcesCore *>(ptrFE)
359}
360
363 return static_cast<FatPrismElementForcesAndSourcesCore *>(ptrFE);
364}
365
366} // namespace MoFEM
367
368#endif //__FATPRISMELEMENTFORCESANDSURCESCORE_HPP__
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:460
@ MOFEM_NOT_IMPLEMENTED
Definition: definitions.h:45
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:453
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:67
VectorShallowArrayAdaptor< double > VectorAdaptor
Definition: Types.hpp:126
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
Deprecated interface functions.
data structure for finite element entity
MatrixDouble & getGaussPtsThroughThickness()
get Gauss pts. through thickness
MatrixDouble & getTangent2AtGaussPtF3()
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 & getTangent2AtGaussPtF4()
if higher order geometry return tangent vector to triangle at Gauss pts.
MatrixDouble & getNormalsAtGaussPtsF4()
if higher order geometry return normals at face F4 at Gauss pts.
MatrixDouble & getNormalsAtGaussPtsF3()
if higher order geometry return normals at face F3 at Gauss pts.
const FatPrismElementForcesAndSourcesCore * getPrismFE()
return pointer to fat prism finite element
MatrixDouble & getHOCoordsAtGaussPtsF3()
coordinate at Gauss points on face 3 (if hierarchical approximation of element geometry)
MatrixDouble & getTangent1AtGaussPtF4()
if higher order geometry return tangent vector to triangle at Gauss pts.
MatrixDouble & getTangent1AtGaussPtF3()
if higher order geometry return tangent vector to triangle at Gauss pts.
virtual MoFEMErrorCode setGaussPtsTrianglesOnly(int order_triangles_only)
virtual MoFEMErrorCode setGaussPtsThroughThickness(int order_thickness)
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
MatrixDouble coordsAtGaussPts
coordinated at gauss points
MatrixDouble gaussPts
Matrix of integration points.
calculate normals at Gauss points of triangle element