v0.10.0
VolumeElementForcesAndSourcesCoreOnSide.cpp
Go to the documentation of this file.
1 /** \file VolumeElementForcesAndSourcesCoreOnSide.cpp
2 
3 \brief Implementation of volume element on side
4 
5 */
6 
7 /* This file is part of MoFEM.
8  * MoFEM is free software: you can redistribute it and/or modify it under
9  * the terms of the GNU Lesser General Public License as published by the
10  * Free Software Foundation, either version 3 of the License, or (at your
11  * option) any later version.
12  *
13  * MoFEM is distributed in the hope that it will be useful, but WITHOUT
14  * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
15  * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
16  * License for more details.
17  *
18  * You should have received a copy of the GNU Lesser General Public
19  * License along with MoFEM. If not, see <http://www.gnu.org/licenses/>. */
20 
21 namespace MoFEM {
22 
24  return -1;
25 };
26 
30 
31  if (!sidePtrFE)
32  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "Side element not set");
33 
34  const EntityHandle face_entity =
36  SideNumber_multiIndex &side_table = const_cast<SideNumber_multiIndex &>(
37  numeredEntFiniteElementPtr->getSideNumberTable());
38  SideNumber_multiIndex::nth_index<0>::type::iterator sit =
39  side_table.get<0>().find(face_entity);
40  if (sit == side_table.get<0>().end())
41  SETERRQ(PETSC_COMM_WORLD, MOFEM_DATA_INCONSISTENCY,
42  "Face can not be found on volume element");
43 
44  auto face_ptr_fe =
45  static_cast<FaceElementForcesAndSourcesCoreBase *>(sidePtrFE);
46 
47  faceSense = (*sit)->sense;
48  faceSideNumber = (*sit)->side_number;
49  fill(tetConnMap.begin(), tetConnMap.end(), -1);
50  for (int nn = 0; nn != 3; ++nn) {
51  faceConnMap[nn] =
52  std::distance(conn, find(conn, &conn[4], face_ptr_fe->conn[nn]));
53  tetConnMap[faceConnMap[nn]] = nn;
54  if (faceConnMap[nn] > 3)
55  SETERRQ(PETSC_COMM_WORLD, MOFEM_DATA_INCONSISTENCY,
56  "No common node on face and element can not be found");
57  }
58 
59  oppositeNode = std::distance(tetConnMap.begin(),
60  find(tetConnMap.begin(), tetConnMap.end(), -1));
61 
62  const int nb_gauss_pts = face_ptr_fe->gaussPts.size2();
63  gaussPts.resize(4, nb_gauss_pts, false);
64  gaussPts.clear();
65  DataForcesAndSourcesCore &dataH1_on_face = *face_ptr_fe->dataOnElement[H1];
66  const MatrixDouble &face_shape_funtions =
67  dataH1_on_face.dataOnEntities[MBVERTEX][0].getN(NOBASE);
68  const double tet_coords[] = {0, 0, 0, 1, 0, 0, 0, 1, 0, 0, 0, 1};
69  for (int gg = 0; gg != nb_gauss_pts; ++gg) {
70  gaussPts(0, gg) =
71  face_shape_funtions(gg, 0) * tet_coords[3 * faceConnMap[0] + 0] +
72  face_shape_funtions(gg, 1) * tet_coords[3 * faceConnMap[1] + 0] +
73  face_shape_funtions(gg, 2) * tet_coords[3 * faceConnMap[2] + 0];
74  gaussPts(1, gg) =
75  face_shape_funtions(gg, 0) * tet_coords[3 * faceConnMap[0] + 1] +
76  face_shape_funtions(gg, 1) * tet_coords[3 * faceConnMap[1] + 1] +
77  face_shape_funtions(gg, 2) * tet_coords[3 * faceConnMap[2] + 1];
78  gaussPts(2, gg) =
79  face_shape_funtions(gg, 0) * tet_coords[3 * faceConnMap[0] + 2] +
80  face_shape_funtions(gg, 1) * tet_coords[3 * faceConnMap[1] + 2] +
81  face_shape_funtions(gg, 2) * tet_coords[3 * faceConnMap[2] + 2];
82  gaussPts(3, gg) = face_ptr_fe->gaussPts(2, gg);
83  }
85 }
86 
89  ForcesAndSourcesCore *ptr) {
91  if(!(ptrFE = dynamic_cast<VolumeElementForcesAndSourcesCoreOnSideBase *>(ptr)))
92  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
93  "User operator and finite element do not work together");
95 }
96 
97 } // namespace MoFEM
MatrixDouble gaussPts
Matrix of integration points.
ForcesAndSourcesCore * sidePtrFE
Element to integrate on the sides.
multi_index_container< boost::shared_ptr< SideNumber >, indexed_by< ordered_unique< member< SideNumber, EntityHandle, &SideNumber::ent > >, ordered_non_unique< composite_key< SideNumber, const_mem_fun< SideNumber, EntityType, &SideNumber::getEntType >, member< SideNumber, char, &SideNumber::side_number > > > > > SideNumber_multiIndex
SideNumber_multiIndex for SideNumber.
data structure for finite element entityIt keeps that about indices of degrees of freedom,...
std::array< boost::ptr_vector< EntData >, MBMAXTYPE > dataOnEntities
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:509
ublas::matrix< double, ublas::row_major, DoubleAllocator > MatrixDouble
Definition: Types.hpp:76
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:485
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:516
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
boost::shared_ptr< const NumeredEntFiniteElement > numeredEntFiniteElementPtr
constexpr int order
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:67
structure to get information form mofem into DataForcesAndSourcesCore
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:415
continuous field
Definition: definitions.h:177