v0.9.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  return static_cast<FlatPrismElementForcesAndSourcesCore *>(ptrFE)
55  ->aRea[0];
56  }
57 
58  inline double getAreaF3() {
59  return static_cast<FlatPrismElementForcesAndSourcesCore *>(ptrFE)
60  ->aRea[0];
61  }
62  inline double getAreaF4() {
63  return static_cast<FlatPrismElementForcesAndSourcesCore *>(ptrFE)
64  ->aRea[1];
65  }
66 
67  /** \brief get triangle normal
68 
69  Normal has 6 elements, first 3 are for face F3 another three for face F4
70 
71  */
72  inline VectorDouble &getNormal() {
73  return static_cast<FlatPrismElementForcesAndSourcesCore *>(ptrFE)->normal;
74  }
75 
77  double *data =
78  &(static_cast<FlatPrismElementForcesAndSourcesCore *>(ptrFE)
79  ->normal[0]);
80  return VectorAdaptor(3, ublas::shallow_array_adaptor<double>(3, data));
81  }
82 
84  double *data =
85  &(static_cast<FlatPrismElementForcesAndSourcesCore *>(ptrFE)
86  ->normal[3]);
87  return VectorAdaptor(3, ublas::shallow_array_adaptor<double>(3, data));
88  }
89 
90  /** \brief get triangle coordinates
91 
92  Vector has 6 elements, i.e. coordinates on face F3 and F4
93 
94  */
95  inline VectorDouble &getCoords() {
96  return static_cast<FlatPrismElementForcesAndSourcesCore *>(ptrFE)->coords;
97  }
98 
99  /** \brief get coordinates at Gauss pts.
100 
101  Matrix has size (nb integration points)x(coordinates on F3 and F4 = 6),
102  i.e. coordinates on face F3 and F4
103 
104  */
106  return static_cast<FlatPrismElementForcesAndSourcesCore *>(ptrFE)
107  ->coordsAtGaussPts;
108  }
109 
110  /** \brief coordinate at Gauss points on face 3 (if hierarchical
111  * approximation of element geometry)
112  */
114  return static_cast<FlatPrismElementForcesAndSourcesCore *>(ptrFE)
115  ->hoCoordsAtGaussPtsF3;
116  }
117 
118  /** \brief coordinate at Gauss points on face 4 (if hierarchical
119  * approximation of element geometry)
120  */
122  return static_cast<FlatPrismElementForcesAndSourcesCore *>(ptrFE)
123  ->hoCoordsAtGaussPtsF4;
124  }
125 
126  /** \brief if higher order geometry return normals at face F3 at Gauss pts.
127  *
128  * Face 3 is top face in canonical triangle numeration, see \cite
129  * tautges2010canonical
130  *
131  */
133  return static_cast<FlatPrismElementForcesAndSourcesCore *>(ptrFE)
134  ->nOrmals_at_GaussPtF3;
135  }
136 
137  /** \brief if higher order geometry return normals at face F4 at Gauss pts.
138  *
139  * Face 4 is top face in canonical triangle numeration, see \cite
140  * tautges2010canonical
141  *
142  */
144  return static_cast<FlatPrismElementForcesAndSourcesCore *>(ptrFE)
145  ->nOrmals_at_GaussPtF4;
146  }
147 
148  /** \brief if higher order geometry return normals at Gauss pts.
149  *
150  * Face 3 is top face in canonical triangle numeration, see \cite
151  * tautges2010canonical
152  *
153  * \param gg gauss point number
154  */
155  inline ublas::matrix_row<MatrixDouble> getNormalsAtGaussPtF3(const int gg) {
156  return ublas::matrix_row<MatrixDouble>(
157  static_cast<FlatPrismElementForcesAndSourcesCore *>(ptrFE)
158  ->nOrmals_at_GaussPtF3,
159  gg);
160  }
161 
162  /** \brief if higher order geometry return normals at Gauss pts.
163  *
164  * Face 3 is top face in canonical triangle numeration, see \cite
165  * tautges2010canonical
166  *
167  * \param gg gauss point number
168  */
169  inline ublas::matrix_row<MatrixDouble> getNormalsAtGaussPtF4(const int gg) {
170  return ublas::matrix_row<MatrixDouble>(
171  static_cast<FlatPrismElementForcesAndSourcesCore *>(ptrFE)
172  ->nOrmals_at_GaussPtF4,
173  gg);
174  }
175 
176  /** \brief if higher order geometry return tangent vector to triangle at
177  * Gauss pts.
178  */
180  return static_cast<FlatPrismElementForcesAndSourcesCore *>(ptrFE)
181  ->tAngent1_at_GaussPtF3;
182  }
183 
184  /** \brief if higher order geometry return tangent vector to triangle at
185  * Gauss pts.
186  */
188  return static_cast<FlatPrismElementForcesAndSourcesCore *>(ptrFE)
189  ->tAngent2_at_GaussPtF3;
190  }
191 
192  /** \brief if higher order geometry return tangent vector to triangle at
193  * Gauss pts.
194  */
196  return static_cast<FlatPrismElementForcesAndSourcesCore *>(ptrFE)
197  ->tAngent1_at_GaussPtF4;
198  }
199 
200  /** \brief if higher order geometry return tangent vector to triangle at
201  * Gauss pts.
202  */
204  return static_cast<FlatPrismElementForcesAndSourcesCore *>(ptrFE)
205  ->tAngent2_at_GaussPtF4;
206  }
207 
208  /** \brief return pointer to triangle finite element object
209  */
212  return static_cast<FlatPrismElementForcesAndSourcesCore *>(ptrFE);
213  }
214  };
215 
216  MoFEMErrorCode operator()();
217 
218  protected:
219 
220  double aRea[2];
224 
226 
236 
237  friend class UserDataOperator;
238 
239 };
240 
241 /// \brief USe FlatPrismElementForcesAndSourcesCore
244 
245 } // namespace MoFEM
246 
247 #endif //__FLATPRISMELEMENTFORCESANDSURCESCORE_HPP__
248 
249 /**
250  * \defgroup mofem_forces_and_sources_prism_element Prism Element
251  * \ingroup mofem_forces_and_sources
252  **/
MatrixDouble & getTangent2AtGaussPtF3()
if higher order geometry return tangent vector to triangle at Gauss pts.
MatrixDouble & getNormalsAtGaussPtF3()
if higher order geometry return normals at face F3 at Gauss pts.
ublas::matrix< double, ublas::row_major, DoubleAllocator > MatrixDouble
Definition: Types.hpp:74
MatrixDouble & getHoCoordsAtGaussPtsF3()
coordinate at Gauss points on face 3 (if hierarchical approximation of element geometry)
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
DEPRECATED typedef FlatPrismElementForcesAndSourcesCore FlatPrismElementForcesAndSurcesCore
USe FlatPrismElementForcesAndSourcesCore.
const FlatPrismElementForcesAndSourcesCore * getFlatPrismElementForcesAndSourcesCore()
return pointer to triangle finite element object
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
ublas::matrix_row< MatrixDouble > getNormalsAtGaussPtF3(const int gg)
if higher order geometry return normals at Gauss pts.
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:66
VectorShallowArrayAdaptor< double > VectorAdaptor
Definition: Types.hpp:101
calculate normals at Gauss points of triangle element
FlatPrism finite elementUser is implementing own operator at Gauss points level, by own object derive...
ForcesAndSourcesCore::UserDataOperator UserDataOperator
MatrixDouble & getTangent1AtGaussPtF3()
if higher order geometry return tangent vector to triangle at Gauss pts.
#define DEPRECATED
Definition: definitions.h:29
MatrixDouble & getTangent1AtGaussPtF4()
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)
structure to get information form mofem into DataForcesAndSourcesCore
MatrixDouble & getTangent2AtGaussPtF4()
if higher order geometry return tangent vector to triangle at Gauss pts.
MatrixDouble & getNormalsAtGaussPtF4()
if higher order geometry return normals at face F4 at Gauss pts.
ublas::vector< double, DoubleAllocator > VectorDouble
Definition: Types.hpp:72
Data operator to do calculations at integration points.Is inherited and implemented by user to do cal...
ublas::matrix_row< MatrixDouble > getNormalsAtGaussPtF4(const int gg)
if higher order geometry return normals at Gauss pts.