static char help[] =
"...\n\n";
int main(
int argc,
char *argv[]) {
try {
moab::Core mb_instance;
moab::Interface &moab = mb_instance;
int rank;
MPI_Comm_rank(PETSC_COMM_WORLD, &rank);
PetscBool flg = PETSC_TRUE;
#if PETSC_VERSION_GE(3, 6, 4)
255, &flg);
#else
CHKERR PetscOptionsGetString(PETSC_NULL, PETSC_NULL,
"-my_file",
#endif
if (flg != PETSC_TRUE) {
SETERRQ(PETSC_COMM_SELF, 1, "*** ERROR -my_file (MESH FILE NEEDED)");
}
#if PETSC_VERSION_GE(3, 6, 4)
CHKERR PetscOptionsGetInt(PETSC_NULL,
"",
"-my_order", &
order, PETSC_NULL);
#else
CHKERR PetscOptionsGetInt(PETSC_NULL, PETSC_NULL,
"-my_order", &
order,
PETSC_NULL);
#endif
const char *option;
option = "";
MB_TAG_DENSE);
MB_TAG_DENSE);
MB_TAG_DENSE);
auto check_axpy = [&]() {
CHKERR fb->setField(+1, MBVERTEX,
"FIELD_A");
CHKERR fb->setField(-2, MBVERTEX,
"FIELD_B");
CHKERR fb->fieldAxpy(+0.5,
"FIELD_B",
"FIELD_A");
CHKERR fb->fieldScale(-0.5,
"FIELD_B");
for (auto dit : *dofs_ptr) {
auto check = [&](const std::string name, const double expected) {
if (dit->getName() == name) {
cout << name << " " << dit->getFieldData() << " " << expected
<< endl;
if (dit->getFieldData() != expected)
"Wrong DOF value 0 != %4.3e for %s", dit->getFieldData(),
boost::lexical_cast<std::string>(*dit).c_str());
}
};
}
};
auto check_set_vertex = [&]() {
auto set_distance = [&](
VectorAdaptor &&field_data,
double *xcoord,
double *ycoord, double *zcoord) {
xcoord, ycoord, zcoord);
field_data[0] = sqrt(t_coord(
i) * t_coord(
i));
};
CHKERR fb->setVertexDofs(set_distance,
"FIELD_A");
TestMethod(moab::Interface &moab) :
EntityMethod(), moab(moab) {}
if (entPtr->getEntType() == MBVERTEX) {
CHKERR moab.get_coords(&
v, 1, &coords[0]);
if (std::abs(entPtr->getEntFieldData()[0] - norm_2(coords)) >
std::numeric_limits<double>::epsilon())
"Wrong vals %3.4e != %3.4e",
entPtr->getEntFieldData()[0], norm_2(coords));
}
}
private:
moab::Interface &moab;
};
TestMethod method(moab);
};
auto check_lambda = [&]() {
if (it->getMeshsetId() != 1)
continue;
CHKERR it->getMeshsetIdEntitiesByType(m_field.
get_moab(), MBTRI, ents,
true);
meshset_ents.merge(nodes);
}
CHKERR fb->setField(+1, MBVERTEX,
"FIELD_A");
CHKERR fb->setField(-1, MBVERTEX,
"FIELD_B");
CHKERR fb->setField(1.5, MBVERTEX,
"FIELD_C");
auto field_axpy = [&](const double val_y, const double val_x) {
return 0.5 * val_x + val_y;
};
auto field_scale = [&](const double val) { return -0.5 * val; };
auto vector_times_scalar_field =
[&](boost::shared_ptr<FieldEntity> ent_ptr_y,
boost::shared_ptr<FieldEntity> ent_ptr_x) {
auto x_data = ent_ptr_x->getEntFieldData();
auto y_data = ent_ptr_y->getEntFieldData();
const auto size_y = y_data.size();
for (size_t dd = 0; dd != size_y; ++dd)
y_data[dd] *= x_data[0];
};
CHKERR fb->fieldLambdaOnValues(field_axpy,
"FIELD_B",
"FIELD_A",
true);
CHKERR fb->fieldLambdaOnValues(field_scale,
"FIELD_B", &meshset_ents);
CHKERR fb->fieldLambdaOnEntities(vector_times_scalar_field,
"FIELD_C",
"FIELD_B", true, &meshset_ents);
constexpr double eps = std::numeric_limits<double>::epsilon();
for (auto dit : *dofs_ptr) {
auto check = [&](const std::string name, const double expected) {
if (dit->getName() == name && dit->getEntType() == MBVERTEX) {
cout << name << " " << dit->getFieldData() << " " << expected
<< endl;
if (abs(dit->getFieldData() - expected) >
eps)
"Wrong DOF value %4.3e != %4.3e for %s",
dit->getFieldData(), expected,
boost::lexical_cast<std::string>(*dit).c_str());
}
};
if (meshset_ents.find(dit->getEnt()) != meshset_ents.end()) {
CHKERR check(
"FIELD_B", 0.75);
}
}
};
}
return 0;
}
#define CATCH_ERRORS
Catch errors.
@ AINSWORTH_LEGENDRE_BASE
Ainsworth Cole (Legendre) approx. base .
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
@ MOFEM_ATOM_TEST_INVALID
#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 ...
virtual const DofEntity_multiIndex * get_dofs() const =0
Get the dofs object.
virtual MoFEMErrorCode build_fields(int verb=DEFAULT_VERBOSITY)=0
virtual MoFEMErrorCode set_field_order(const EntityHandle meshset, const EntityType type, const std::string &name, const ApproximationOrder order, int verb=DEFAULT_VERBOSITY)=0
Set order approximation of the entities in the field.
virtual MoFEMErrorCode add_ents_to_field_by_type(const Range &ents, const EntityType type, const std::string &name, int verb=DEFAULT_VERBOSITY)=0
Add entities to field meshset.
virtual MoFEMErrorCode loop_entities(const std::string field_name, EntityMethod &method, Range const *const ents=nullptr, int verb=DEFAULT_VERBOSITY)=0
Loop over field entities.
#define _IT_CUBITMESHSETS_BY_BCDATA_TYPE_FOR_LOOP_(MESHSET_MANAGER, CUBITBCTYPE, IT)
Iterator that loops over a specific Cubit MeshSet in a moFEM field.
FTensor::Index< 'i', SPACE_DIM > i
const double v
phase velocity of light in medium (cm/ns)
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
VectorBoundedArray< double, 3 > VectorDouble3
VectorShallowArrayAdaptor< double > VectorAdaptor
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
implementation of Data Operators for Forces and Sources
virtual moab::Interface & get_moab()=0
virtual MoFEMErrorCode add_field(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_EXCL, int verb=DEFAULT_VERBOSITY)=0
Add field.
static MoFEMErrorCode Initialize(int *argc, char ***args, const char file[], const char help[])
Initializes the MoFEM database PETSc, MOAB and MPI.
static MoFEMErrorCode Finalize()
Checks for options to be called at the conclusion of the program.
Deprecated interface functions.
Data structure to exchange data between mofem and User Loop Methods on entities.
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface refernce to pointer of interface.