12 meshPositionsFieldName(
"MESH_NODE_POSITIONS"), aRea(elementMeasure) {}
29 auto get_ftensor_n_diff = [&]() {
37 &
m(0, 0), &
m(0, 1), &
m(0, 2));
42 const size_t nb_gauss_pts =
gaussPts.size2();
58 for (
int gg = 0; gg != nb_gauss_pts; ++gg) {
67 }
else if (type == MBQUAD) {
75 const size_t nb_gauss_pts =
gaussPts.size2();
90 auto t_diff = get_ftensor_n_diff();
91 for (
int gg = 0; gg != nb_gauss_pts; ++gg) {
92 auto t_coords = get_ftensor_from_vec_3d(
coords);
94 t_t1(
i) += t_coords(
i) * t_diff(
N0);
95 t_t2(
i) += t_coords(
i) * t_diff(
N1);
107 "Element type not implemented");
124 auto calc_normal = [&](
const double *diff_ptr) {
135 diff_ptr, &diff_ptr[1]);
146 for (
int nn = 0; nn !=
num_nodes; ++nn) {
147 t_t1(
i) += t_coords(
i) * t_diff(
N0);
148 t_t2(
i) += t_coords(
i) * t_diff(
N1);
153 aRea = sqrt(t_normal(
i) * t_normal(
i));
157 const double *diff_ptr;
161 CHKERR calc_normal(diff_ptr);
167 CHKERR calc_normal(diff_ptr);
171 "Element type not implemented");
186 auto get_rule_by_type = [&]() {
189 return getRule(order_row + 1, order_col + 1, order_data + 1);
191 return getRule(order_row, order_col, order_data);
195 const int rule = get_rule_by_type();
197 auto set_integration_pts_for_tri = [&]() {
208 gaussPts.resize(3, nb_gauss_pts,
false);
209 cblas_dcopy(nb_gauss_pts, &
QUAD_2D_TABLE[rule]->points[1], 3,
211 cblas_dcopy(nb_gauss_pts, &
QUAD_2D_TABLE[rule]->points[2], 3,
219 cblas_dcopy(3 * nb_gauss_pts,
QUAD_2D_TABLE[rule]->points, 1, shape_ptr,
232 auto calc_base_for_tri = [&]() {
234 const size_t nb_gauss_pts =
gaussPts.size2();
237 base.resize(nb_gauss_pts, 3,
false);
238 diff_base.resize(3, 2,
false);
247 auto calc_base_for_quad = [&]() {
249 const size_t nb_gauss_pts =
gaussPts.size2();
252 base.resize(nb_gauss_pts, 4,
false);
253 diff_base.resize(nb_gauss_pts, 8,
false);
254 for (
int gg = 0; gg != nb_gauss_pts; ++gg) {
276 CHKERR set_integration_pts_for_tri();
281 CHKERR calc_base_for_quad();
285 "Element type not implemented: %d", type);
291 const size_t nb_gauss_pts =
gaussPts.size2();
295 CHKERR calc_base_for_tri();
298 CHKERR calc_base_for_quad();
302 "Element type not implemented: %d", type);
378 const size_t nb_nodes =
380 double *shape_functions =
382 const size_t nb_gauss_pts =
gaussPts.size2();
384 for (
int gg = 0; gg != nb_gauss_pts; ++gg)
385 for (
int dd = 0; dd != 3; ++dd)
387 nb_nodes, &shape_functions[nb_nodes * gg], 1, &
coords[dd], 3);
397 "User operator and finite element do not work together");
442FaceElementForcesAndSourcesCore::UserDataOperator::loopSideVolumes(
444 return loopSide(fe_name, &fe_method, 3);
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
@ L2
field with C-1 continuity
@ HCURL
field with continuous tangents
@ HDIV
field with continuous normal traction
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
@ MOFEM_DATA_INCONSISTENCY
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define CHKERR
Inline error check.
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
FTensor::Index< 'm', SPACE_DIM > m
FTensor::Index< 'i', SPACE_DIM > i
const double v
phase velocity of light in medium (cm/ns)
FTensor::Index< 'j', 3 > j
FTensor::Index< 'k', 3 > k
constexpr std::enable_if<(Dim0<=2 &&Dim1<=2), Tensor2_Expr< Levi_Civita< T >, T, Dim0, Dim1, i, j > >::type levi_civita(const Index< i, Dim0 > &, const Index< j, Dim1 > &)
levi_civita functions to make for easy adhoc use
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
implementation of Data Operators for Forces and Sources
double zeta
Viscous hardening.
#define QUAD_2D_TABLE_SIZE
static QUAD *const QUAD_2D_TABLE[]
virtual moab::Interface & get_moab()=0
Deprecated interface functions.
std::array< std::bitset< LASTSPACE >, MBMAXTYPE > spacesOnEntities
spaces on entity types
std::array< boost::ptr_vector< EntData >, MBMAXTYPE > dataOnEntities
boost::shared_ptr< const NumeredEntFiniteElement > numeredEntFiniteElementPtr
MatrixDouble tangentOneAtGaussPts
virtual MoFEMErrorCode calculateAreaAndNormal()
Calculate element area and normal of the face.
MoFEMErrorCode operator()()
function is run for every finite element
FaceElementForcesAndSourcesCore(Interface &m_field)
MatrixDouble normalsAtGaussPts
virtual MoFEMErrorCode getSpaceBaseAndOrderOnElement()
Determine approximation space and order of base functions.
virtual MoFEMErrorCode calculateAreaAndNormalAtIntegrationPts()
Calculate element area and normal of the face at integration points.
virtual MoFEMErrorCode setIntegrationPts()
Set integration points.
const EntityHandle * conn
virtual MoFEMErrorCode calculateCoordinatesAtGaussPts()
Calculate coordinate at integration points.
MatrixDouble tangentTwoAtGaussPts
ForcesAndSourcesCore * ptrFE
structure to get information form mofem into EntitiesFieldData
MoFEMErrorCode loopOverOperators()
Iterate user data operators.
int getMaxRowOrder() const
Get max order of approximation for field in rows.
EntitiesFieldData & dataHdiv
MoFEMErrorCode getSpacesAndBaseOnEntities(EntitiesFieldData &data) const
Get field approximation space and base on entities.
virtual MoFEMErrorCode setGaussPts(int order_row, int order_col, int order_data)
set user specific integration rule
EntitiesFieldData & dataHcurl
auto & getElementPolynomialBase()
Get the Entity Polynomial Base object.
EntitiesFieldData & dataH1
MatrixDouble coordsAtGaussPts
coordinated at gauss points
virtual int getRule(int order_row, int order_col, int order_data)
another variant of getRule
MoFEMErrorCode calHierarchicalBaseFunctionsOnElement()
Calculate base functions.
MatrixDouble gaussPts
Matrix of integration points.
MoFEMErrorCode calBernsteinBezierBaseFunctionsOnElement()
Calculate Bernstein-Bezier base.
int getMaxColOrder() const
Get max order of approximation for field in columns.
EntityType lastEvaluatedElementEntityType
Last evaluated type of element entity.
int getMaxDataOrder() const
Get max order of approximation for data fields.
MoFEMErrorCode createDataOnElement(EntityType type)
Create a entity data on element object.
EntitiesFieldData & dataL2
Calculate base functions on triangle.
Calculate base functions on triangle.
Base volume element used to integrate on skeleton.