23 #ifndef __VOLUMEELEMENTFORCESANDSOURCESCORE_ONSIDE_HPP__
24 #define __VOLUMEELEMENTFORCESANDSOURCESCORE_ONSIDE_HPP__
26 using namespace boost::numeric;
46 inline const std::array<int, 4> &getFaceConnMap()
const;
56 inline const std::array<int, 8> &getTetConnMap()
const;
66 inline int getOppositeNode()
const;
73 inline int getFaceSense()
const;
80 inline int getFaceSideNumber()
const;
87 int getRule(
int order);
114 inline int getFaceSense()
const;
120 inline int getFaceSideNumber()
const;
122 inline bool getEdgeFace(
const int ee)
const;
140 inline auto getFTensor1Normal();
144 inline auto getFTensor1Tangent1();
148 inline auto getFTensor1Tangent2();
162 inline ublas::matrix_row<MatrixDouble> getNormalsAtGaussPts(
const int gg);
178 inline auto getFTensor1NormalsAtGaussPts();
200 template <
int SWITCH>
204 using VolumeElementForcesAndSourcesCoreOnSideBase::
205 VolumeElementForcesAndSourcesCoreOnSideBase;
220 const std::array<int, 4> &
221 VolumeElementForcesAndSourcesCoreOnSideBase::getFaceConnMap()
const {
225 const std::array<int, 8> &
226 VolumeElementForcesAndSourcesCoreOnSideBase::getTetConnMap()
const {
230 int VolumeElementForcesAndSourcesCoreOnSideBase::getOppositeNode()
const {
234 int VolumeElementForcesAndSourcesCoreOnSideBase::getFaceSense()
const {
238 int VolumeElementForcesAndSourcesCoreOnSideBase::getFaceSideNumber()
const {
239 return faceSideNumber;
243 VolumeElementForcesAndSourcesCoreOnSideBase::UserDataOperator::getFaceFE()
246 getVolumeFE()->sidePtrFE);
250 VolumeElementForcesAndSourcesCoreOnSideBase::UserDataOperator::getNormal() {
251 return getFaceFE()->
nOrmal;
255 VolumeElementForcesAndSourcesCoreOnSideBase::UserDataOperator::getTangent1() {
256 return getFaceFE()->tangentOne;
260 VolumeElementForcesAndSourcesCoreOnSideBase::UserDataOperator::getTangent2() {
261 return getFaceFE()->tangentTwo;
265 VolumeElementForcesAndSourcesCoreOnSideBase::UserDataOperator::getVolumeFE()
270 int VolumeElementForcesAndSourcesCoreOnSideBase::UserDataOperator::
271 getFaceSense()
const {
275 auto VolumeElementForcesAndSourcesCoreOnSideBase::UserDataOperator::
276 getFTensor1Normal() {
277 double *ptr = &*getNormal().data().begin();
281 auto VolumeElementForcesAndSourcesCoreOnSideBase::UserDataOperator::
282 getFTensor1Tangent1() {
283 double *ptr = &*getTangent1().data().begin();
287 auto VolumeElementForcesAndSourcesCoreOnSideBase::UserDataOperator::
288 getFTensor1Tangent2() {
289 double *ptr = &*getTangent2().data().begin();
293 inline auto VolumeElementForcesAndSourcesCoreOnSideBase::UserDataOperator::
294 getFTensor1NormalsAtGaussPts() {
295 double *ptr = &*getNormalsAtGaussPts().data().begin();
300 int VolumeElementForcesAndSourcesCoreOnSideBase::UserDataOperator::
301 getFaceSideNumber()
const {
302 return getVolumeFE()->faceSideNumber;
305 MatrixDouble &VolumeElementForcesAndSourcesCoreOnSideBase::UserDataOperator::
306 getNormalsAtGaussPts() {
307 return getFaceFE()->normalsAtGaussPts;
310 MatrixDouble &VolumeElementForcesAndSourcesCoreOnSideBase::UserDataOperator::
311 getFaceCoordsAtGaussPts() {
312 return getFaceFE()->coordsAtGaussPts;
315 bool VolumeElementForcesAndSourcesCoreOnSideBase::UserDataOperator::getEdgeFace(
316 const int ee)
const {
317 constexpr
bool edges_on_faces[6][4] = {{
true,
false,
false,
true},
318 {
false,
true,
false,
true},
319 {
false,
false,
true,
true},
320 {
true,
false,
true,
false},
321 {
true,
true,
false,
false},
322 {
false,
true,
true,
false}};
323 return edges_on_faces[ee][getFaceSideNumber()];
326 ublas::matrix_row<MatrixDouble> VolumeElementForcesAndSourcesCoreOnSideBase::
327 UserDataOperator::getNormalsAtGaussPts(
const int gg) {
328 return ublas::matrix_row<MatrixDouble>(getNormalsAtGaussPts(), gg);
331 template <
int SWITCH>
334 return opSwitch<SWITCH>();
ForcesAndSourcesCore::UserDataOperator UserDataOperator
VolumeElementForcesAndSourcesCoreSwitch< 0 > VolumeElementForcesAndSourcesCore
Volume finite element default.
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
UBlasMatrix< double > MatrixDouble
UBlasVector< double > VectorDouble
implementation of Data Operators for Forces and Sources
structure to get information form mofem into EntitiesFieldData
default operator for TET element
Base volume element used to integrate on skeleton.
int faceSideNumber
Face side number.
std::array< int, 4 > faceConnMap
int faceSense
Sense of face, could be 1 or -1.
std::array< int, 8 > tetConnMap
Volume side finite element with switches.
Volume finite element with switches.