v0.13.1
VolumeElementForcesAndSourcesCore.hpp
Go to the documentation of this file.
1/** \file VolumeElementForcesAndSourcesCore.hpp
2 \brief Volume element.
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#ifndef __VOLUMEELEMENTFORCESANDSOURCESCORE_HPP__
24#define __VOLUMEELEMENTFORCESANDSOURCESCORE_HPP__
25
26using namespace boost::numeric;
27
28namespace MoFEM {
29
30/** \brief Volume finite element base
31 \ingroup mofem_forces_and_sources_volume_element
32
33 User is implementing own operator at Gauss point level, by class
34 derived from VolumeElementForcesAndSourcesCore::UserDataOperator. Arbitrary
35 number of operator can be added by pushing objects to OpPtrVector
36
37 */
39
41 const EntityType type = MBTET);
42
43
44 std::string meshPositionsFieldName; ///< \deprecated DO NOT USE!
45
46 /** \brief default operator for TET element
47 * \ingroup mofem_forces_and_sources_volume_element
48 */
49 struct UserDataOperator;
50
51 enum Switches {
52 };
53
55
56protected:
57
58 // Note that functions below could be overloaded by user to change default
59 // behavior of the element.
60
61 /**
62 * \brief Set integration points
63 * @return Error code
64 */
66
67 /**
68 * \brief Calculate element volume and Jacobian
69 *
70 * Note that at that point is assumed that geometry is exclusively defined by
71 * corner nodes.
72 *
73 * @return Error code
74 */
76
77 /**
78 * \brief Calculate coordinate at integration points
79 * @return Error code
80 */
82
83 /**
84 * \brief Determine approximation space and order of base functions
85 * @return Error code
86 */
88
89 /**
90 * \brief Transform base functions based on geometric element Jacobian.
91 *
92 * This function apply transformation to base functions and its derivatives.
93 * For example when base functions for H-div are present the
94 * Piola-Transformarion is applied to base functions and their derivatives.
95 *
96 * @return Error code
97 */
99
100
104
109
110 double &vOlume;
111
113 const EntityHandle *conn;
116
117 friend class UserDataOperator;
118};
119
122
124
125 /** \brief get element number of nodes
126 */
127 inline int getNumNodes();
128
129 /** \brief get element connectivity
130 */
131 inline const EntityHandle *getConn();
132
133 /** \brief element volume (linear geometry)
134 */
135 inline double getVolume() const;
136
137 /** \brief element volume (linear geometry)
138 */
139 inline double &getVolume();
140
141 /**
142 * \brief get element Jacobian
143 */
145
146 /**
147 * \brief get element inverse Jacobian
148 */
150
151 /** \brief nodal coordinates
152 */
153 inline VectorDouble &getCoords();
154
155 /** \brief return pointer to Generic Volume Finite Element object
156 */
158
159protected:
161};
162
164 return static_cast<VolumeElementForcesAndSourcesCore *>(ptrFE)->num_nodes;
165}
166
167const EntityHandle *
169 return static_cast<VolumeElementForcesAndSourcesCore *>(ptrFE)->conn;
170}
171
172double
174 return static_cast<VolumeElementForcesAndSourcesCore *>(ptrFE)->vOlume;
175}
176
178 return static_cast<VolumeElementForcesAndSourcesCore *>(ptrFE)->vOlume;
179}
180
183 return static_cast<VolumeElementForcesAndSourcesCore *>(ptrFE)->tJac;
184}
185
188 return static_cast<VolumeElementForcesAndSourcesCore *>(ptrFE)->tInvJac;
189}
190
193 return static_cast<VolumeElementForcesAndSourcesCore *>(ptrFE)->coords;
194}
195
198 return static_cast<VolumeElementForcesAndSourcesCore *>(ptrFE);
199}
200
201/**
202 * @deprecated use VolumeElementForcesAndSourcesCore
203 *
204 */
207
208/**
209 * @deprecated obsolete template. it will be removed in next versions
210 */
211template <int SWITCH>
214
217};
218
219} // namespace MoFEM
220
221#endif //__VOLUMEELEMENTFORCESANDSOURCESCORE_HPP__
222
223/**
224 * \defgroup mofem_forces_and_sources_volume_element Volume Element
225 * \brief Implementation of general volume element.
226 *
227 * \ingroup mofem_forces_and_sources
228 **/
#define DEPRECATED
Definition: definitions.h:30
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:67
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
DEPRECATED typedef VolumeElementForcesAndSourcesCore VolumeElementForcesAndSourcesCoreBase
Deprecated interface functions.
UserDataOperator(const FieldSpace space, const char type=OPSPACE, const bool symm=true)
structure to get information form mofem into EntitiesFieldData
apply contravariant (Piola) transfer to Hdiv space
apply covariant transfer to Hcurl space
Transform local reference derivatives of shape function to global derivatives.
brief Transform local reference derivatives of shape function to global derivatives
VolumeElementForcesAndSourcesCore * getVolumeFE() const
return pointer to Generic Volume Finite Element object
FTensor::Tensor2< double *, 3, 3 > & getInvJac()
get element inverse Jacobian
FTensor::Tensor2< double *, 3, 3 > & getJac()
get element Jacobian
virtual MoFEMErrorCode transformBaseFunctions()
Transform base functions based on geometric element Jacobian.
virtual MoFEMErrorCode setIntegrationPts()
Set integration points.
MoFEMErrorCode operator()()
function is run for every finite element
virtual MoFEMErrorCode getSpaceBaseAndOrderOnElement()
Determine approximation space and order of base functions.
virtual MoFEMErrorCode calculateVolumeAndJacobian()
Calculate element volume and Jacobian.
virtual MoFEMErrorCode calculateCoordinatesAtGaussPts()
Calculate coordinate at integration points.
VolumeElementForcesAndSourcesCore(Interface &m_field, const EntityType type=MBTET)