v0.13.2
Loading...
Searching...
No Matches
EdgeElementForcesAndSourcesCore.hpp
Go to the documentation of this file.
1/** \file EdgeElementForcesAndSourcesCore.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
11
12#ifndef __EDGEELEMENTFORCESANDSURCESCORE_HPP__
13#define __EDGEELEMENTFORCESANDSURCESCORE_HPP__
14
15using namespace boost::numeric;
16
17namespace MoFEM {
18
19struct FaceElementForcesAndSourcesCoreOnSide;
20
21/** \brief Edge finite element
22 * \ingroup mofem_forces_and_sources_edge_element
23 *
24 * User is implementing own operator at Gauss points level, by own object
25 * derived from EdgeElementForcesAndSourcesCoreL::UserDataOperator. Arbitrary
26 * number of operator added pushing objects to rowOpPtrVector and
27 * rowColOpPtrVector.
28 *
29 */
32
34
35 /** \brief default operator for EDGE element
36 \ingroup mofem_forces_and_sources_edge_element
37 */
38 struct UserDataOperator;
39
40 enum Switches {};
41
43
44protected:
45
46
48
49 double lEngth;
50
55
59
61};
62
63/** \brief default operator for EDGE element
64 \ingroup mofem_forces_and_sources_edge_element
65 */
68
69 using ForcesAndSourcesCore::UserDataOperator::UserDataOperator;
70
71 /** \brief get element connectivity
72 */
73 inline const EntityHandle *getConn();
74
75 /**
76 * \brief get edge length
77 */
78 inline double getLength();
79
80 /**
81 * \brief get measure of element
82 * @return length of face
83 */
84 inline double getMeasure();
85
86 /**
87 * \brief get edge direction
88 */
89 inline VectorDouble &getDirection();
90
91 /**
92 * \brief get edge normal
93 * NOTE: it should be used only in 2D analysis
94 */
95 inline auto getFTensor1Normal();
96
97 /**
98 * @brief get ftensor1 edge normal
99 *
100 * @param vec vector in third direction
101 * @return auto
102 */
103 inline auto getFTensor1Normal(const FTensor::Tensor1<double, 3> &vec);
104
105 /**
106 * \brief get edge node coordinates
107 */
108 inline VectorDouble &getCoords();
109
110 /**
111 * \brief get tangent vector to edge curve at integration points
112 */
114
116 return getTangentAtGaussPts();
117 }
118
119 /**
120 * \brief get pointer to this finite element
121 */
123
125
126 /**
127 * \brief get get coords at gauss points
128
129 \code
130 FTensor::Index<'i',3> i;
131 auto t_center;
132 auto t_coords = getFTensor1Coords();
133 t_center(i) = 0;
134 for(int nn = 0;nn!=2;nn++) {
135 t_center(i) += t_coords(i);
136 ++t_coords;
137 }
138 t_center(i) /= 2;
139 \endcode
140
141 */
143
144 /**
145 * @deprecated use getFTensor1Coords
146 */
149 return getFTensor1Coords();
150 }
151
152 /**
153 * @brief Get tangent vector at integration points
154 *
155 * @return FTensor::Tensor1<FTensor::PackPtr<double *, 3>, DIM>
156 */
157 template <int DIM = 3>
160
162 loopSideFaces(const string fe_name,
164
165protected:
167};
168
169const EntityHandle *
170EdgeElementForcesAndSourcesCore::UserDataOperator::getConn() {
171 return static_cast<EdgeElementForcesAndSourcesCore *>(ptrFE)->cOnn;
172}
173
174double EdgeElementForcesAndSourcesCore::UserDataOperator::getLength() {
175 return static_cast<EdgeElementForcesAndSourcesCore *>(ptrFE)->lEngth;
176}
177
178double EdgeElementForcesAndSourcesCore::UserDataOperator::getMeasure() {
179 return getLength();
180}
181
183EdgeElementForcesAndSourcesCore::UserDataOperator::getDirection() {
184 return static_cast<EdgeElementForcesAndSourcesCore *>(ptrFE)->dIrection;
185}
186
188EdgeElementForcesAndSourcesCore::UserDataOperator::getCoords() {
189 return static_cast<EdgeElementForcesAndSourcesCore *>(ptrFE)->cOords;
190}
191
193EdgeElementForcesAndSourcesCore::UserDataOperator::getTangentAtGaussPts() {
194 return static_cast<EdgeElementForcesAndSourcesCore *>(ptrFE)
196}
197
199EdgeElementForcesAndSourcesCore::UserDataOperator::getEdgeFE() {
200 return static_cast<EdgeElementForcesAndSourcesCore *>(ptrFE);
201}
202
204EdgeElementForcesAndSourcesCore::UserDataOperator::getFTensor1Direction() {
205 return FTensor::Tensor1<double, 3>(getDirection()[0], getDirection()[1],
206 getDirection()[2]);
207}
208
210EdgeElementForcesAndSourcesCore::UserDataOperator::getFTensor1Coords() {
211 double *ptr = &*getCoords().data().begin();
213 &ptr[2]);
214}
215
216template <>
218EdgeElementForcesAndSourcesCore::UserDataOperator::getFTensor1TangentAtGaussPts<
219 3>() {
220 double *ptr = &*getTangentAtGaussPts().data().begin();
222 &ptr[2]);
223}
224
225template <>
227EdgeElementForcesAndSourcesCore::UserDataOperator::getFTensor1TangentAtGaussPts<
228 2>() {
229 double *ptr = &*getTangentAtGaussPts().data().begin();
230 return FTensor::Tensor1<FTensor::PackPtr<double *, 3>, 2>(ptr, &ptr[1]);
231}
232
233auto EdgeElementForcesAndSourcesCore::UserDataOperator::getFTensor1Normal(
234 const FTensor::Tensor1<double, 3> &vec) {
239 auto t_dir = getFTensor1Direction();
240 t_normal(i) = FTensor::levi_civita(i, j, k) * t_dir(j) * vec(k);
241 return t_normal;
242}
243
244auto EdgeElementForcesAndSourcesCore::UserDataOperator::
245 getFTensor1Normal() {
246 FTensor::Tensor1<double, 3> t_normal{0., 0., 1.};
247 return getFTensor1Normal(t_normal);
248}
249
250/**
251 * @deprecated use EdgeElementForcesAndSourcesCore
252 */
255
256} // namespace MoFEM
257
258#endif //__EDGEELEMENTFORCESANDSURCESCORE_HPP__
259
260/**
261 * \defgroup mofem_forces_and_sources_edge_element Edge Element
262 *
263 * \brief Implementation of edge element.
264 *
265 * \ingroup mofem_forces_and_sources
266 */
T data[Tensor_Dim]
#define DEPRECATED
Definition: definitions.h:17
FTensor::Index< 'i', SPACE_DIM > i
FTensor::Index< 'j', 3 > j
FTensor::Index< 'k', 3 > k
constexpr std::enable_if<(Dim0<=2 &&Dim1<=2), Tensor2_Expr< Levi_Civita< T >, T, Dim0, Dim1, i, j > >::type levi_civita(const Index< i, Dim0 > &, const Index< j, Dim1 > &)
levi_civita functions to make for easy adhoc use
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:10
DEPRECATED typedef EdgeElementForcesAndSourcesCore EdgeElementForcesAndSourcesCoreBase
Deprecated interface functions.
DEPRECATED FTensor::Tensor1< FTensor::PackPtr< double *, 3 >, 3 > getTensor1Coords()
auto getFTensor1Normal()
get edge normal NOTE: it should be used only in 2D analysis
MatrixDouble & getTangentAtGaussPts()
get tangent vector to edge curve at integration points
MoFEMErrorCode loopSideFaces(const string fe_name, FaceElementForcesAndSourcesCoreOnSide &fe_side)
const EdgeElementForcesAndSourcesCore * getEdgeFE()
get pointer to this finite element
FTensor::Tensor1< FTensor::PackPtr< double *, 3 >, DIM > getFTensor1TangentAtGaussPts()
Get tangent vector at integration points.
FTensor::Tensor1< FTensor::PackPtr< double *, 3 >, 3 > getFTensor1Coords()
get get coords at gauss points
MoFEMErrorCode operator()()
function is run for every finite element
Base face element used to integrate on skeleton.
structure to get information form mofem into EntitiesFieldData