|
| v0.14.0
|
Go to the documentation of this file.
13 using namespace MoFEM;
15 static char help[] =
"...\n\n";
40 if (
type == MBEDGE && side == getEdgeSideNumber()) {
42 MatrixDouble diff = getCoordsAtGaussPts() - getEdgeCoordsAtGaussPts();
44 const double eps = 1e-12;
45 if (norm_inf(diff) >
eps)
47 "coordinates at integration pts are different");
49 const size_t nb_dofs = data.
getN().size2() / 3;
50 const size_t nb_integration_pts = data.
getN().size1();
52 auto &dir = getDirection();
57 if (getFEMethod()->nInTheLoop == 0)
62 dot_elem_data.resize(nb_integration_pts, nb_dofs,
false);
64 for (
size_t gg = 0; gg != nb_integration_pts; ++gg) {
65 for (
size_t bb = 0; bb != nb_dofs; ++bb) {
66 dot_elem_data(gg, bb) = t_normal(
i) * t_hdiv_base(
i);
80 faceSideFe(m_field), elemData(elem_data) {
82 auto jac_ptr = boost::make_shared<MatrixDouble>();
99 const size_t nb_dofs = data.
getN().size2() / 3;
100 const size_t nb_integration_pts = data.
getN().size1();
102 auto &dir = getDirection();
106 elemData.
dotEdge.resize(nb_integration_pts, nb_dofs,
false);
107 elemData.
dotEleLeft.resize(nb_integration_pts, 0,
false);
108 elemData.
dotEleRight.resize(nb_integration_pts, 0,
false);
110 for (
size_t gg = 0; gg != nb_integration_pts; ++gg) {
111 for (
size_t bb = 0; bb != nb_dofs; ++bb) {
112 elemData.
dotEdge(gg, bb) = t_normal(
i) * t_hdiv_base(
i);
117 CHKERR loopSideFaces(
"dFE", faceSideFe);
119 auto check_continuity_of_base = [&](
auto &vol_dot_data) {
122 if (vol_dot_data.size1() != elemData.
dotEdge.size1())
124 "Inconsistent number of integration points");
126 if (vol_dot_data.size2() != elemData.
dotEdge.size2())
128 "Inconsistent number of base functions");
129 const double eps = 1e-12;
130 for (
size_t gg = 0; gg != vol_dot_data.size1(); ++gg)
131 for (
size_t bb = 0; bb != vol_dot_data.size2(); ++bb) {
133 std::abs(vol_dot_data(gg, bb) - elemData.
dotEdge(gg, bb));
136 "Inconsistency (%d, %d) %3.4e != %3.4e", gg, bb,
137 vol_dot_data(gg, bb), elemData.
dotEdge(gg, bb));
152 int main(
int argc,
char *argv[]) {
168 DMType dm_name =
"DMMOFEM";
181 enum bases { AINSWORTH, DEMKOWICZ, LASTBASEOP };
182 const char *list_bases[] = {
"ainsworth",
"demkowicz"};
184 PetscInt choice_base_value = AINSWORTH;
186 LASTBASEOP, &choice_base_value, &flg);
187 if (flg == PETSC_TRUE) {
189 if (choice_base_value == AINSWORTH)
191 else if (choice_base_value == DEMKOWICZ)
199 auto base = get_base();
207 auto dm = simple_interface->
getDM();
211 boost::shared_ptr<EdgeEle> skeleton_fe =
212 boost::shared_ptr<EdgeEle>(
new EdgeEle(m_field));
214 skeleton_fe->getOpPtrVector().push_back(
216 skeleton_fe->getOpPtrVector().push_back(
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface reference to pointer of interface.
Data on single entity (This is passed as argument to DataOperator::doWork)
OpSetContravariantPiolaTransformOnFace2DImpl< 2 > OpSetContravariantPiolaTransformOnFace2D
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
UBlasMatrix< double > MatrixDouble
const std::string getSkeletonFEName() const
Get the Skeleton FE Name.
MoFEMErrorCode loadFile(const std::string options, const std::string mesh_file_name, LoadFileFunc loadFunc=defaultLoadFileFunc)
Load mesh file.
MoFEM::EdgeElementForcesAndSourcesCore EdgeEle
static MoFEMErrorCode Finalize()
Checks for options to be called at the conclusion of the program.
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
Simple interface for fast problem set-up.
Deprecated interface functions.
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
MoFEMErrorCode addSkeletonField(const std::string &name, const FieldSpace space, const FieldApproximationBase base, const FieldCoefficientsNumber nb_of_coefficients, const TagType tag_type=MB_TAG_SPARSE, const enum MoFEMTypes bh=MF_ZERO, int verb=-1)
Add field on skeleton.
DeprecatedCoreInterface Interface
MoFEMErrorCode getOptions()
get options
MoFEMErrorCode getDM(DM *dm)
Get DM.
#define CHKERR
Inline error check.
friend class UserDataOperator
implementation of Data Operators for Forces and Sources
int main(int argc, char *argv[])
MoFEMErrorCode addDomainField(const std::string &name, const FieldSpace space, const FieldApproximationBase base, const FieldCoefficientsNumber nb_of_coefficients, const TagType tag_type=MB_TAG_SPARSE, const enum MoFEMTypes bh=MF_ZERO, int verb=-1)
Add field on domain.
PetscErrorCode DMRegister_MoFEM(const char sname[])
Register MoFEM problem.
OpFaceSide(CommonData &elem_data)
friend class UserDataOperator
FTensor::Index< 'i', SPACE_DIM > i
MoFEMErrorCode setFieldOrder(const std::string field_name, const int order, const Range *ents=NULL)
Set field order.
default operator for EDGE element
SkeletonFE(MoFEM::Interface &m_field, CommonData &elem_data)
static MoFEMErrorCode Initialize(int *argc, char ***args, const char file[], const char help[])
Initializes the MoFEM database PETSc, MOAB and MPI.
#define CATCH_ERRORS
Catch errors.
MatrixDouble & getN(const FieldApproximationBase base)
get base functions this return matrix (nb. of rows is equal to nb. of Gauss pts, nb....
FTensor::Tensor1< FTensor::PackPtr< double *, Tensor_Dim >, Tensor_Dim > getFTensor1N(FieldApproximationBase base)
Get base functions for Hdiv/Hcurl spaces.
@ AINSWORTH_LEGENDRE_BASE
Ainsworth Cole (Legendre) approx. base .
boost::ptr_deque< UserDataOperator > & getOpPtrVector()
Use to push back operator for row operator.
@ HCURL
field with continuous tangents
PetscErrorCode PetscOptionsGetEList(PetscOptions *, const char pre[], const char name[], const char *const *list, PetscInt next, PetscInt *value, PetscBool *set)
OpCalculateHOJacForFaceImpl< 2 > OpCalculateHOJacForFace
default operator for Face element
FieldApproximationBase
approximation base
@ MOFEM_ATOM_TEST_INVALID
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Make Hdiv space from Hcurl space in 2d.
MoFEMErrorCode setUp(const PetscBool is_partitioned=PETSC_TRUE)
Setup problem.
PetscErrorCode DMoFEMLoopFiniteElements(DM dm, const char fe_name[], MoFEM::FEMethod *method, CacheTupleWeakPtr cache_ptr=CacheTupleSharedPtr())
Executes FEMethod for finite elements in DM.
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
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 ...