v0.10.0
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 /* This file is part of MoFEM.
10  * MoFEM is free software: you can redistribute it and/or modify it under
11  * the terms of the GNU Lesser General Public License as published by the
12  * Free Software Foundation, either version 3 of the License, or (at your
13  * option) any later version.
14  *
15  * MoFEM is distributed in the hope that it will be useful, but WITHOUT
16  * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
17  * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
18  * License for more details.
19  *
20  * You should have received a copy of the GNU Lesser General Public
21  * License along with MoFEM. If not, see <http://www.gnu.org/licenses/>. */
22 
23 using namespace boost::numeric;
24 
25 #ifndef __FLATPRISMELEMENTFORCESANDSURCESCORE_HPP__
26 #define __FLATPRISMELEMENTFORCESANDSURCESCORE_HPP__
27 
28 namespace MoFEM {
29 
30 /** \brief FlatPrism finite element
31  \ingroup mofem_forces_and_sources_prism_element
32 
33  User is implementing own operator at Gauss points level, by own object
34  derived from FlatPrismElementForcesAndSourcesCoreL::UserDataOperator. Arbitrary
35  number of operator added pushing objects to rowOpPtrVector and
36  rowColOpPtrVector.
37 
38  */
40 
42 
43  /** \brief default operator for Flat Prism element
44  * \ingroup mofem_forces_and_sources_prism_element
45  */
47 
49 
50  /** \brief get face aRea
51  \param dd if dd == 0 it is for face F3 if dd == 1 is for face F4
52  */
53  inline double getArea(const int dd);
54 
55  inline double getAreaF3();
56 
57  inline double getAreaF4();
58 
59  /** \brief get triangle normal
60 
61  Normal has 6 elements, first 3 are for face F3 another three for face F4
62 
63  */
64  inline VectorDouble &getNormal();
65 
66  inline VectorAdaptor getNormalF3();
67 
68  inline VectorAdaptor getNormalF4();
69 
70  /** \brief get triangle coordinates
71 
72  Vector has 6 elements, i.e. coordinates on face F3 and F4
73 
74  */
75  inline VectorDouble &getCoords();
76 
77  /** \brief get coordinates at Gauss pts.
78 
79  Matrix has size (nb integration points)x(coordinates on F3 and F4 = 6),
80  i.e. coordinates on face F3 and F4
81 
82  */
83  inline MatrixDouble &getCoordsAtGaussPts();
84 
85  /** \brief coordinate at Gauss points on face 3 (if hierarchical
86  * approximation of element geometry)
87  */
88  inline MatrixDouble &getHoCoordsAtGaussPtsF3();
89 
90  /** \brief coordinate at Gauss points on face 4 (if hierarchical
91  * approximation of element geometry)
92  */
93  inline MatrixDouble &getHoCoordsAtGaussPtsF4();
94 
95  /** \brief if higher order geometry return normals at face F3 at Gauss pts.
96  *
97  * Face 3 is top face in canonical triangle numeration, see \cite
98  * tautges2010canonical
99  *
100  */
101  inline MatrixDouble &getNormalsAtGaussPtF3();
102 
103  /** \brief if higher order geometry return normals at face F4 at Gauss pts.
104  *
105  * Face 4 is top face in canonical triangle numeration, see \cite
106  * tautges2010canonical
107  *
108  */
109  inline MatrixDouble &getNormalsAtGaussPtF4();
110 
111  /** \brief if higher order geometry return normals at Gauss pts.
112  *
113  * Face 3 is top face in canonical triangle numeration, see \cite
114  * tautges2010canonical
115  *
116  * \param gg gauss point number
117  */
118  inline ublas::matrix_row<MatrixDouble> getNormalsAtGaussPtF3(const int gg);
119 
120  /** \brief if higher order geometry return normals at Gauss pts.
121  *
122  * Face 3 is top face in canonical triangle numeration, see \cite
123  * tautges2010canonical
124  *
125  * \param gg gauss point number
126  */
127  inline ublas::matrix_row<MatrixDouble> getNormalsAtGaussPtF4(const int gg);
128 
129  /** \brief if higher order geometry return tangent vector to triangle at
130  * Gauss pts.
131  */
132  inline MatrixDouble &getTangent1AtGaussPtF3();
133 
134  /** \brief if higher order geometry return tangent vector to triangle at
135  * Gauss pts.
136  */
137  inline MatrixDouble &getTangent2AtGaussPtF3();
138 
139  /** \brief if higher order geometry return tangent vector to triangle at
140  * Gauss pts.
141  */
142  inline MatrixDouble &getTangent1AtGaussPtF4();
143 
144  /** \brief if higher order geometry return tangent vector to triangle at
145  * Gauss pts.
146  */
147  inline MatrixDouble &getTangent2AtGaussPtF4();
148 
149  /** \brief return pointer to triangle finite element object
150  */
152  getFlatPrismElementForcesAndSourcesCore();
153 
154  protected:
155  MoFEMErrorCode setPtrFE(ForcesAndSourcesCore *ptr);
156 
157  };
158 
159  MoFEMErrorCode operator()();
160 
161 protected:
162  double aRea[2];
166 
168 
178 
179  friend class UserDataOperator;
180 };
181 
182 /// \brief USe FlatPrismElementForcesAndSourcesCore
185 
186 inline double
187 FlatPrismElementForcesAndSourcesCore::UserDataOperator::getArea(const int dd) {
188  return static_cast<FlatPrismElementForcesAndSourcesCore *>(ptrFE)->aRea[0];
189 }
190 
191 inline double
192 FlatPrismElementForcesAndSourcesCore::UserDataOperator::getAreaF3() {
193  return static_cast<FlatPrismElementForcesAndSourcesCore *>(ptrFE)->aRea[0];
194 }
195 inline double
196 FlatPrismElementForcesAndSourcesCore::UserDataOperator::getAreaF4() {
197  return static_cast<FlatPrismElementForcesAndSourcesCore *>(ptrFE)->aRea[1];
198 }
199 
200 inline VectorDouble &
201 FlatPrismElementForcesAndSourcesCore::UserDataOperator::getNormal() {
202  return static_cast<FlatPrismElementForcesAndSourcesCore *>(ptrFE)->normal;
203 }
204 
205 inline VectorAdaptor
206 FlatPrismElementForcesAndSourcesCore::UserDataOperator::getNormalF3() {
207  double *data =
208  &(static_cast<FlatPrismElementForcesAndSourcesCore *>(ptrFE)->normal[0]);
209  return VectorAdaptor(3, ublas::shallow_array_adaptor<double>(3, data));
210 }
211 
212 inline VectorAdaptor
213 FlatPrismElementForcesAndSourcesCore::UserDataOperator::getNormalF4() {
214  double *data =
215  &(static_cast<FlatPrismElementForcesAndSourcesCore *>(ptrFE)->normal[3]);
216  return VectorAdaptor(3, ublas::shallow_array_adaptor<double>(3, data));
217 }
218 
219 inline VectorDouble &
220 FlatPrismElementForcesAndSourcesCore::UserDataOperator::getCoords() {
221  return static_cast<FlatPrismElementForcesAndSourcesCore *>(ptrFE)->coords;
222 }
223 
224 inline MatrixDouble &
225 FlatPrismElementForcesAndSourcesCore::UserDataOperator::getCoordsAtGaussPts() {
226  return static_cast<FlatPrismElementForcesAndSourcesCore *>(ptrFE)
227  ->coordsAtGaussPts;
228 }
229 
230 inline MatrixDouble &FlatPrismElementForcesAndSourcesCore::UserDataOperator::
231  getHoCoordsAtGaussPtsF3() {
232  return static_cast<FlatPrismElementForcesAndSourcesCore *>(ptrFE)
233  ->hoCoordsAtGaussPtsF3;
234 }
235 
236 inline MatrixDouble &FlatPrismElementForcesAndSourcesCore::UserDataOperator::
237  getHoCoordsAtGaussPtsF4() {
238  return static_cast<FlatPrismElementForcesAndSourcesCore *>(ptrFE)
239  ->hoCoordsAtGaussPtsF4;
240 }
241 
242 inline MatrixDouble &FlatPrismElementForcesAndSourcesCore::UserDataOperator::
243  getNormalsAtGaussPtF3() {
244  return static_cast<FlatPrismElementForcesAndSourcesCore *>(ptrFE)
245  ->nOrmals_at_GaussPtF3;
246 }
247 
248 inline MatrixDouble &FlatPrismElementForcesAndSourcesCore::UserDataOperator::
249  getNormalsAtGaussPtF4() {
250  return static_cast<FlatPrismElementForcesAndSourcesCore *>(ptrFE)
251  ->nOrmals_at_GaussPtF4;
252 }
253 
254 inline ublas::matrix_row<MatrixDouble>
255 FlatPrismElementForcesAndSourcesCore::UserDataOperator::getNormalsAtGaussPtF3(
256  const int gg) {
257  return ublas::matrix_row<MatrixDouble>(
258  static_cast<FlatPrismElementForcesAndSourcesCore *>(ptrFE)
259  ->nOrmals_at_GaussPtF3,
260  gg);
261 }
262 
263 inline ublas::matrix_row<MatrixDouble>
264 FlatPrismElementForcesAndSourcesCore::UserDataOperator::getNormalsAtGaussPtF4(
265  const int gg) {
266  return ublas::matrix_row<MatrixDouble>(
267  static_cast<FlatPrismElementForcesAndSourcesCore *>(ptrFE)
268  ->nOrmals_at_GaussPtF4,
269  gg);
270 }
271 
272 inline MatrixDouble &FlatPrismElementForcesAndSourcesCore::UserDataOperator::
273  getTangent1AtGaussPtF3() {
274  return static_cast<FlatPrismElementForcesAndSourcesCore *>(ptrFE)
275  ->tAngent1_at_GaussPtF3;
276 }
277 
278 inline MatrixDouble &FlatPrismElementForcesAndSourcesCore::UserDataOperator::
279  getTangent2AtGaussPtF3() {
280  return static_cast<FlatPrismElementForcesAndSourcesCore *>(ptrFE)
281  ->tAngent2_at_GaussPtF3;
282 }
283 
284 inline MatrixDouble &FlatPrismElementForcesAndSourcesCore::UserDataOperator::
285  getTangent1AtGaussPtF4() {
286  return static_cast<FlatPrismElementForcesAndSourcesCore *>(ptrFE)
287  ->tAngent1_at_GaussPtF4;
288 }
289 
290 inline MatrixDouble &FlatPrismElementForcesAndSourcesCore::UserDataOperator::
291  getTangent2AtGaussPtF4() {
292  return static_cast<FlatPrismElementForcesAndSourcesCore *>(ptrFE)
293  ->tAngent2_at_GaussPtF4;
294 }
295 
297 FlatPrismElementForcesAndSourcesCore::UserDataOperator::
298  getFlatPrismElementForcesAndSourcesCore() {
299  return static_cast<FlatPrismElementForcesAndSourcesCore *>(ptrFE);
300 }
301 
302 } // namespace MoFEM
303 
304 #endif //__FLATPRISMELEMENTFORCESANDSURCESCORE_HPP__
305 
306 /**
307  * \defgroup mofem_forces_and_sources_prism_element Prism Element
308  * \ingroup mofem_forces_and_sources
309  **/
Deprecated interface functions.
ublas::matrix< double, ublas::row_major, DoubleAllocator > MatrixDouble
Definition: Types.hpp:76
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
DEPRECATED typedef FlatPrismElementForcesAndSourcesCore FlatPrismElementForcesAndSurcesCore
USe FlatPrismElementForcesAndSourcesCore.
ForcesAndSourcesCore::UserDataOperator UserDataOperator
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:109
calculate normals at Gauss points of triangle element
FlatPrism finite elementUser is implementing own operator at Gauss points level, by own object derive...
#define DEPRECATED
Definition: definitions.h:29
structure to get information form mofem into DataForcesAndSourcesCore
ublas::vector< double, DoubleAllocator > VectorDouble
Definition: Types.hpp:74
Data operator to do calculations at integration points.Is inherited and implemented by user to do cal...