![Logo](MoFEMLogo.png) |
| v0.14.0
|
Go to the documentation of this file.
10 *iface =
const_cast<Simple *
>(
this);
16 static_assert(DIM == 2 || DIM == 3,
"not implemented");
26 boost::make_shared<ParentFiniteElementAdjacencyFunctionSkeleton<1>>(
40 boost::make_shared<ParentFiniteElementAdjacencyFunctionSkeleton<2>>(
59 return setSkeletonAdjacency<2>(fe_name);
61 return setSkeletonAdjacency<3>(fe_name);
68 static_assert(DIM == 2 || DIM == 3,
"not implemented");
77 boost::make_shared<ParentFiniteElementAdjacencyFunction<3>>(
80 boost::make_shared<ParentFiniteElementAdjacencyFunction<2>>(
108 boost::make_shared<ParentFiniteElementAdjacencyFunction<2>>(
111 boost::make_shared<ParentFiniteElementAdjacencyFunction<1>>(
134 return setParentAdjacency<2>();
136 return setParentAdjacency<3>();
153 return setSkeletonAdjacency<2>(fe_name);
155 return setSkeletonAdjacency<3>(fe_name);
164 bitLevelMask(
BitRefLevel().set()), meshSet(0), boundaryMeshset(0),
165 skeletonMeshset(0), nameOfProblem(
"SimpleProblem"), domainFE(
"dFE"),
166 boundaryFE(
"bFE"), skeletonFE(
"sFE"), dIm(-1), addSkeletonFE(false),
167 addBoundaryFE(false), addParentAdjacencies(false),
173 PetscLogEventRegister(
"SimpBuildFEs", 0,
181 PetscBool flg = PETSC_TRUE;
183 ierr = PetscOptionsBegin(PETSC_COMM_WORLD,
"",
"Simple interface options",
186 ierr = PetscOptionsString(
"-file_name",
"file name",
"",
"mesh.h5m",
189 ierr = PetscOptionsEnd();
213 if (nb_ents_3d > 0) {
219 if (nb_ents_2d > 0) {
241 MOFEM_LOG(
"WORLD", Sev::warning) <<
"BitRefLevel is none and not set";
257 "PARALLEL_RESOLVE_SHARED_ENTS;"
258 "PARTITION=PARALLEL_PARTITION;",
267 const TagType tag_type,
const enum MoFEMTypes bh,
272 CHKERR m_field.
add_field(name, space, base, nb_of_cooficients, tag_type, bh,
285 const TagType tag_type,
const enum MoFEMTypes bh,
303 const TagType tag_type,
const enum MoFEMTypes bh,
307 CHKERR m_field.
add_field(name, space, base, nb_of_cooficients, tag_type, bh,
313 "NOFIELD space for boundary filed not implemented in Simple interface");
321 const TagType tag_type,
const enum MoFEMTypes bh,
326 CHKERR m_field.
add_field(name, space, base, nb_of_cooficients, tag_type, bh,
332 "NOFIELD space for boundary filed not implemented in Simple interface");
341 const TagType tag_type,
const enum MoFEMTypes bh,
346 CHKERR m_field.
add_field(name, space, base, nb_of_cooficients, tag_type, bh,
358 auto remove_field_from_list = [&](
auto &vec) {
359 auto it = std::find(vec.begin(), vec.end(), name);
373 auto remove_field_from_list = [&](
auto &vec) {
374 auto it = std::find(vec.begin(), vec.end(), name);
387 auto remove_field_from_list = [&](
auto &vec) {
388 auto it = std::find(vec.begin(), vec.end(), name);
402 auto clear_rows_and_cols = [&](
auto &fe_name) {
406 ->get<FiniteElement_name_mi_tag>();
407 auto it_fe = fe_by_name.find(fe_name);
408 if (it_fe != fe_by_name.end()) {
412 "modification unsuccessful");
416 "modification unsuccessful");
427 auto add_fields = [&](
auto &fe_name,
auto &fields) {
429 for (
auto &field : fields) {
437 auto add_data_fields = [&](
auto &fe_name,
auto &fields) {
439 for (
auto &field : fields)
496 ents == NULL ?
false :
true);
509 auto add_ents_to_field = [&](
auto meshset,
auto dim,
auto &fields) {
511 for (
auto &field : fields) {
513 CHKERR comm_interface_ptr->synchroniseFieldEntities(field, 0);
518 auto make_no_field_ents = [&](
auto &fields) {
520 for (
auto &field : fields) {
521 std::array<double, 6> coords = {0, 0, 0, 0, 0, 0};
523 CHKERR comm_interface_ptr->makeFieldEntitiesMultishared(field, 0);
534 std::set<std::string> nofield_fields;
536 nofield_fields.insert(
f);
538 nofield_fields.insert(
f);
540 CHKERR make_no_field_ents(nofield_fields);
546 const auto f = std::get<0>(
t);
547 const auto order = std::get<1>(
t);
550 <<
"Set order to field " <<
f <<
" order " <<
order;
552 if (std::get<3>(
t)) {
554 <<
"To ents: " << std::endl
555 << std::get<2>(
t) << std::endl;
559 if (std::get<3>(
t)) {
564 auto f_ptr = get_field_ptr(
f);
566 if (f_ptr->getSpace() ==
H1) {
574 for (
auto d = 1;
d <=
dIm; ++
d) {
575 for (EntityType
t = CN::TypeDimensionMap[
d].first;
576 t <= CN::TypeDimensionMap[
d].second; ++
t) {
607 for (std::vector<std::string>::iterator fit =
otherFEs.begin();
678 const auto problem_name = problem_ptr->
getName();
729 const std::string col_field)
const {
740 ParallelComm *pcomm =
743 auto get_skin = [&](
auto meshset) {
749 CHKERR pcomm->filter_pstatus(domain_ents,
750 PSTATUS_SHARED | PSTATUS_MULTISHARED,
751 PSTATUS_NOT, -1,
nullptr);
755 CHKERR skin.find_skin(0, domain_ents,
false, domain_skin);
756 CHKERR pcomm->filter_pstatus(domain_skin,
757 PSTATUS_SHARED | PSTATUS_MULTISHARED,
758 PSTATUS_NOT, -1,
nullptr);
762 auto create_boundary_meshset = [&](
auto &&domain_skin) {
773 moab::Interface::UNION);
788 ParallelComm *pcomm =
791 auto create_skeleton_meshset = [&](
auto meshset) {
797 Range boundary_ents, skeleton_ents;
799 dIm - 1, boundary_ents);
804 domain_ents,
dIm - 1,
false, skeleton_ents, moab::Interface::UNION);
805 skeleton_ents = subtract(skeleton_ents, boundary_ents);
806 CHKERR pcomm->filter_pstatus(skeleton_ents, PSTATUS_NOT_OWNED, PSTATUS_NOT,
821 MOFEM_LOG(
"WORLD", Sev::verbose) <<
"Exchange ghost cells";
823 ParallelComm *pcomm =
834 for (
auto d =
dIm - 1;
d >= 0; --
d) {
836 moab::Interface::UNION);
838 CHKERR pcomm->filter_pstatus(shared, PSTATUS_SHARED | PSTATUS_MULTISHARED,
839 PSTATUS_OR, -1, &shared);
840 Tag part_tag = pcomm->part_tag();
841 CHKERR pcomm->exchange_tags(part_tag, shared);
#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.
virtual MoFEMErrorCode modify_finite_element_adjacency_table(const std::string &fe_name, const EntityType type, ElementAdjacencyFunct function)=0
modify finite element table, only for advanced user
std::vector< std::string > boundaryFields
boundary fields
boost::shared_ptr< ParentFiniteElementAdjacencyFunction< 1 > > parentAdjFunctionDim1
MoFEMErrorCode setSkeletonAdjacency(int dim=-1, std::string fe_name="")
Set the skeleton adjacency object.
#define MYPCOMM_INDEX
default communicator number PCOMM
std::vector< std::string > otherFEs
Other finite elements.
bool addParentAdjacencies
If set to true parent adjacencies are build.
bool addBoundaryFE
Add boundary FE.
#define MOFEM_LOG_SEVERITY_SYNC(comm, severity)
Synchronise "SYNC" on curtain severity level.
MoFEMErrorCode exchangeGhostCells()
std::list< std::tuple< std::string, int, Range, bool > > fieldsOrder
virtual MoFEMErrorCode modify_problem_ref_level_set_bit(const std::string &name_problem, const BitRefLevel &bit)=0
set ref level for problem
#define MOFEM_LOG_CHANNEL(channel)
Set and reset channel.
Problem manager is used to build and partition problems.
Simple(const MoFEM::Core &core)
boost::shared_ptr< ParentFiniteElementAdjacencyFunction< 2 > > parentAdjFunctionDim2
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
std::vector< std::string > noFieldFields
NOFIELD field name.
EntityHandle meshSet
domain meshset
virtual MPI_Comm & get_comm() const =0
virtual MoFEMErrorCode modify_finite_element_add_field_row(const std::string &fe_name, const std::string name_row)=0
set field row which finite element use
SmartPetscObj< DM > getDM()
Return smart DM object.
MoFEMErrorCode createBoundaryMeshset()
virtual const Field * get_field_structure(const std::string &name, enum MoFEMTypes bh=MF_EXIST) const =0
get field structure
char meshFileName[255]
mesh file name
PetscLogEvent MOFEM_EVENT_SimpleKSPSolve
MoFEMErrorCode buildProblem()
Build problem.
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
virtual MoFEMErrorCode create_vertices_and_add_to_field(const std::string name, const double coords[], int size, int verb=DEFAULT_VERBOSITY)=0
Create a vertices and add to field object.
PetscLogEvent MOFEM_EVENT_SimpleBuildFields
PetscObject getPetscObject(T obj)
MoFEMErrorCode loadFile(const std::string options, const std::string mesh_file_name, LoadFileFunc loadFunc=defaultLoadFileFunc)
Load mesh file.
PetscErrorCode DMSetUp_MoFEM(DM dm)
int dIm
dimension of problem
virtual const FiniteElement_multiIndex * get_finite_elements() const =0
Get the finite elements object.
bool addSkeletonFE
Add skeleton FE.
std::string nameOfProblem
problem name
boost::shared_ptr< ParentFiniteElementAdjacencyFunctionSkeleton< 2 > > parentAdjSkeletonFunctionDim2
Simple interface for fast problem set-up.
virtual MoFEMErrorCode modify_problem_mask_ref_level_set_bit(const std::string &name_problem, const BitRefLevel &bit)=0
set dof mask ref level for problem
#define THROW_MESSAGE(msg)
Throw MoFEM exception.
virtual MoFEMErrorCode add_ents_to_field_by_dim(const Range &ents, const int dim, const std::string &name, int verb=DEFAULT_VERBOSITY)=0
Add entities to field meshset.
PetscErrorCode DMMoFEMAddElement(DM dm, std::string fe_name)
add element to dm
Deprecated interface functions.
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.
MoFEMErrorCode removeDomainField(const std::string &name)
Remove field form domain.
PetscLogEvent MOFEM_EVENT_SimpleLoadMesh
virtual MoFEMErrorCode add_broken_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.
MoFEMErrorCode deleteDM()
Delete dm.
MoFEMErrorCode getOptions()
get options
virtual MoFEMErrorCode rebuild_database(int verb=DEFAULT_VERBOSITY)=0
Clear database and initialize it once again.
BitRefLevel bitLevelMask
BitRefLevel of the problem.
std::vector< std::string > dataFields
Data fields.
FieldSpace
approximation spaces
#define CHKERR
Inline error check.
BitRefLevel bitAdjEntMask
bit ref level for parent parent
virtual MoFEMErrorCode add_finite_element(const std::string &fe_name, enum MoFEMTypes bh=MF_EXCL, int verb=DEFAULT_VERBOSITY)=0
add finite element
virtual moab::Interface & get_moab()=0
implementation of Data Operators for Forces and Sources
virtual MoFEMErrorCode modify_finite_element_add_field_col(const std::string &fe_name, const std::string name_row)=0
set field col which finite element use
MoFEMErrorCode deleteFiniteElements()
Delete finite elements.
auto createDM(MPI_Comm comm, const std::string dm_type_name)
Creates smart DM object.
std::string domainFE
domain finite element
virtual MoFEMErrorCode build_finite_elements(int verb=DEFAULT_VERBOSITY)=0
Build finite elements.
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.
int FieldCoefficientsNumber
Number of field coefficients.
PetscLogEvent MOFEM_EVENT_SimpleBuildFiniteElements
boost::function< MoFEMErrorCode(Interface &, const char *, const char *)> LoadFileFunc
MoFEMErrorCode createSkeletonMeshset()
std::string boundaryFE
boundary finite element
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.
virtual MoFEMErrorCode delete_finite_element(const std::string name, int verb=DEFAULT_VERBOSITY)=0
delete finite element from mofem database
virtual int get_comm_size() const =0
boost::shared_ptr< ParentFiniteElementAdjacencyFunction< 3 > > parentAdjFunctionDim3
BitRefLevel bitAdjParentMask
bit ref level for parent parent
std::string skeletonFE
skeleton finite element
MoFEMErrorCode removeBoundaryField(const std::string &name)
Remove field form boundary.
@ MOFEM_OPERATION_UNSUCCESSFUL
BitRefLevel bitAdjEnt
bit ref level for parent
virtual MoFEMErrorCode add_ents_to_finite_element_by_dim(const EntityHandle entities, const int dim, const std::string &name, const bool recursive=true)=0
add entities to finite element
constexpr double t
plate stiffness
PetscErrorCode DMMoFEMCreateMoFEM(DM dm, MoFEM::Interface *m_field_ptr, const char problem_name[], const MoFEM::BitRefLevel bit_level, const MoFEM::BitRefLevel bit_mask=MoFEM::BitRefLevel().set())
Must be called by user to set MoFEM data structures.
MoFEMErrorCode addDataField(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.
constexpr auto field_name
@ AINSWORTH_BERNSTEIN_BEZIER_BASE
MoFEMErrorCode setFieldOrder(const std::string field_name, const int order, const Range *ents=NULL)
Set field order.
MoFEMErrorCode buildFiniteElements()
Build finite elements.
base class for all interface classes
const Tensor2_symmetric_Expr< const ddTensor0< T, Dim, i, j >, typename promote< T, double >::V, Dim, i, j > dd(const Tensor0< T * > &a, const Index< i, Dim > index1, const Index< j, Dim > index2, const Tensor1< int, Dim > &d_ijk, const Tensor1< double, Dim > &d_xyz)
#define MOFEM_TAG_AND_LOG(channel, severity, tag)
Tag and log in channel.
#define MOFEM_LOG(channel, severity)
Log.
std::vector< std::string > domainFields
domain fields
MoFEMErrorCode addDomainBrokenField(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 broken field on domain.
std::vector< std::string > skeletonFields
fields on the skeleton
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
PetscErrorCode DMMoFEMGetProblemPtr(DM dm, const MoFEM::Problem **problem_ptr)
Get pointer to problem data structure.
BitRefLevel bitLevel
BitRefLevel of the problem.
MoFEMErrorCode defineFiniteElements()
Define finite elements.
int getDim() const
Get the problem dimension.
MoFEMErrorCode buildFields()
Build fields.
virtual MoFEMErrorCode build_fields(int verb=DEFAULT_VERBOSITY)=0
BitRefLevel bitAdjParent
bit ref level for parent
std::vector< std::string > noFieldDataFields
NOFIELD field name.
virtual MoFEMErrorCode modify_finite_element_add_field_data(const std::string &fe_name, const std::string name_filed)=0
set finite element field data
FieldApproximationBase
approximation base
MoFEMTypes
Those types control how functions respond on arguments, f.e. error handling.
MoFEMErrorCode addBoundaryField(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 boundary.
MoFEMErrorCode reSetUp(bool only_dm=false)
Rebuild internal MoFEM data structures.
multi_index_container< boost::shared_ptr< FiniteElement >, indexed_by< hashed_unique< tag< FiniteElement_Meshset_mi_tag >, member< FiniteElement, EntityHandle, &FiniteElement::meshset > >, hashed_unique< tag< BitFEId_mi_tag >, const_mem_fun< FiniteElement, BitFEId, &FiniteElement::getId >, HashBit< BitFEId >, EqBit< BitFEId > >, ordered_unique< tag< FiniteElement_name_mi_tag >, const_mem_fun< FiniteElement, boost::string_ref, &FiniteElement::getNameRef > > > > FiniteElement_multiIndex
MultiIndex for entities for FiniteElement.
EntityHandle boundaryMeshset
meshset with boundary
boost::shared_ptr< ParentFiniteElementAdjacencyFunctionSkeleton< 1 > > parentAdjSkeletonFunctionDim1
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
virtual MoFEMErrorCode build_adjacencies(const Range &ents, int verb=DEFAULT_VERBOSITY)=0
build adjacencies
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
MoFEMErrorCode defineProblem(const PetscBool is_partitioned=PETSC_TRUE)
define problem
keeps basic data about problem
virtual bool check_finite_element(const std::string &name) const =0
Check if finite element is in database.
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.
SmartPetscObj< DM > dM
Discrete manager (interface to PETSc/MoFEM functions)
MoFEMErrorCode setParentAdjacency()
MoFEMErrorCode setUp(const PetscBool is_partitioned=PETSC_TRUE)
Setup problem.
const std::string getProblemName() const
Get the Problem Name.
MoFEMErrorCode addFieldToEmptyFieldBlocks(const std::string row_field, const std::string col_field) const
add empty block to problem
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
EntityHandle skeletonMeshset
skeleton meshset with boundary
PetscLogEvent MOFEM_EVENT_SimpleSetUP
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
@ NOFIELD
scalar or vector of scalars describe (no true field)
PetscErrorCode DMMoFEMSetIsPartitioned(DM dm, PetscBool is_partitioned)
PetscLogEvent MOFEM_EVENT_SimpleBuildProblem
MoFEMErrorCode removeSkeletonField(const std::string &name)
Remove field form skeleton.