10#ifndef __ENTITIES_FIELD_DATA_HPP__
11#define __ENTITIES_FIELD_DATA_HPP__
13using namespace boost::numeric;
23using VectorDofs = ublas::vector<FEDofEntity *, DofsAllocator>;
45 std::bitset<LASTBASE>
bAse;
49 std::array<std::bitset<LASTSPACE>, MBMAXTYPE>
51 std::array<std::bitset<LASTBASE>, MBMAXTYPE>
53 std::array<std::bitset<LASTBASE>,
LASTSPACE>
55 std::array<boost::ptr_vector<EntData>, MBMAXTYPE>
114 const boost::shared_ptr<EntitiesFieldData> &data_ptr);
118 const boost::shared_ptr<EntitiesFieldData>
dataPtr;
142 EntData(
const bool allocate_base_matrices =
true);
216 template <
int Tensor_Dim>
232 template <
int Tensor_Dim0,
int Tensor_Dim1>
234 Tensor_Dim0, Tensor_Dim1>
249 template <
int Tensor_Dim>
285 virtual boost::shared_ptr<MatrixDouble> &
292 virtual boost::shared_ptr<MatrixDouble> &
298 virtual boost::shared_ptr<MatrixDouble> &
425 const int gg,
const int nb_base_functions);
454 const int nb_base_functions);
468 const int nb_base_functions);
502 template <
int DIM0,
int DIM1>
512 template <
int DIM0,
int DIM1>
522 template <
int DIM0,
int DIM1>
524 const int dof,
const int gg);
532 template <
int DIM0,
int DIM1>
611 template <
int Tensor_Dim>
631 template <
int Tensor_Dim>
653 template <
int Tensor_Dim>
676 template <
int Tensor_Dim>
705 template <
int Tensor_Dim>
733 template <
int Tensor_Dim0,
int Tensor_Dim1>
735 Tensor_Dim0, Tensor_Dim1>
741 template <
int Tensor_Dim0,
int Tensor_Dim1>
743 Tensor_Dim0, Tensor_Dim1>
748 template <
int Tensor_Dim0,
int Tensor_Dim1>
750 Tensor_Dim0, Tensor_Dim1>
752 return getFTensor2DiffN<Tensor_Dim0, Tensor_Dim1>(
bAse);
758 template <
int Tensor_Dim0,
int Tensor_Dim1>
760 Tensor_Dim0, Tensor_Dim1>
762 return getFTensor2DiffN<Tensor_Dim0, Tensor_Dim1>(
bAse, gg, bb);
767 template <
int Tensor_Dim0,
int Tensor_Dim1,
int Tensor_Dim2>
770 Tensor_Dim0, Tensor_Dim1, Tensor_Dim2>
775 template <
int Tensor_Dim0,
int Tensor_Dim1,
int Tensor_Dim2>
778 Tensor_Dim0, Tensor_Dim1, Tensor_Dim2>
780 return getFTensor3Diff2N<Tensor_Dim0, Tensor_Dim1, Tensor_Dim2>(
bAse);
797 template <
int Tensor_Dim>
815 template <
int Tensor_Dim>
844 template <
int Tensor_Dim0,
int Tensor_Dim1>
846 Tensor_Dim0, Tensor_Dim1>
871 template <
int Tensor_Dim0,
int Tensor_Dim1>
auto getFTensor2N();
899 template <
int Tensor_Dim0,
int Tensor_Dim1>
901 Tensor_Dim0, Tensor_Dim1>
925 template <
int Tensor_Dim0,
int Tensor_Dim1>
934 friend std::ostream &
operator<<(std::ostream &os,
963 virtual boost::shared_ptr<MatrixInt> &
969 virtual boost::shared_ptr<MatrixDouble> &
975 virtual const boost::shared_ptr<MatrixDouble> &
981 virtual boost::shared_ptr<MatrixDouble> &
987 virtual const boost::shared_ptr<MatrixDouble> &
990 virtual std::map<std::string, boost::shared_ptr<MatrixInt>> &
998 virtual std::map<std::string, boost::shared_ptr<MatrixDouble>> &
getBBNMap();
1006 virtual std::map<std::string, boost::shared_ptr<MatrixDouble>> &
1015 virtual boost::shared_ptr<MatrixInt> &
1024 virtual boost::shared_ptr<MatrixDouble> &
1033 virtual boost::shared_ptr<MatrixDouble> &
1077 std::array<std::array<boost::shared_ptr<MatrixDouble>,
LASTBASE>,
1082 std::array<boost::shared_ptr<MatrixDouble>,
LASTBASE>
1087 std::map<std::string, boost::shared_ptr<MatrixDouble>>
bbN;
1088 std::map<std::string, boost::shared_ptr<MatrixDouble>>
bbDiffN;
1089 std::map<std::string, boost::shared_ptr<MatrixInt>>
1131 const boost::shared_ptr<EntitiesFieldData::EntData> &ent_data_ptr);
1138 boost::shared_ptr<MatrixDouble> &
1142 boost::shared_ptr<MatrixDouble> &
1145 boost::shared_ptr<MatrixDouble> &
1148 const boost::shared_ptr<MatrixDouble> &
1151 const boost::shared_ptr<MatrixDouble> &
1154 inline boost::shared_ptr<MatrixDouble> &
1157 inline boost::shared_ptr<MatrixDouble> &
1160 boost::shared_ptr<MatrixInt> &
1166 boost::shared_ptr<MatrixDouble> &
1172 const boost::shared_ptr<MatrixDouble> &
1178 boost::shared_ptr<MatrixDouble> &
1184 const boost::shared_ptr<MatrixDouble> &
1207 unsigned int size = 0;
1208 if (
auto dof = dOfs[0]) {
1209 size = dof->getOrderNbDofs(
order) * dof->getNbOfCoeffs();
1210 size = size < iNdices.size() ? size : iNdices.size();
1212 int *data = &*iNdices.data().begin();
1213 return VectorIntAdaptor(size, ublas::shallow_array_adaptor<int>(size, data));
1217 return localIndices;
1222 unsigned int size = 0;
1223 if (
auto dof = dOfs[0]) {
1224 size = dof->getOrderNbDofs(
order) * dof->getNbOfCoeffs();
1225 size = size < localIndices.size() ? size : localIndices.size();
1227 int *data = &*localIndices.data().begin();
1228 return VectorIntAdaptor(size, ublas::shallow_array_adaptor<int>(size, data));
1238 return localIndices;
1247 unsigned int size = 0;
1248 if (
auto dof = dOfs[0]) {
1249 size = dof->getOrderNbDofs(
order) * dof->getNbOfCoeffs();
1250 size = size < fieldData.size() ? size : fieldData.size();
1252 double *data = &*fieldData.data().begin();
1265 return fieldEntities;
1270 return fieldEntities;
1273template <
int Tensor_Dim>
1276 std::stringstream s;
1277 s <<
"Not implemented for this dimension dim = " << Tensor_Dim;
1281template <
int Tensor_Dim0,
int Tensor_Dim1>
1283 Tensor_Dim0, Tensor_Dim1>
1285 std::stringstream s;
1286 s <<
"Not implemented for this dimension dim0 = " << Tensor_Dim0;
1287 s <<
" and dim1 " << Tensor_Dim1;
1291template <
int Tensor_Dim>
1293 FTensor::PackPtr<
double *, (Tensor_Dim * (Tensor_Dim + 1)) / 2>, Tensor_Dim>
1295 std::stringstream s;
1296 s <<
"Not implemented for this dimension dim = " << Tensor_Dim;
1306 return *(getNSharedPtr(base));
1317 return *(getDiffNSharedPtr(base));
1324 if (!getNSharedPtr(base, derivative)) {
1326 "Ptr to base %s functions derivative %d is null",
1331 return *(getNSharedPtr(base, derivative));
1343 return getN(
bAse, derivative);
1349 int size = getN(base).size2();
1350 double *data = &getN(base)(gg, 0);
1351 return VectorAdaptor(size, ublas::shallow_array_adaptor<double>(size, data));
1355 return getN(
bAse, gg);
1364 if (getN(base).size1() == getDiffN(base).size1()) {
1365 int size = getN(base).size2();
1366 int dim = getDiffN(base).size2() / size;
1367 double *data = &getDiffN(base)(gg, 0);
1369 getN(base).size2(),
dim,
1370 ublas::shallow_array_adaptor<double>(getDiffN(base).size2(), data));
1375 getN(base).size1(), getN(base).size2(),
1376 ublas::shallow_array_adaptor<double>(getDiffN(base).data().size(),
1377 &getDiffN(base).data()[0]));
1382 return getDiffN(
bAse, gg);
1387 const int gg,
const int nb_base_functions) {
1388 (void)getN()(gg, nb_base_functions -
1390 double *data = &getN(base)(gg, 0);
1391 return VectorAdaptor(nb_base_functions, ublas::shallow_array_adaptor<double>(
1392 nb_base_functions, data));
1397 return getN(
bAse, gg, nb_base_functions);
1403 const int nb_base_functions) {
1407 if (getN(base).size1() == getDiffN(base).size1()) {
1408 (void)getN(base)(gg,
1411 int dim = getDiffN(base).size2() / getN(base).size2();
1412 double *data = &getDiffN(base)(gg, 0);
1414 nb_base_functions,
dim,
1415 ublas::shallow_array_adaptor<double>(
dim * nb_base_functions, data));
1420 getN(base).size1(), getN(base).size2(),
1421 ublas::shallow_array_adaptor<double>(getDiffN(base).data().size(),
1422 &getDiffN(base).data()[0]));
1428 const int nb_base_functions) {
1429 return getDiffN(
bAse, gg, nb_base_functions);
1436 if (PetscUnlikely(getN(base).size2() % DIM)) {
1440 const int nb_base_functions = getN(base).size2() / DIM;
1441 double *data = &getN(base)(gg, 0);
1443 nb_base_functions, DIM,
1444 ublas::shallow_array_adaptor<double>(DIM * nb_base_functions, data));
1449 return getVectorN<DIM>(
bAse, gg);
1452template <
int DIM0,
int DIM1>
1456 if (PetscUnlikely(getDiffN(base).size2() % (DIM0 * DIM1))) {
1460 const int nb_base_functions = getN(base).size2() / (DIM0 * DIM1);
1461 double *data = &getN(base)(gg, 0);
1463 ublas::shallow_array_adaptor<double>(
1464 DIM0 * DIM1 * nb_base_functions, data));
1467template <
int DIM0,
int DIM1>
1469 return getVectorDiffN<DIM0, DIM1>(
bAse, gg);
1472template <
int DIM0,
int DIM1>
1475 const int dof,
const int gg) {
1479 ublas::shallow_array_adaptor<double>(DIM0 * DIM1, data));
1482template <
int DIM0,
int DIM1>
1485 return getVectorDiffN<DIM0, DIM1>(
bAse, dof, gg);
1490 double *ptr = &*getN(base).data().begin();
1496 return getFTensor0N(
bAse);
1501 const int gg,
const int bb) {
1502 double *ptr = &getN(base)(gg, bb);
1508 return getFTensor0N(
bAse, gg, bb);
1512 return getFTensor1N<Tensor_Dim>(
bAse);
1515template <
int Tensor_Dim>
1517 return getFTensor1N<Tensor_Dim>(
bAse, gg, bb);
1520template <
int Tensor_Dim0,
int Tensor_Dim1>
1522 return getFTensor2N<Tensor_Dim0, Tensor_Dim1>(
bAse);
1525template <
int Tensor_Dim0,
int Tensor_Dim1>
1527 return getFTensor2N<Tensor_Dim0, Tensor_Dim1>(
bAse, gg, bb);
1537 return *getBBAlphaIndicesSharedPtr(bbFieldName);
1546boost::shared_ptr<MatrixDouble> &
1552boost::shared_ptr<MatrixDouble> &
1575template <
typename T = EntityStorage>
1578 const double *ptr, InsertMode iora) {
1579 static_assert(!std::is_same<T, T>::value,
1580 "VecSetValues value for this data storage is not implemented");
1587 const double *ptr, InsertMode iora) {
1607template <
typename T = EntityStorage>
1611 return VecSetValues<T>(V, data, &*vec.data().begin(), iora);
1630template <
typename T = EntityStorage>
1634 const double *ptr, InsertMode iora) {
1635 static_assert(!std::is_same<T, T>::value,
1636 "MatSetValues value for this data storage is not implemented");
1656template <
typename T = EntityStorage>
1661 return MatSetValues<T>(
M, row_data, col_data, &*mat.data().begin(), iora);
1668 const double *ptr, InsertMode iora) {
1685 const int gg,
const int bb);
1699EntitiesFieldData::EntData::getFTensor1DiffN<3>(
1703EntitiesFieldData::EntData::getFTensor1DiffN<3>();
1707EntitiesFieldData::EntData::getFTensor1DiffN<2>(
1711EntitiesFieldData::EntData::getFTensor1DiffN<2>();
1719 const int gg,
const int bb);
1733EntitiesFieldData::EntData::getFTensor1FieldData<3>();
1737EntitiesFieldData::EntData::getFTensor1FieldData<2>();
1741EntitiesFieldData::EntData::getFTensor1FieldData<1>();
1745EntitiesFieldData::EntData::getFTensor2FieldData<3, 3>();
1749EntitiesFieldData::EntData::getFTensor2SymmetricFieldData<3>();
1753EntitiesFieldData::EntData::getFTensor2SymmetricFieldData<2>();
#define MOFEM_LOG_C(channel, severity, format,...)
FieldApproximationBase
approximation base
FieldSpace
approximation spaces
@ LASTSPACE
FieldSpace in [ 0, LASTSPACE )
static const char *const ApproximationBaseNames[]
#define BITFEID_SIZE
max number of finite elements
#define THROW_MESSAGE(msg)
Throw MoFEM exception.
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
VectorShallowArrayAdaptor< int > VectorIntAdaptor
int ApproximationOrder
Approximation on the entity.
VectorShallowArrayAdaptor< double > VectorAdaptor
UBlasMatrix< double > MatrixDouble
UBlasVector< double > VectorDouble
MatrixShallowArrayAdaptor< double > MatrixAdaptor
Matrix adaptor.
UBlasMatrix< int > MatrixInt
UBlasVector< int > VectorInt
implementation of Data Operators for Forces and Sources
auto getVectorAdaptor(T1 ptr, const size_t n)
Get Vector adaptor.
DEPRECATED typedef DerivedEntitiesFieldData DerivedDataForcesAndSourcesCore
ublas::vector< FieldEntity *, FieldEntAllocator > VectorFieldEntities
DEPRECATED typedef EntitiesFieldData DataForcesAndSourcesCore
MoFEMErrorCode MatSetValues(Mat M, const EntitiesFieldData::EntData &row_data, const EntitiesFieldData::EntData &col_data, const double *ptr, InsertMode iora)
Assemble PETSc matrix.
ublas::unbounded_array< FieldEntity *, std::allocator< FieldEntity * > > FieldEntAllocator
MoFEMErrorCode VecSetValues< EntityStorage >(Vec V, const EntitiesFieldData::EntData &data, const double *ptr, InsertMode iora)
MoFEMErrorCode MatSetValues< EntityStorage >(Mat M, const EntitiesFieldData::EntData &row_data, const EntitiesFieldData::EntData &col_data, const double *ptr, InsertMode iora)
ublas::unbounded_array< FEDofEntity *, std::allocator< FEDofEntity * > > DofsAllocator
MoFEMErrorCode VecSetValues(Vec V, const EntitiesFieldData::EntData &data, const double *ptr, InsertMode iora)
Assemble PETSc vector.
ublas::vector< FEDofEntity *, DofsAllocator > VectorDofs
constexpr auto field_name
Derived ata on single entity (This is passed as argument to DataOperator::doWork)
boost::shared_ptr< MatrixDouble > & getBBNSharedPtr(const std::string &field_name)
boost::shared_ptr< MatrixDouble > & getDerivedDiffNSharedPtr(const FieldApproximationBase base)
const boost::shared_ptr< EntitiesFieldData::EntData > entDataPtr
MoFEMErrorCode baseSwap(const std::string &field_name, const FieldApproximationBase base)
Used by Bernstein base to keep temporally pointer.
std::vector< BitRefLevel > & getEntDataBitRefLevel()
int getSense() const
get entity sense, need to calculate base functions with conforming approximation fields
boost::shared_ptr< MatrixInt > & getBBAlphaIndicesSharedPtr(const std::string &field_name)
boost::shared_ptr< MatrixDouble > & getBBDiffNSharedPtr(const std::string &field_name)
boost::shared_ptr< MatrixDouble > & getDiffNSharedPtr(const FieldApproximationBase base)
boost::shared_ptr< MatrixDouble > & getNSharedPtr(const FieldApproximationBase base, const BaseDerivatives derivative)
boost::shared_ptr< MatrixDouble > & getDerivedNSharedPtr(const FieldApproximationBase base)
this class derive data form other data structure
MoFEMErrorCode setElementType(const EntityType type)
const boost::shared_ptr< EntitiesFieldData > dataPtr
Data on single entity (This is passed as argument to DataOperator::doWork)
std::array< boost::shared_ptr< MatrixDouble >, MaxBernsteinBezierOrder > bbDiffNByOrder
BB base functions direvatives by order.
boost::shared_ptr< MatrixDouble > swapBaseDiffNPtr
Used by Bernstein base to keep temporally pointer.
virtual std::array< boost::shared_ptr< MatrixDouble >, MaxBernsteinBezierOrder > & getBBDiffNByOrderArray()
FieldSpace sPace
Entity space.
const VectorIntAdaptor getIndicesUpToOrder(int order)
get global indices of dofs on entity up to given order
virtual std::array< boost::shared_ptr< MatrixDouble >, MaxBernsteinBezierOrder > & getBBNByOrderArray()
int sEnse
Entity sense (orientation)
auto getFTensor2N()
Get base functions for Hdiv space.
VectorDouble fieldData
Field data on entity.
FTensor::Tensor1< FTensor::PackPtr< double *, Tensor_Dim >, Tensor_Dim > getFTensor1DiffN(const int gg, const int bb)
Get derivatives of base functions (Loop by integration points)
MatrixDouble & getN()
get base functions this return matrix (nb. of rows is equal to nb. of Gauss pts, nb....
FieldApproximationBase & getBase()
Get approximation base.
boost::shared_ptr< MatrixDouble > swapBaseNPtr
Used by Bernstein base to keep temporally pointer.
ApproximationOrder getOrder() const
get approximation order
std::map< std::string, boost::shared_ptr< MatrixInt > > bbAlphaIndices
Indices for Bernstein-Bezier (BB) base.
friend std::ostream & operator<<(std::ostream &os, const EntitiesFieldData::EntData &e)
VectorInt localIndices
Local indices on entity.
const MatrixAdaptor getVectorDiffN(FieldApproximationBase base, const int gg)
get DiffHdiv of base functions at Gauss pts
const VectorFieldEntities & getFieldEntities() const
get field entities
virtual std::array< boost::shared_ptr< MatrixInt >, MaxBernsteinBezierOrder > & getBBAlphaIndicesByOrderArray()
virtual std::map< std::string, boost::shared_ptr< MatrixInt > > & getBBAlphaIndicesMap()
const VectorIntAdaptor getLocalIndicesUpToOrder(int order)
get local indices of dofs on entity up to given order
const VectorAdaptor getFieldDataUpToOrder(int order)
get dofs values up to given order
FTensor::Tensor1< FTensor::PackPtr< double *, Tensor_Dim >, Tensor_Dim > getFTensor1DiffN()
Get derivatives of base functions.
std::vector< BitRefLevel > entDataBitRefLevel
Bit ref level in entity.
FTensor::Tensor2< FTensor::PackPtr< double *, Tensor_Dim0 *Tensor_Dim1 >, Tensor_Dim0, Tensor_Dim1 > getFTensor2DiffN(const int gg, const int bb)
Get derivatives of base functions for Hdiv space at integration pts.
static constexpr size_t MaxBernsteinBezierOrder
VectorInt iNdices
Global indices on entity.
virtual boost::shared_ptr< MatrixInt > & getBBAlphaIndicesSharedPtr(const std::string &field_name)
const VectorDouble & getFieldData() const
get dofs values
VectorInt bbNodeOrder
order of nodes
virtual boost::shared_ptr< MatrixDouble > & getBBDiffNByOrderSharedPtr(const size_t o)
get BB base derivative by order
FTensor::Tensor3< FTensor::PackPtr< double *, Tensor_Dim0 *Tensor_Dim1 *Tensor_Dim2 >, Tensor_Dim0, Tensor_Dim1, Tensor_Dim2 > getFTensor3Diff2N(FieldApproximationBase base)
Get second derivatives of base functions for Hvec space.
virtual std::map< std::string, boost::shared_ptr< MatrixDouble > > & getBBNMap()
get hash map of base function for BB base, key is a field name
const MatrixAdaptor getVectorN(const FieldApproximationBase base, const int gg)
get Hdiv of base functions at Gauss pts
virtual MoFEMErrorCode baseSwap(const std::string &field_name, const FieldApproximationBase base)
Swap bases functions.
const VectorInt & getLocalIndices() const
get local indices of dofs on entity
ApproximationOrder oRder
Entity order.
FTensor::Tensor1< FTensor::PackPtr< double *, Tensor_Dim >, Tensor_Dim > getFTensor1FieldData()
Return FTensor of rank 1, i.e. vector from filed data coeffinects.
virtual std::map< std::string, boost::shared_ptr< MatrixDouble > > & getBBDiffNMap()
get hash map of direvarives base function for BB base, key is a field name
virtual std::vector< BitRefLevel > & getEntDataBitRefLevel()
virtual int getSense() const
get entity sense, need to calculate base functions with conforming approximation fields
VectorFieldEntities fieldEntities
Field entities.
virtual ~EntData()=default
FTensor::Tensor0< FTensor::PackPtr< double *, 1 > > getFTensor0N()
Get base function as Tensor0.
FTensor::Tensor0< FTensor::PackPtr< double *, 1 > > getFTensor0FieldData()
Resturn scalar files as a FTensor of rank 0.
MatrixInt & getBBAlphaIndices()
Get file BB indices.
std::array< boost::shared_ptr< MatrixInt >, MaxBernsteinBezierOrder > bbAlphaIndicesByOrder
BB alpha indices by order.
MatrixDouble & getDiffN()
get derivatives of base functions
VectorInt & getBBNodeOrder()
Get orders at the nodes.
std::string bbFieldName
field name
virtual boost::shared_ptr< MatrixInt > & getBBAlphaIndicesByOrderSharedPtr(const size_t o)
get ALpha indices for BB base by order
std::array< boost::shared_ptr< MatrixDouble >, LASTBASE > & diffN
Derivatives of base functions.
std::map< std::string, boost::shared_ptr< MatrixDouble > > bbDiffN
MoFEMErrorCode resetFieldDependentData()
const VectorDofs & getFieldDofs() const
get dofs data stature FEDofEntity
FTensor::Tensor2_symmetric< FTensor::PackPtr< double *,(Tensor_Dim *(Tensor_Dim+1))/2 >, Tensor_Dim > getFTensor2SymmetricFieldData()
Return symmetric FTensor rank 2, i.e. matrix from filed data coeffinects.
FTensor::Tensor3< FTensor::PackPtr< double *, Tensor_Dim0 *Tensor_Dim1 *Tensor_Dim2 >, Tensor_Dim0, Tensor_Dim1, Tensor_Dim2 > getFTensor3Diff2N()
Get second derivatives of base functions for Hvec space.
FieldSpace & getSpace()
Get field space.
auto getFTensor1N()
Get base functions for Hdiv space.
virtual boost::shared_ptr< MatrixDouble > & getBBNSharedPtr(const std::string &field_name)
std::array< boost::shared_ptr< MatrixDouble >, MaxBernsteinBezierOrder > bbNByOrder
BB base functions by order.
std::array< std::array< boost::shared_ptr< MatrixDouble >, LASTBASE >, LastDerivative > baseFunctionsAndBaseDerivatives
FTensor::Tensor2< FTensor::PackPtr< double *, Tensor_Dim0 *Tensor_Dim1 >, Tensor_Dim0, Tensor_Dim1 > getFTensor2DiffN()
Get derivatives of base functions for Hdiv space.
std::map< std::string, boost::shared_ptr< MatrixDouble > > bbN
VectorDofs dOfs
DoFs on entity.
virtual boost::shared_ptr< MatrixDouble > & getDiffNSharedPtr(const FieldApproximationBase base)
FTensor::Tensor2< FTensor::PackPtr< double *, Tensor_Dim0 *Tensor_Dim1 >, Tensor_Dim0, Tensor_Dim1 > getFTensor2FieldData()
Return FTensor rank 2, i.e. matrix from filed data coeffinects.
virtual boost::shared_ptr< MatrixDouble > & getBBNByOrderSharedPtr(const size_t o)
get BB base by order
virtual boost::shared_ptr< MatrixDouble > & getNSharedPtr(const FieldApproximationBase base, const BaseDerivatives direvatie)
const VectorInt & getIndices() const
Get global indices of dofs on entity.
FieldApproximationBase bAse
Field approximation base.
virtual boost::shared_ptr< MatrixDouble > & getBBDiffNSharedPtr(const std::string &field_name)
std::array< boost::shared_ptr< MatrixDouble >, LASTBASE > & N
Base functions.
data structure for finite element entity
friend std::ostream & operator<<(std::ostream &os, const EntitiesFieldData &e)
MoFEMErrorCode resetFieldDependentData()
std::bitset< LASTSPACE > sPace
spaces on element
MatrixInt facesNodes
nodes on finite element faces
std::array< std::bitset< LASTSPACE >, MBMAXTYPE > spacesOnEntities
spaces on entity types
std::bitset< LASTBASE > bAse
bases on element
std::array< boost::ptr_vector< EntData >, MBMAXTYPE > dataOnEntities
std::array< std::bitset< LASTBASE >, LASTSPACE > basesOnSpaces
base on spaces
virtual ~EntitiesFieldData()=default
virtual MoFEMErrorCode baseSwap(const std::string &field_name, const FieldApproximationBase base)
Swap approximation base.
virtual MoFEMErrorCode setElementType(const EntityType type)
MatrixInt facesNodesOrder
order of face nodes on element
std::array< std::bitset< LASTBASE >, MBMAXTYPE > basesOnEntities
bases on entity types
keeps information about indexed dofs for the finite element
Struct keeps handle to entity in the field.
Operator to project base functions from parent entity to child.