v0.9.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 
27 using namespace boost::numeric;
28 
29 namespace 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 
86  inline VectorAdaptor getNormalF3();
87 
88  inline VectorAdaptor getNormalF4();
89 
90  /** \brief get Gauss pts. in the prism
91  */
92  inline MatrixDouble &getGaussPts();
93 
94  /** \brief get Gauss pts. on triangles
95  */
96  inline MatrixDouble &getGaussPtsTrianglesOnly();
97 
98  /** \brief get Gauss pts. through thickness
99  */
100  inline MatrixDouble &getGaussPtsThroughThickness();
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  */
108  inline MatrixDouble &getCoordsAtGaussPts();
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  */
116  inline MatrixDouble &getCoordsAtGaussPtsTrianglesOnly();
117 
118  /** \brief coordinate at Gauss points on face 3 (if hierarchical
119  * approximation of element geometry)
120  */
121  inline MatrixDouble &getHoCoordsAtGaussPtsF3();
122 
123  /** \brief coordinate at Gauss points on face 4 (if hierarchical
124  * approximation of element geometry)
125  */
126  inline MatrixDouble &getHoCoordsAtGaussPtsF4();
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  */
134  inline MatrixDouble &getNormalsAtGaussPtF3();
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  */
142  inline MatrixDouble &getNormalsAtGaussPtF4();
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> getNormalsAtGaussPtF3(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> getNormalsAtGaussPtF4(const int gg);
161 
162  /** \brief if higher order geometry return tangent vector to triangle at
163  * Gauss pts.
164  */
165  inline MatrixDouble &getTangent1AtGaussPtF3();
166 
167  /** \brief if higher order geometry return tangent vector to triangle at
168  * Gauss pts.
169  */
170  inline MatrixDouble &getTangent2AtGaussPtF3();
171 
172  /** \brief if higher order geometry return tangent vector to triangle at
173  * Gauss pts.
174  */
175  inline MatrixDouble &getTangent1AtGaussPtF4();
176 
177  /** \brief if higher order geometry return tangent vector to triangle at
178  * Gauss pts.
179  */
180  inline MatrixDouble &getTangent2AtGaussPtF4();
181 
182  inline DataForcesAndSourcesCore &getTrianglesOnlyDataStructure();
183 
184  inline DataForcesAndSourcesCore &getTroughThicknessDataStructure();
185 
186  /** \brief return pointer to fat prism finite element
187  */
188  inline const FatPrismElementForcesAndSourcesCore *getPrismFE();
189 
190  protected:
191  MoFEMErrorCode setPtrFE(ForcesAndSourcesCore *ptr);
192 
193  };
194 
195  MoFEMErrorCode operator()();
196 
197 protected:
198  double aRea[2];
200 
204 
207 
217 
218  friend class UserDataOperator;
219 };
220 
221 inline double
222 FatPrismElementForcesAndSourcesCore::UserDataOperator::getArea(const int dd) {
223  return static_cast<FatPrismElementForcesAndSourcesCore *>(ptrFE)->aRea[0];
224 }
225 
226 inline double
227 FatPrismElementForcesAndSourcesCore::UserDataOperator::getAreaF3() {
228  return static_cast<FatPrismElementForcesAndSourcesCore *>(ptrFE)->aRea[0];
229 }
230 inline double
231 FatPrismElementForcesAndSourcesCore::UserDataOperator::getAreaF4() {
232  return static_cast<FatPrismElementForcesAndSourcesCore *>(ptrFE)->aRea[1];
233 }
234 
235 inline VectorDouble &
236 FatPrismElementForcesAndSourcesCore::UserDataOperator::getNormal() {
237  return static_cast<FatPrismElementForcesAndSourcesCore *>(ptrFE)->normal;
238 }
239 
240 inline VectorAdaptor
241 FatPrismElementForcesAndSourcesCore::UserDataOperator::getNormalF3() {
242  double *data =
243  &(static_cast<FatPrismElementForcesAndSourcesCore *>(ptrFE)->normal[0]);
244  return VectorAdaptor(3, ublas::shallow_array_adaptor<double>(3, data));
245 }
246 
247 inline VectorAdaptor
248 FatPrismElementForcesAndSourcesCore::UserDataOperator::getNormalF4() {
249  double *data =
250  &(static_cast<FatPrismElementForcesAndSourcesCore *>(ptrFE)->normal[3]);
251  return VectorAdaptor(3, ublas::shallow_array_adaptor<double>(3, data));
252 }
253 
254 inline MatrixDouble &
255 FatPrismElementForcesAndSourcesCore::UserDataOperator::getGaussPts() {
256  return static_cast<FatPrismElementForcesAndSourcesCore *>(ptrFE)->gaussPts;
257 }
258 
259 inline MatrixDouble &FatPrismElementForcesAndSourcesCore::UserDataOperator::
260  getGaussPtsTrianglesOnly() {
261  return static_cast<FatPrismElementForcesAndSourcesCore *>(ptrFE)
262  ->gaussPtsTrianglesOnly;
263 }
264 
265 inline MatrixDouble &FatPrismElementForcesAndSourcesCore::UserDataOperator::
266  getGaussPtsThroughThickness() {
267  return static_cast<FatPrismElementForcesAndSourcesCore *>(ptrFE)
268  ->gaussPtsThroughThickness;
269 }
270 
271 inline MatrixDouble &
272 FatPrismElementForcesAndSourcesCore::UserDataOperator::getCoordsAtGaussPts() {
273  return static_cast<FatPrismElementForcesAndSourcesCore *>(ptrFE)
274  ->coordsAtGaussPts;
275 }
276 
277 inline MatrixDouble &FatPrismElementForcesAndSourcesCore::UserDataOperator::
278  getCoordsAtGaussPtsTrianglesOnly() {
279  return static_cast<FatPrismElementForcesAndSourcesCore *>(ptrFE)
280  ->coordsAtGaussPtsTrianglesOnly;
281 }
282 
283 inline MatrixDouble &FatPrismElementForcesAndSourcesCore::UserDataOperator::
284  getHoCoordsAtGaussPtsF3() {
285  return static_cast<FatPrismElementForcesAndSourcesCore *>(ptrFE)
286  ->hoCoordsAtGaussPtsF3;
287 }
288 
289 inline MatrixDouble &FatPrismElementForcesAndSourcesCore::UserDataOperator::
290  getHoCoordsAtGaussPtsF4() {
291  return static_cast<FatPrismElementForcesAndSourcesCore *>(ptrFE)
292  ->hoCoordsAtGaussPtsF4;
293 }
294 
295 inline MatrixDouble &
296 FatPrismElementForcesAndSourcesCore::UserDataOperator::getNormalsAtGaussPtF3() {
297  return static_cast<FatPrismElementForcesAndSourcesCore *>(ptrFE)
298  ->nOrmals_at_GaussPtF3;
299 }
300 
301 inline MatrixDouble &
302 FatPrismElementForcesAndSourcesCore::UserDataOperator::getNormalsAtGaussPtF4() {
303  return static_cast<FatPrismElementForcesAndSourcesCore *>(ptrFE)
304  ->nOrmals_at_GaussPtF4;
305 }
306 
307 inline ublas::matrix_row<MatrixDouble>
308 FatPrismElementForcesAndSourcesCore::UserDataOperator::getNormalsAtGaussPtF3(
309  const int gg) {
310  return ublas::matrix_row<MatrixDouble>(
311  static_cast<FatPrismElementForcesAndSourcesCore *>(ptrFE)
312  ->nOrmals_at_GaussPtF3,
313  gg);
314 }
315 
316 inline ublas::matrix_row<MatrixDouble>
317 FatPrismElementForcesAndSourcesCore::UserDataOperator::getNormalsAtGaussPtF4(
318  const int gg) {
319  return ublas::matrix_row<MatrixDouble>(
320  static_cast<FatPrismElementForcesAndSourcesCore *>(ptrFE)
321  ->nOrmals_at_GaussPtF4,
322  gg);
323 }
324 
325 inline MatrixDouble &FatPrismElementForcesAndSourcesCore::UserDataOperator::
326  getTangent1AtGaussPtF3() {
327  return static_cast<FatPrismElementForcesAndSourcesCore *>(ptrFE)
328  ->tAngent1_at_GaussPtF3;
329 }
330 
331 inline MatrixDouble &FatPrismElementForcesAndSourcesCore::UserDataOperator::
332  getTangent2AtGaussPtF3() {
333  return static_cast<FatPrismElementForcesAndSourcesCore *>(ptrFE)
334  ->tAngent2_at_GaussPtF3;
335 }
336 
337 inline MatrixDouble &FatPrismElementForcesAndSourcesCore::UserDataOperator::
338  getTangent1AtGaussPtF4() {
339  return static_cast<FatPrismElementForcesAndSourcesCore *>(ptrFE)
340  ->tAngent1_at_GaussPtF4;
341 }
342 
343 inline MatrixDouble &FatPrismElementForcesAndSourcesCore::UserDataOperator::
344  getTangent2AtGaussPtF4() {
345  return static_cast<FatPrismElementForcesAndSourcesCore *>(ptrFE)
346  ->tAngent2_at_GaussPtF4;
347 }
348 
349 inline DataForcesAndSourcesCore &FatPrismElementForcesAndSourcesCore::
350  UserDataOperator::getTrianglesOnlyDataStructure() {
351  return static_cast<FatPrismElementForcesAndSourcesCore *>(ptrFE)
352  ->dataH1TrianglesOnly;
353 }
354 
355 inline DataForcesAndSourcesCore &FatPrismElementForcesAndSourcesCore::
356  UserDataOperator::getTroughThicknessDataStructure() {
357  return static_cast<FatPrismElementForcesAndSourcesCore *>(ptrFE)
358  ->dataH1TroughThickness;
359 }
360 
362 FatPrismElementForcesAndSourcesCore::UserDataOperator::getPrismFE() {
363  return static_cast<FatPrismElementForcesAndSourcesCore *>(ptrFE);
364 }
365 
366 } // namespace MoFEM
367 
368 #endif //__FATPRISMELEMENTFORCESANDSURCESCORE_HPP__
Deprecated interface functions.
data structure for finite element entityIt keeps that about indices of degrees of freedom,...
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:506
ublas::matrix< double, ublas::row_major, DoubleAllocator > MatrixDouble
Definition: Types.hpp:74
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:513
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
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
constexpr int order
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:66
VectorShallowArrayAdaptor< double > VectorAdaptor
Definition: Types.hpp:104
calculate normals at Gauss points of triangle element
virtual MoFEMErrorCode setGaussPtsThroughThickness(int order_thickness)
ForcesAndSourcesCore::UserDataOperator UserDataOperator
virtual MoFEMErrorCode setGaussPtsTrianglesOnly(int order_triangles_only)
structure to get information form mofem into DataForcesAndSourcesCore
FatPrism finite elementUser is implementing own operator at Gauss points level, by own object derived...
ublas::vector< double, DoubleAllocator > VectorDouble
Definition: Types.hpp:72