|
| v0.14.0
|
Go to the documentation of this file.
28 meshPositionsFieldName(
"MESH_NODE_POSITIONS"), lEngth(elementMeasure) {
32 "Problem with creation data on element");
53 int rule =
getRule(order_row, order_col, order_data);
66 gaussPts.resize(2, nb_gauss_pts,
false);
75 cblas_dcopy(2 * nb_gauss_pts,
QUAD_1D_TABLE[rule]->points, 1, shape_ptr,
89 for (
int gg = 0; gg != nb_gauss_pts; gg++) {
104 const auto nb_gauss_pts =
gaussPts.size2();
121 t_dir(
i) = t_coords_node1(
i) - t_coords_node0(
i);
123 for (
int gg = 0; gg != nb_gauss_pts; ++gg) {
124 t_coords(
i) =
N_MBEDGE0(t_ksi) * t_coords_node0(
i) +
126 t_tangent(
i) = t_dir(
i);
140 "User operator and finite element do not work together");
147 return loopSide(fe_name, &fe_side, 2);
165 CHKERR getEntityDataOrder<MBEDGE>(data_l2,
L2);
166 data_l2.dataOnEntities[MBEDGE][0].getSense() =
168 data_l2.spacesOnEntities[MBEDGE].set(
L2);
174 CHKERR getEntityDataOrder<MBEDGE>(data_curl,
HCURL);
175 data_curl.dataOnEntities[MBEDGE][0].getSense() =
177 data_curl.spacesOnEntities[MBEDGE].set(
HCURL);
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
MoFEMErrorCode calBernsteinBezierBaseFunctionsOnElement()
Calculate Bernstein-Bezier base.
MoFEMErrorCode setIntegrationPts()
int getMaxRowOrder() const
Get max order of approximation for field in rows.
EdgeElementForcesAndSourcesCore(Interface &m_field)
#define CHK_THROW_MESSAGE(err, msg)
Check and throw MoFEM exception.
@ L2
field with C-1 continuity
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
MatrixDouble tangentAtGaussPts
virtual MoFEMErrorCode setGaussPts(int order_row, int order_col, int order_data)
set user specific integration rule
int getMaxDataOrder() const
Get max order of approximation for data fields.
const EntityHandle * cOnn
MoFEMErrorCode calculateEdgeDirection()
auto & getElementPolynomialBase()
Get the Entity Polynomial Base object.
int getMaxColOrder() const
Get max order of approximation for field in columns.
MoFEMErrorCode calHierarchicalBaseFunctionsOnElement()
Calculate base functions.
std::array< std::bitset< LASTSPACE >, MBMAXTYPE > spacesOnEntities
spaces on entity types
Deprecated interface functions.
static QUAD *const QUAD_1D_TABLE[]
#define CHKERR
Inline error check.
virtual moab::Interface & get_moab()=0
implementation of Data Operators for Forces and Sources
MatrixDouble coordsAtGaussPts
coordinated at gauss points
const std::array< boost::shared_ptr< EntitiesFieldData >, LASTSPACE > dataOnElement
Entity data on element entity rows fields.
MoFEMErrorCode createDataOnElement(EntityType type)
Create a entity data on element object.
virtual int getRule(int order_row, int order_col, int order_data)
another variant of getRule
MoFEMErrorCode loopSideFaces(const string fe_name, FaceElementForcesAndSourcesCoreOnSide &fe_side)
FTensor::Index< 'i', SPACE_DIM > i
MoFEMErrorCode getSpacesAndBaseOnEntities(EntitiesFieldData &data) const
Get field approximation space and base on entities.
MoFEMErrorCode operator()()
function is run for every finite element
static FTensor::Tensor1< double, 3 > tFaceOrientation
structure to get information form mofem into EntitiesFieldData
EntitiesFieldData & dataH1
Calculate base functions on tetrahedral.
#define QUAD_1D_TABLE_SIZE
@ HCURL
field with continuous tangents
@ MOFEM_DATA_INCONSISTENCY
MatrixDouble gaussPts
Matrix of integration points.
std::array< boost::ptr_vector< EntData >, MBMAXTYPE > dataOnEntities
MoFEMErrorCode loopOverOperators()
Iterate user data operators.
MoFEMErrorCode calculateCoordsAtIntegrationPts()
boost::shared_ptr< const NumeredEntFiniteElement > numeredEntFiniteElementPtr
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
ForcesAndSourcesCore * ptrFE
Base face element used to integrate on skeleton.
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
MoFEMErrorCode setPtrFE(ForcesAndSourcesCore *ptr)