v0.13.1
Files | Classes
BitRefManager

Managing BitRefLevels. More...

Collaboration diagram for BitRefManager:

Files

file  BitRefManager.hpp
 Interface managing BitRefLevels.
 

Classes

struct  MoFEM::BitRefManager
 Managing BitRefLevels. More...
 
struct  MoFEM::CommInterface
 Managing BitRefLevels. More...
 

Setting and shifting bits

MoFEMErrorCode MoFEM::BitRefManager::setBitRefLevel (const Range &ents, const BitRefLevel bit, const bool only_tets=true, int verb=0, Range *adj_ents_ptr=nullptr) const
 add entities to database and set bit ref level More...
 
MoFEMErrorCode MoFEM::BitRefManager::setElementsBitRefLevel (const Range &ents, const BitRefLevel bit=BitRefLevel(), int verb=QUIET) const
 add entities to database and set bit ref level More...
 
MoFEMErrorCode MoFEM::BitRefManager::setEntitiesBitRefLevel (const Range &ents, const BitRefLevel bit=BitRefLevel(), int verb=QUIET) const
 add entities to database and set bit ref level More...
 
MoFEMErrorCode MoFEM::BitRefManager::setBitLevelToMeshset (const EntityHandle meshset, const BitRefLevel bit, int verb=0) const
 
MoFEMErrorCode MoFEM::BitRefManager::addBitRefLevel (const Range &ents, const BitRefLevel bit, int verb=QUIET) const
 add bit ref level to ref entity More...
 
MoFEMErrorCode MoFEM::BitRefManager::setNthBitRefLevel (const int n, const bool b, int verb=QUIET) const
 Set nth bit ref level to all entities in databse. More...
 
MoFEMErrorCode MoFEM::BitRefManager::shiftLeftBitRef (const int shift, const BitRefLevel mask=BitRefLevel().set(), int verb=DEFAULT_VERBOSITY) const
 left shift bit ref levelthis results of deletion of entities on far left side More...
 
MoFEMErrorCode MoFEM::BitRefManager::shiftRightBitRef (const int shift, const BitRefLevel mask=BitRefLevel().set(), int verb=DEFAULT_VERBOSITY, MoFEMTypes mf=MF_ZERO) const
 right shift bit ref level More...
 
MoFEMErrorCode MoFEM::BitRefManager::setFieldEntitiesBitRefLevel (const std::string field_name, const BitRefLevel bit=BitRefLevel(), int verb=QUIET) const
 Set the bit ref level to entities in the field meshset. More...
 
MoFEMErrorCode MoFEM::BitRefManager::setBitRefLevelByDim (const EntityHandle meshset, const int dim, const BitRefLevel bit, int verb=QUIET) const
 Set the Bit Ref Level By Dim object. More...
 
MoFEMErrorCode MoFEM::BitRefManager::setBitRefLevelByType (const EntityHandle meshset, const EntityType type, const BitRefLevel bit, int verb=QUIET) const
 Set the Bit Ref Level By Type object. More...
 
MoFEMErrorCode MoFEM::BitRefManager::addToDatabaseBitRefLevelByType (const EntityType type, const BitRefLevel bit, const BitRefLevel mask=BitRefLevel().set(), int verb=QUIET) const
 Add entities which exist in MoAB database, and have set appropiate BitRef level tag, to multi-indices in MoFEM. More...
 
MoFEMErrorCode MoFEM::BitRefManager::addToDatabaseBitRefLevelByDim (const int dim, const BitRefLevel bit, const BitRefLevel mask=BitRefLevel().set(), int verb=QUIET) const
 Add entities which exist in MoAB database, and have set appropiate BitRef level tag, to multi-indices in MoFEM. More...
 
MoFEMErrorCode MoFEM::BitRefManager::addBitRefLevelByDim (const EntityHandle meshset, const int dim, const BitRefLevel bit, int verb=QUIET) const
 add bit ref level by dimension More...
 
MoFEMErrorCode MoFEM::BitRefManager::setNthBitRefLevel (const Range &ents, const int n, const bool b, int verb=QUIET) const
 Set nth bit ref level. More...
 

Entity handlers by bit ref level

MoFEMErrorCode MoFEM::BitRefManager::getEntitiesByTypeAndRefLevel (const BitRefLevel bit, const BitRefLevel mask, const EntityType type, const EntityHandle meshset, int verb=0) const
 add all ents from ref level given by bit to meshset More...
 
MoFEMErrorCode MoFEM::BitRefManager::getEntitiesByTypeAndRefLevel (const BitRefLevel bit, const BitRefLevel mask, const EntityType type, Range &ents, int verb=0) const
 add all ents from ref level given by bit to meshset More...
 
MoFEMErrorCode MoFEM::BitRefManager::getEntitiesByDimAndRefLevel (const BitRefLevel bit, const BitRefLevel mask, const int dim, const EntityHandle meshset, int verb=0) const
 add all ents from ref level given by bit to meshset More...
 
MoFEMErrorCode MoFEM::BitRefManager::getEntitiesByDimAndRefLevel (const BitRefLevel bit, const BitRefLevel mask, const int dim, Range &ents, int verb=0) const
 add all ents from ref level given by bit to meshset More...
 
MoFEMErrorCode MoFEM::BitRefManager::getEntitiesByRefLevel (const BitRefLevel bit, const BitRefLevel mask, const EntityHandle meshset, const int verb=QUIET) const
 add all ents from ref level given by bit to meshset More...
 
MoFEMErrorCode MoFEM::BitRefManager::getEntitiesByRefLevel (const BitRefLevel bit, const BitRefLevel mask, Range &ents, const int verb=QUIET) const
 add all ents from ref level given by bit to meshset More...
 

Get adjacencies bit ref level

virtual MoFEMErrorCode MoFEM::BitRefManager::getAdjacenciesEquality (const EntityHandle from_entity, const int to_dimension, Range &adj_entities) const
 Get the adjacencies associated with a entity to entities of a specified dimension. More...
 
virtual MoFEMErrorCode MoFEM::BitRefManager::getAdjacenciesAny (const EntityHandle from_entity, const int to_dimension, Range &adj_entities) const
 Get the adjacencies associated with a entity to entities of a specified dimension. More...
 
virtual MoFEMErrorCode MoFEM::BitRefManager::getAdjacencies (const Problem *problem_ptr, const EntityHandle *from_entities, const int num_entities, const int to_dimension, Range &adj_entities, const int operation_type=moab::Interface::INTERSECT, const int verb=0) const
 Get the adjacencies associated with a entity to entities of a specified dimension. More...
 
virtual MoFEMErrorCode MoFEM::BitRefManager::getAdjacencies (const BitRefLevel bit, const EntityHandle *from_entities, const int num_entities, const int to_dimension, Range &adj_entities, const int operation_type=moab::Interface::INTERSECT, const int verb=0) const
 Get the adjacencies associated with a entity to entities of a specified dimension. More...
 

Update meshsets and ranges by children

MoFEMErrorCode MoFEM::BitRefManager::updateMeshsetByEntitiesChildren (const EntityHandle parent, const BitRefLevel &parent_bit, const BitRefLevel &parent_mask, const BitRefLevel &child_bit, const BitRefLevel &child_mask, const EntityHandle child, EntityType child_type, const bool recursive=false, int verb=0)
 Get child entities form meshset containing parent entities. More...
 
MoFEMErrorCode MoFEM::BitRefManager::updateMeshsetByEntitiesChildren (const EntityHandle parent, const BitRefLevel &child_bit, const EntityHandle child, EntityType child_type, const bool recursive=false, int verb=0)
 Get child entities form meshset containing parent entities. More...
 
MoFEMErrorCode MoFEM::BitRefManager::updateFieldMeshsetByEntitiesChildren (const BitRefLevel &child_bit, int verb=0)
 update fields meshesets by child entities More...
 
MoFEMErrorCode MoFEM::BitRefManager::updateFieldMeshsetByEntitiesChildren (const std::string name, const BitRefLevel &child_bit, int verb=0)
 update field meshset by child entities More...
 
MoFEMErrorCode MoFEM::BitRefManager::updateFiniteElementMeshsetByEntitiesChildren (const std::string name, const BitRefLevel &child_bit, const EntityType fe_ent_type, int verb=0)
 update finite element meshset by child entities More...
 
MoFEMErrorCode MoFEM::BitRefManager::updateRangeByChildren (const Range &parent, Range &child, MoFEMTypes bh=MF_ZERO)
 Update range by childresn. More...
 
MoFEMErrorCode MoFEM::BitRefManager::updateRangeByParent (const Range &parent, Range &child, MoFEMTypes bh=MF_ZERO)
 Update range by prents. More...
 
DEPRECATED MoFEMErrorCode MoFEM::BitRefManager::updateRange (const Range &parent, Range &child, MoFEMTypes bh=MF_ZERO)
 

Detailed Description

Managing BitRefLevels.

Function Documentation

◆ addBitRefLevel()

MoFEMErrorCode MoFEM::BitRefManager::addBitRefLevel ( const Range &  ents,
const BitRefLevel  bit,
int  verb = QUIET 
) const

add bit ref level to ref entity

Parameters
entsrange of entities
bitbit ref level
verbverbosity level
Returns
error code
Examples
hanging_node_approx.cpp.

Definition at line 538 of file BitRefManager.cpp.

540 {
541 MoFEM::Interface &m_field = cOre;
543 std::vector<const BitRefLevel *> ents_bits_vec;
544 CHKERR RefEntity::getBitRefLevel(m_field.get_moab(), ents, ents_bits_vec);
545 for (auto it : ents_bits_vec)
546 const_cast<BitRefLevel &>(*it) |= bit;
548}
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:429
#define CHKERR
Inline error check.
Definition: definitions.h:548
auto bit
set bit
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Types.hpp:51
virtual moab::Interface & get_moab()=0
Deprecated interface functions.
const BitRefLevel & getBitRefLevel() const
Get entity ref bit refinement signature.

◆ addBitRefLevelByDim()

MoFEMErrorCode MoFEM::BitRefManager::addBitRefLevelByDim ( const EntityHandle  meshset,
const int  dim,
const BitRefLevel  bit,
int  verb = QUIET 
) const

add bit ref level by dimension

Parameters
meshset
dimdimension of entities
bitadded bit
verbverbosity level
Returns
MoFEMErrorCode

Definition at line 550 of file BitRefManager.cpp.

553 {
554 MoFEM::Interface &m_field = cOre;
555 moab::Interface &moab = m_field.get_moab();
556 Range ents, adj;
558 CHKERR moab.get_entities_by_dimension(meshset, dim, ents, true);
559 for (int dd = dim - 1; dd >= 0; dd--) {
560 CHKERR moab.get_adjacencies(ents, dd, false, adj, moab::Interface::UNION);
561 }
562 ents.merge(adj);
563 if (verb == VERY_NOISY) {
565 MOFEM_LOG("BitRefSelf", Sev::noisy) << "Add add bit ref level by dim ";
566 }
567
568 CHKERR addBitRefLevel(ents, bit, verb);
570}
@ VERY_NOISY
Definition: definitions.h:225
const int dim
MoFEMErrorCode addBitRefLevel(const Range &ents, const BitRefLevel bit, int verb=QUIET) const
add bit ref level to ref entity
#define MOFEM_LOG(channel, severity)
Log.
Definition: LogManager.hpp:311
#define MOFEM_LOG_FUNCTION()
Set scope.
Definition: LogManager.hpp:328
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)
Definition: ddTensor0.hpp:33
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1965

◆ addToDatabaseBitRefLevelByDim()

MoFEMErrorCode MoFEM::BitRefManager::addToDatabaseBitRefLevelByDim ( const int  dim,
const BitRefLevel  bit,
const BitRefLevel  mask = BitRefLevel().set(),
int  verb = QUIET 
) const

Add entities which exist in MoAB database, and have set appropiate BitRef level tag, to multi-indices in MoFEM.

Note
Every entity, used for create DoFS, or elements has to have set BitRefLevel and be added to MoFEM database.
This functions add lower dimension entities by calling setEntitiesBitRefLevel
Parameters
dimdimension of entity
bitbit ref level
verbverbosity level
Returns
MoFEMErrorCode

Definition at line 477 of file BitRefManager.cpp.

479 {
481 Range ents;
483 CHKERR setBitRefLevel(ents, BitRefLevel(), false, verb);
485}
MoFEMErrorCode getEntitiesByDimAndRefLevel(const BitRefLevel bit, const BitRefLevel mask, const int dim, const EntityHandle meshset, int verb=0) const
add all ents from ref level given by bit to meshset
MoFEMErrorCode setBitRefLevel(const Range &ents, const BitRefLevel bit, const bool only_tets=true, int verb=0, Range *adj_ents_ptr=nullptr) const
add entities to database and set bit ref level

◆ addToDatabaseBitRefLevelByType()

MoFEMErrorCode MoFEM::BitRefManager::addToDatabaseBitRefLevelByType ( const EntityType  type,
const BitRefLevel  bit,
const BitRefLevel  mask = BitRefLevel().set(),
int  verb = QUIET 
) const

Add entities which exist in MoAB database, and have set appropiate BitRef level tag, to multi-indices in MoFEM.

Note
Every entity, used for create DoFS, or elements has to have set BitRefLevel and be added to MoFEM database.
This functions add lower dimension entities by calling setEntitiesBitRefLevel
Parameters
typeof entity
bitbit ref level
mask
verbverbosity level
Returns
MoFEMErrorCode

Definition at line 467 of file BitRefManager.cpp.

469 {
471 Range ents;
473 CHKERR setBitRefLevel(ents, BitRefLevel(), false, verb);
475}
MoFEMErrorCode getEntitiesByTypeAndRefLevel(const BitRefLevel bit, const BitRefLevel mask, const EntityType type, const EntityHandle meshset, int verb=0) const
add all ents from ref level given by bit to meshset

◆ getAdjacencies() [1/2]

MoFEMErrorCode MoFEM::BitRefManager::getAdjacencies ( const BitRefLevel  bit,
const EntityHandle *  from_entities,
const int  num_entities,
const int  to_dimension,
Range &  adj_entities,
const int  operation_type = moab::Interface::INTERSECT,
const int  verb = 0 
) const
virtual

Get the adjacencies associated with a entity to entities of a specified dimension.

bit ref level of adjacent entities is equal to bit ref level of adjacent entities

Definition at line 995 of file BitRefManager.cpp.

998 {
999 MoFEM::Interface &m_field = cOre;
1000 moab::Interface &moab(m_field.get_moab());
1002
1003 if (verb > VERBOSE) {
1005 MOFEM_LOG("BitRef", Sev::noisy) << "from: " << bit;
1006 }
1007
1008 CHKERR moab.get_adjacencies(from_entities, num_entities, to_dimension, false,
1009 adj_entities, operation_type);
1010 std::vector<BitRefLevel> bit_levels(adj_entities.size());
1011 CHKERR moab.tag_get_data(cOre.get_th_RefBitLevel(), adj_entities,
1012 &*bit_levels.begin());
1013 std::vector<BitRefLevel>::iterator b_it = bit_levels.begin();
1014 for (Range::iterator eit = adj_entities.begin(); eit != adj_entities.end();
1015 b_it++) {
1016
1017 if (verb > VERBOSE) {
1018 RefEntity adj_entity(m_field.get_basic_entity_data_ptr(), *eit);
1020 MOFEM_LOG("BitRef", Sev::noisy)
1021 << "to: " << adj_entity.getBitRefLevel() << " : " << adj_entity;
1022 }
1023
1024 if (!((*b_it) & bit).any()) {
1025 eit = adj_entities.erase(eit);
1026 } else {
1027 eit++;
1028 }
1029 }
1030 if (b_it != bit_levels.end()) {
1031 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "Data inconsistency");
1032 }
1034}
@ VERBOSE
Definition: definitions.h:222
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:44
RefEntityTmp< 0 > RefEntity
virtual boost::shared_ptr< BasicEntityData > & get_basic_entity_data_ptr()=0
Get pointer to basic entity data.
Tag get_th_RefBitLevel() const
Definition: Core.hpp:208

◆ getAdjacencies() [2/2]

MoFEMErrorCode MoFEM::BitRefManager::getAdjacencies ( const Problem problem_ptr,
const EntityHandle *  from_entities,
const int  num_entities,
const int  to_dimension,
Range &  adj_entities,
const int  operation_type = moab::Interface::INTERSECT,
const int  verb = 0 
) const
virtual

Get the adjacencies associated with a entity to entities of a specified dimension.

bit ref level of adjacent entities is equal to bit ref level of adjacent entities

Definition at line 984 of file BitRefManager.cpp.

987 {
989 BitRefLevel bit = problem_ptr->getBitRefLevel();
990 CHKERR getAdjacencies(bit, from_entities, num_entities, to_dimension,
991 adj_entities, operation_type);
993}
virtual MoFEMErrorCode getAdjacencies(const Problem *problem_ptr, const EntityHandle *from_entities, const int num_entities, const int to_dimension, Range &adj_entities, const int operation_type=moab::Interface::INTERSECT, const int verb=0) const
Get the adjacencies associated with a entity to entities of a specified dimension.

◆ getAdjacenciesAny()

MoFEMErrorCode MoFEM::BitRefManager::getAdjacenciesAny ( const EntityHandle  from_entity,
const int  to_dimension,
Range &  adj_entities 
) const
virtual

Get the adjacencies associated with a entity to entities of a specified dimension.

bit ref level of adjacent entities is any of bit ref level of adjacent entities

Definition at line 958 of file BitRefManager.cpp.

960 {
961 MoFEM::Interface &m_field = cOre;
962 moab::Interface &moab(m_field.get_moab());
964 BitRefLevel bit_from_entity;
965 CHKERR moab.tag_get_data(cOre.get_th_RefBitLevel(), &from_entity, 1,
966 &bit_from_entity);
967 CHKERR moab.get_adjacencies(&from_entity, 1, to_dimension, false,
968 adj_entities);
969 std::vector<BitRefLevel> bit_levels(adj_entities.size());
970 CHKERR moab.tag_get_data(cOre.get_th_RefBitLevel(), adj_entities,
971 &*bit_levels.begin());
972 std::vector<BitRefLevel>::iterator b_it = bit_levels.begin();
973 Range::iterator eit = adj_entities.begin();
974 for (; eit != adj_entities.end(); b_it++) {
975 if (!(bit_from_entity & (*b_it)).any()) {
976 eit = adj_entities.erase(eit);
977 } else {
978 eit++;
979 }
980 }
982}
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:460
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:453

◆ getAdjacenciesEquality()

MoFEMErrorCode MoFEM::BitRefManager::getAdjacenciesEquality ( const EntityHandle  from_entity,
const int  to_dimension,
Range &  adj_entities 
) const
virtual

Get the adjacencies associated with a entity to entities of a specified dimension.

bit ref level of adjacent entities is equal to bit ref level of adjacent entities

Definition at line 932 of file BitRefManager.cpp.

934 {
935 MoFEM::Interface &m_field = cOre;
936 moab::Interface &moab(m_field.get_moab());
938 BitRefLevel bit_from_entity;
939 CHKERR moab.tag_get_data(cOre.get_th_RefBitLevel(), &from_entity, 1,
940 &bit_from_entity);
941 CHKERR moab.get_adjacencies(&from_entity, 1, to_dimension, false,
942 adj_entities);
943 std::vector<BitRefLevel> bit_levels(adj_entities.size());
944 CHKERR moab.tag_get_data(cOre.get_th_RefBitLevel(), adj_entities,
945 &*bit_levels.begin());
946 auto b_it = bit_levels.begin();
947 auto eit = adj_entities.begin();
948 for (; eit != adj_entities.end(); b_it++) {
949 if (bit_from_entity != *b_it) {
950 eit = adj_entities.erase(eit);
951 } else {
952 eit++;
953 }
954 }
956}

◆ getEntitiesByDimAndRefLevel() [1/2]

MoFEMErrorCode MoFEM::BitRefManager::getEntitiesByDimAndRefLevel ( const BitRefLevel  bit,
const BitRefLevel  mask,
const int  dim,
const EntityHandle  meshset,
int  verb = 0 
) const

add all ents from ref level given by bit to meshset

Note
Entities NOT have to be added to MoFEM database
Parameters
BitRefLevelbitLevel
BitRefLevelmask
EntityTypedimension of entities
Return values
ents
Examples
hanging_node_approx.cpp.

Definition at line 829 of file BitRefManager.cpp.

831 {
832 MoFEM::Interface &m_field = cOre;
833 moab::Interface &moab(m_field.get_moab());
835 Range ents;
836 CHKERR getEntitiesByDimAndRefLevel(bit, mask, dim, ents, verb);
837 CHKERR moab.add_entities(meshset, ents);
839}

◆ getEntitiesByDimAndRefLevel() [2/2]

MoFEMErrorCode MoFEM::BitRefManager::getEntitiesByDimAndRefLevel ( const BitRefLevel  bit,
const BitRefLevel  mask,
const int  dim,
Range &  ents,
int  verb = 0 
) const

add all ents from ref level given by bit to meshset

Note
Entities NOT have to be added to MoFEM database
Parameters
BitRefLevelbitLevel
BitRefLevelmask
EntityTypedimension of entities
Return values
ents

Definition at line 841 of file BitRefManager.cpp.

843 {
844 MoFEM::Interface &m_field = cOre;
845 moab::Interface &moab(m_field.get_moab());
847 CHKERR moab.get_entities_by_dimension(0, dim, ents, false);
848 CHKERR filterEntitiesByRefLevel(bit, mask, ents, verb);
850}
MoFEMErrorCode filterEntitiesByRefLevel(const BitRefLevel bit, const BitRefLevel mask, Range &ents, int verb=QUIET) const
filter entities by bit ref level

◆ getEntitiesByRefLevel() [1/2]

MoFEMErrorCode MoFEM::BitRefManager::getEntitiesByRefLevel ( const BitRefLevel  bit,
const BitRefLevel  mask,
const EntityHandle  meshset,
const int  verb = QUIET 
) const

add all ents from ref level given by bit to meshset

Note
Entities NOT have to be added to MoFEM database
Parameters
BitRefLevelbitLevel
BitRefLevelmask
EntityHandlemeshset

Definition at line 852 of file BitRefManager.cpp.

855 {
856 MoFEM::Interface &m_field = cOre;
857 moab::Interface &moab(m_field.get_moab());
859 Range ents;
860 CHKERR getEntitiesByRefLevel(bit, mask, ents, verb);
861 CHKERR moab.add_entities(meshset, ents);
863}
MoFEMErrorCode getEntitiesByRefLevel(const BitRefLevel bit, const BitRefLevel mask, const EntityHandle meshset, const int verb=QUIET) const
add all ents from ref level given by bit to meshset

◆ getEntitiesByRefLevel() [2/2]

MoFEMErrorCode MoFEM::BitRefManager::getEntitiesByRefLevel ( const BitRefLevel  bit,
const BitRefLevel  mask,
Range &  ents,
const int  verb = QUIET 
) const

add all ents from ref level given by bit to meshset

Note
Entities NOT have to be added to MoFEM database
Parameters
BitRefLevelbitLevel
BitRefLevelmask
Return values
ents

Definition at line 865 of file BitRefManager.cpp.

868 {
869 MoFEM::Interface &m_field = cOre;
870 moab::Interface &moab(m_field.get_moab());
872 Range meshset_ents;
873 CHKERR moab.get_entities_by_type(0, MBENTITYSET, meshset_ents, false);
874 CHKERR moab.get_entities_by_handle(0, ents, false);
875 ents.merge(meshset_ents);
876 CHKERR filterEntitiesByRefLevel(bit, mask, ents, verb);
878}

◆ getEntitiesByTypeAndRefLevel() [1/2]

MoFEMErrorCode MoFEM::BitRefManager::getEntitiesByTypeAndRefLevel ( const BitRefLevel  bit,
const BitRefLevel  mask,
const EntityType  type,
const EntityHandle  meshset,
int  verb = 0 
) const

add all ents from ref level given by bit to meshset

Note
Entities NOT have to be added to MoFEM database
Parameters
BitRefLevelbitLevel
BitRefLevelmask
EntityTypetype of entities
Return values
EntityHandlemeshset
Examples
mesh_cut.cpp.

Definition at line 741 of file BitRefManager.cpp.

743 {
744 MoFEM::Interface &m_field = cOre;
745 moab::Interface &moab(m_field.get_moab());
747 Range ents;
748 CHKERR getEntitiesByTypeAndRefLevel(bit, mask, type, ents, verb);
749 CHKERR moab.add_entities(meshset, ents);
751}

◆ getEntitiesByTypeAndRefLevel() [2/2]

MoFEMErrorCode MoFEM::BitRefManager::getEntitiesByTypeAndRefLevel ( const BitRefLevel  bit,
const BitRefLevel  mask,
const EntityType  type,
Range &  ents,
int  verb = 0 
) const

add all ents from ref level given by bit to meshset

Note
Entities NOT have to be added to MoFEM database
Parameters
BitRefLevelbitLevel
BitRefLevelmask
EntityTypetype of entities
Return values
ents

Definition at line 818 of file BitRefManager.cpp.

820 {
821 MoFEM::Interface &m_field = cOre;
822 moab::Interface &moab(m_field.get_moab());
824 CHKERR moab.get_entities_by_type(0, type, ents, false);
825 CHKERR filterEntitiesByRefLevel(bit, mask, ents, verb);
827}

◆ setBitLevelToMeshset()

MoFEMErrorCode MoFEM::BitRefManager::setBitLevelToMeshset ( const EntityHandle  meshset,
const BitRefLevel  bit,
int  verb = 0 
) const

brief add meshset and set bit ref level

Parameters
EntityHandleMeshSet
BitRefLevelbitLevel

Definition at line 487 of file BitRefManager.cpp.

489 {
490 MoFEM::Interface &m_field = cOre;
491 auto ref_ents_ptr = m_field.get_ref_ents();
492 auto ref_fe_ptr = m_field.get_ref_finite_elements();
494 // Add ref entity
495 std::pair<RefEntity_multiIndex::iterator, bool> p_ent =
496 const_cast<RefEntity_multiIndex *>(ref_ents_ptr)
497 ->insert(boost::shared_ptr<RefEntity>(
498 new RefEntity(m_field.get_basic_entity_data_ptr(), meshset)));
499 *(const_cast<RefEntity *>(p_ent.first->get())->getBitRefLevelPtr()) |= bit;
500 // Add ref element
501 boost::shared_ptr<RefElement> fe_ptr =
502 boost::shared_ptr<RefElement>(new RefElement_MESHSET(*p_ent.first));
503 std::pair<RefElement_multiIndex::iterator, bool> p_fe =
504 const_cast<RefElement_multiIndex *>(ref_fe_ptr)->insert(fe_ptr);
505
507 MOFEM_LOG("BitRefSelf", Sev::noisy)
508 << "Add meshset as ref_ent " << **p_fe.first;
509
511}
multi_index_container< boost::shared_ptr< RefEntity >, indexed_by< ordered_unique< tag< Ent_mi_tag >, const_mem_fun< RefEntity, EntityHandle, &RefEntity::getEnt > >, ordered_non_unique< tag< Ent_Ent_mi_tag >, const_mem_fun< RefEntity, EntityHandle, &RefEntity::getParentEnt > >, ordered_non_unique< tag< Composite_EntType_and_ParentEntType_mi_tag >, composite_key< RefEntity, const_mem_fun< RefEntity, EntityType, &RefEntity::getEntType >, const_mem_fun< RefEntity, EntityType, &RefEntity::getParentEntType > > >, ordered_non_unique< tag< Composite_ParentEnt_And_EntType_mi_tag >, composite_key< RefEntity, const_mem_fun< RefEntity, EntityType, &RefEntity::getEntType >, const_mem_fun< RefEntity, EntityHandle, &RefEntity::getParentEnt > > > > > RefEntity_multiIndex
multi_index_container< boost::shared_ptr< RefElement >, indexed_by< ordered_unique< tag< Ent_mi_tag >, const_mem_fun< RefElement::interface_type_RefEntity, EntityHandle, &RefElement::getEnt > > > > RefElement_multiIndex
virtual const RefEntity_multiIndex * get_ref_ents() const =0
Get the ref ents object.
virtual const RefElement_multiIndex * get_ref_finite_elements() const =0
Get the ref finite elements object.

◆ setBitRefLevel()

MoFEMErrorCode MoFEM::BitRefManager::setBitRefLevel ( const Range &  ents,
const BitRefLevel  bit,
const bool  only_tets = true,
int  verb = 0,
Range *  adj_ents_ptr = nullptr 
) const

add entities to database and set bit ref level

This function set bit ref level, add entries to core database and create ref finite elements. Finite elements are create of entities in function argument, whereas all lower dimension entities are added as a field entities

Example:

EntityHandle meshset1; //contains ent1,ent2,ent3
BitRefLevel myLevel0;
myLevel0.set(0);
m_field.getInterface<BitRefManager>()->setBitRefLevelByDim(meshset1,3,myLevel0);
//refine meshset1 into meshset2 and get new ents which are ent4, ent5
EntityHandle meshset2; //contains ent1,ent2,ent3,ent4,ent5
BitRefLevel myLevel1;
myLevel1.set(1);
m_field.getInterface<BitRefManager>()->setBitRefLevelByDim(meshset2,3,myLevel1);
MoFEMErrorCode setBitRefLevelByDim(const EntityHandle meshset, const int dim, const BitRefLevel bit, int verb=QUIET) const
Set the Bit Ref Level By Dim object.
Managing BitRefLevels.

So entities 1,2,3 would be assigned to bit level 0 and 1
ent1[1,1,0,0,0,0,0], ent2[1,1,0,0,0,0,0], ent3[1,1,0,0,0,0,0],
and entities 4 and 5 are assigned to bit level 1 only
ent4[0,1,0,0,0,0,0], ent5[0,1,0,0,0,0,0]

Parameters
entsentities to set
bitbit refinement level
only_tetsonly add entities on tetrahedral (obsolete need to be fixed)
verbverbosity level
adj_ents_ptrif pointer is given, it is used to get adj entities by dimension/type
Returns
error code

Definition at line 293 of file BitRefManager.cpp.

296 {
297 MoFEM::Interface &m_field = cOre;
298 auto ref_ents_ptr = m_field.get_ref_ents();
299 auto ref_fe_ptr = m_field.get_ref_finite_elements();
301
303 MOFEM_LOG_C("BitRefSelf", Sev::noisy, "Number of entities to add %d",
304 ents.size());
305
306 CHKERR setElementsBitRefLevel(ents, bit, verb);
307
308 if (!ents.empty()) {
309 for (int d = 3; d >= 1; --d) {
310 Range dim_ents;
311 if (only_tets && d == 3) {
312 dim_ents = ents.subset_by_type(MBTET);
313 } else {
314 dim_ents = ents.subset_by_dimension(d);
315 }
316
318 MOFEM_LOG_C("BitRefSelf", Sev::noisy,
319 "\tNumber of dim %d entities to add %d", d, dim_ents.size());
320
321 if (!dim_ents.empty()) {
322 for (int dd = 0; dd < d; ++dd) {
323 Range adj_ents;
324
325 if (dd == 0) {
326 rval = m_field.get_moab().get_connectivity(ents, adj_ents, true);
327 // rval = m_field.get_moab().get_adjacencies(
328 // dim_ents, dd, true, adj_ents, moab::Interface::UNION);
329
330 } else {
331 if (adj_ents_ptr) {
332 if (dd == 1) {
333 adj_ents = adj_ents_ptr->subset_by_dimension(MBEDGE);
334 } else if (dd == 2) {
335 adj_ents = adj_ents_ptr->subset_by_dimension(MBTRI);
336 }
337 } else {
338 rval = m_field.get_moab().get_adjacencies(
339 dim_ents, dd, true, adj_ents, moab::Interface::UNION);
340 }
341 }
342
343 // rval = m_field.get_moab().get_adjacencies(
344 // dim_ents, dd, true, adj_ents, moab::Interface::UNION);
345
347 MOFEM_LOG_C("BitRefSelf", Sev::noisy,
348 "\tNumber of dim %d adj entities for dim %d to add %d", d,
349 dd, adj_ents.size());
350
351 if (rval == MB_MULTIPLE_ENTITIES_FOUND) {
352 auto log_message = [&](const auto sev) {
355 MOFEM_LOG("BitRefSelf", sev)
356 << "When get adjacencies moab return MB_MULTIPLE_ENTITIES_ "
357 "FOUND for dim = "
358 << dd << " and dim of entities " << d;
359 MOFEM_LOG_CHANNEL("BitRefSelf"); // reset channel
360 };
361
362 if (verb <= QUIET)
363 log_message(Sev::noisy);
364 else
365 log_message(Sev::warning);
366
367 rval = MB_SUCCESS;
368 }
370 for (Range::pair_iterator pit = adj_ents.pair_begin();
371 pit != adj_ents.pair_end(); ++pit) {
372 Range seed_ents_range;
373 // get first and last element of range
374 EntityHandle f = pit->first;
375 EntityHandle s = pit->second;
376 CHKERR SetBitRefLevelTool(m_field, bit, ref_ents_ptr, ref_fe_ptr)
377 .findEntsToAdd(f, s, seed_ents_range);
378 if (!seed_ents_range.empty())
379 CHKERR SetBitRefLevelTool(m_field, bit, ref_ents_ptr, ref_fe_ptr)
380 .addEntsToDatabase(seed_ents_range);
381 }
382 }
383 }
384 }
385 }
386
388}
#define MOFEM_LOG_C(channel, severity, format,...)
Definition: LogManager.hpp:314
@ QUIET
Definition: definitions.h:221
#define MOAB_THROW(err)
Check error code of MoAB function and throw MoFEM exception.
Definition: definitions.h:554
MoFEMErrorCode setElementsBitRefLevel(const Range &ents, const BitRefLevel bit=BitRefLevel(), int verb=QUIET) const
add entities to database and set bit ref level
#define MOFEM_LOG_CHANNEL(channel)
Set and reset channel.
Definition: LogManager.hpp:287
#define MOFEM_LOG_ATTRIBUTES(channel, bit)
Add attributes to channel.
Definition: LogManager.hpp:299
const Tensor1_Expr< const dTensor0< T, Dim, i >, typename promote< T, double >::V, Dim, i > d(const Tensor0< T * > &a, const Index< i, Dim > index, const Tensor1< int, Dim > &d_ijk, const Tensor1< double, Dim > &d_xyz)
Definition: dTensor0.hpp:27
static MoFEMErrorCodeGeneric< moab::ErrorCode > rval
Definition: Exceptions.hpp:85

◆ setBitRefLevelByDim()

MoFEMErrorCode MoFEM::BitRefManager::setBitRefLevelByDim ( const EntityHandle  meshset,
const int  dim,
const BitRefLevel  bit,
int  verb = QUIET 
) const

Set the Bit Ref Level By Dim object.

Note
In THIS variant only entities in range are added. And DO NOT create elements.
Parameters
meshset
dim
bit
verb
Returns
MoFEMErrorCode

Definition at line 513 of file BitRefManager.cpp.

516 {
517 MoFEM::Interface &m_field = cOre;
519 Range ents;
520 CHKERR m_field.get_moab().get_entities_by_dimension(meshset, dim, ents,
521 false);
522 CHKERR setBitRefLevel(ents, bit, false, verb);
524}

◆ setBitRefLevelByType()

MoFEMErrorCode MoFEM::BitRefManager::setBitRefLevelByType ( const EntityHandle  meshset,
const EntityType  type,
const BitRefLevel  bit,
int  verb = QUIET 
) const

Set the Bit Ref Level By Type object.

Note
In THIS variant only entities in range are added. And DO NOT create elements.
Parameters
meshset
type
bit
verb
Returns
MoFEMErrorCode
Examples
mesh_cut.cpp.

Definition at line 526 of file BitRefManager.cpp.

529 {
530 MoFEM::Interface &m_field = cOre;
532 Range ents;
533 CHKERR m_field.get_moab().get_entities_by_type(meshset, type, ents, false);
534 CHKERR setBitRefLevel(ents, bit, false, verb);
536}

◆ setElementsBitRefLevel()

MoFEMErrorCode MoFEM::BitRefManager::setElementsBitRefLevel ( const Range &  ents,
const BitRefLevel  bit = BitRefLevel(),
int  verb = QUIET 
) const

add entities to database and set bit ref level

Note
In THIS variant only entities in range are added and ref finite elements reated.
Parameters
ents
bit
only_tets
verb
Returns
MoFEMErrorCode setBitRefLevel

Definition at line 390 of file BitRefManager.cpp.

392 {
393 MoFEM::Interface &m_field = cOre;
394 auto ref_ents_ptr = m_field.get_ref_ents();
395 auto ref_fe_ptr = m_field.get_ref_finite_elements();
397
398 for (Range::const_pair_iterator pit = ents.pair_begin();
399 pit != ents.pair_end(); pit++) {
400 // get first and last element of range
401 EntityHandle f = pit->first;
402 EntityHandle s = pit->second;
403 Range seed_ents_range; // entities seeded not in database
404 // find ents to add
405 CHKERR SetBitRefLevelTool(m_field, bit, ref_ents_ptr, ref_fe_ptr)
406 .findEntsToAdd(f, s, seed_ents_range);
407 // add elements
408 if (!seed_ents_range.empty())
409 CHKERR SetBitRefLevelTool(m_field, bit, ref_ents_ptr, ref_fe_ptr)
410 .addEntsToDatabase(seed_ents_range);
411
412 Range seed_fe_range;
413 CHKERR SetBitRefLevelTool(m_field, bit, ref_ents_ptr, ref_fe_ptr)
414 .findElementsToAdd(f, s, seed_fe_range);
415 if (!seed_fe_range.empty()) {
416 CHKERR SetBitRefLevelTool(m_field, bit, ref_ents_ptr, ref_fe_ptr)
417 .addElementsToDatabase(seed_fe_range);
418 }
419 }
420
422 MOFEM_LOG("BitRefSelf", Sev::noisy)
423 << "Number of entities in databse " << ref_ents_ptr->size();
424 MOFEM_LOG("BitRefSelf", Sev::noisy)
425 << "Number of finite element entities in databse " << ref_fe_ptr->size();
426
428}

◆ setEntitiesBitRefLevel()

MoFEMErrorCode MoFEM::BitRefManager::setEntitiesBitRefLevel ( const Range &  ents,
const BitRefLevel  bit = BitRefLevel(),
int  verb = QUIET 
) const

add entities to database and set bit ref level

Note
In THIS variant only entities in range are added. And DO NOT create elements.
Parameters
ents
bit
only_tets
verb
Returns
MoFEMErrorCode setBitRefLevel

Definition at line 430 of file BitRefManager.cpp.

432 {
433 MoFEM::Interface &m_field = cOre;
434 auto ref_ents_ptr = m_field.get_ref_ents();
435 auto ref_fe_ptr = m_field.get_ref_finite_elements();
437
438 for (Range::const_pair_iterator pit = ents.pair_begin();
439 pit != ents.pair_end(); pit++) {
440 // get first and last element of range
441 EntityHandle f = pit->first;
442 EntityHandle s = pit->second;
443 Range seed_ents_range; // entities seeded not in database
444 // find ents to add
445 CHKERR SetBitRefLevelTool(m_field, bit, ref_ents_ptr, ref_fe_ptr)
446 .findEntsToAdd(f, s, seed_ents_range);
447 // add elements
448 if (!seed_ents_range.empty())
449 CHKERR SetBitRefLevelTool(m_field, bit, ref_ents_ptr, ref_fe_ptr)
450 .addEntsToDatabase(seed_ents_range);
451 }
453}

◆ setFieldEntitiesBitRefLevel()

MoFEMErrorCode MoFEM::BitRefManager::setFieldEntitiesBitRefLevel ( const std::string  field_name,
const BitRefLevel  bit = BitRefLevel(),
int  verb = QUIET 
) const

Set the bit ref level to entities in the field meshset.

Parameters
field_name
bit
verb
Returns
MoFEMErrorCode

Definition at line 455 of file BitRefManager.cpp.

456 {
457 MoFEM::Interface &m_field = cOre;
459 EntityHandle field_meshset = m_field.get_field_meshset(field_name);
460 Range field_ents;
461 CHKERR m_field.get_moab().get_entities_by_handle(field_meshset, field_ents,
462 true);
463 CHKERR setEntitiesBitRefLevel(field_ents, bit, verb);
465}
MoFEMErrorCode setEntitiesBitRefLevel(const Range &ents, const BitRefLevel bit=BitRefLevel(), int verb=QUIET) const
add entities to database and set bit ref level
constexpr auto field_name
virtual EntityHandle get_field_meshset(const std::string name) const =0
get field meshset

◆ setNthBitRefLevel() [1/2]

MoFEMErrorCode MoFEM::BitRefManager::setNthBitRefLevel ( const int  n,
const bool  b,
int  verb = QUIET 
) const

Set nth bit ref level to all entities in databse.

Parameters
nnth bit
bvalue to set
Returns
error code

Definition at line 587 of file BitRefManager.cpp.

588 {
589 MoFEM::Interface &m_field = cOre;
590 auto ref_ents_ptr = m_field.get_ref_ents();
592 auto hi_dit = ref_ents_ptr->end();
593 for (auto dit = ref_ents_ptr->begin(); dit != hi_dit; ++dit) {
594 (*const_cast<RefEntity *>(dit->get())->getBitRefLevelPtr())[n] = b;
595 if (verb >= VERY_VERBOSE) {
596 MOFEM_LOG("BitRefSelf", Sev::noisy) << "Set bit to " << **dit;
597 }
598 }
600}
static Index< 'n', 3 > n
@ VERY_VERBOSE
Definition: definitions.h:223

◆ setNthBitRefLevel() [2/2]

MoFEMErrorCode MoFEM::BitRefManager::setNthBitRefLevel ( const Range &  ents,
const int  n,
const bool  b,
int  verb = QUIET 
) const

Set nth bit ref level.

Note
This function modify bits only on entities in RefEntity_multiindex
Parameters
entsentities to set bit ref level
nnth bit
bvalue to set
Returns
error code

Definition at line 572 of file BitRefManager.cpp.

573 {
574 MoFEM::Interface &m_field = cOre;
576 std::vector<const BitRefLevel *> ents_bits_vec;
577 CHKERR RefEntity::getBitRefLevel(m_field.get_moab(), ents, ents_bits_vec);
578 for (auto it = ents_bits_vec.begin(); it != ents_bits_vec.end(); ++it) {
579 const_cast<BitRefLevel &>(**it)[n] = b;
580 }
581 if (verb == VERY_NOISY) {
582 MOFEM_LOG("BitRefSelf", Sev::noisy) << "Set bit to " << ents;
583 }
585}

◆ shiftLeftBitRef()

MoFEMErrorCode MoFEM::BitRefManager::shiftLeftBitRef ( const int  shift,
const BitRefLevel  mask = BitRefLevel().set(),
int  verb = DEFAULT_VERBOSITY 
) const

left shift bit ref levelthis results of deletion of entities on far left side

Note
Not implemented

Definition at line 602 of file BitRefManager.cpp.

604 {
606 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
608}
@ MOFEM_NOT_IMPLEMENTED
Definition: definitions.h:45

◆ shiftRightBitRef()

MoFEMErrorCode MoFEM::BitRefManager::shiftRightBitRef ( const int  shift,
const BitRefLevel  mask = BitRefLevel().set(),
int  verb = DEFAULT_VERBOSITY,
MoFEMTypes  mf = MF_ZERO 
) const

right shift bit ref level

Definition at line 610 of file BitRefManager.cpp.

612 {
613 MoFEM::Interface &m_field = cOre;
614 auto ref_ents_ptr = m_field.get_ref_ents();
616 RefEntity_change_right_shift right_shift(1, mask);
617 for (int ii = 0; ii < shift; ii++) {
618 // delete bits on the right which are shifted to zero
619 BitRefLevel delete_bits = BitRefLevel().set(0) & mask;
620 if (delete_bits.any()) {
621 CHKERR m_field.delete_ents_by_bit_ref(delete_bits, delete_bits, true,
622 verb, mf);
623 }
624 for (RefEntity_multiIndex::iterator ent_it = ref_ents_ptr->begin();
625 ent_it != ref_ents_ptr->end(); ent_it++) {
626 if (verb >= NOISY) {
628 MOFEM_LOG("BitRefSelf", Sev::noisy)
629 << (*ent_it)->getBitRefLevel() << " : ";
630 }
631 right_shift(const_cast<boost::shared_ptr<RefEntity> &>(*ent_it));
632 if (verb == VERY_NOISY) {
634 MOFEM_LOG("BitRefSelf", Sev::noisy) << (*ent_it)->getBitRefLevel();
635 }
636 }
637 }
639}
@ NOISY
Definition: definitions.h:224
virtual MoFEMErrorCode delete_ents_by_bit_ref(const BitRefLevel bit, const BitRefLevel mask, const bool remove_parent=false, int verb=DEFAULT_VERBOSITY, MoFEMTypes mf=MF_ZERO)=0
delete entities form mofem and moab database

◆ updateFieldMeshsetByEntitiesChildren() [1/2]

MoFEMErrorCode MoFEM::BitRefManager::updateFieldMeshsetByEntitiesChildren ( const BitRefLevel child_bit,
int  verb = 0 
)

update fields meshesets by child entities

Note
This calls updateMeshsetByEntitiesChildren for all entity types.

Definition at line 1072 of file BitRefManager.cpp.

1073 {
1074 MoFEM::Interface &m_field = cOre;
1075 moab::Interface &moab = m_field.get_moab();
1076 auto fields_ptr = m_field.get_fields();
1077 auto ref_ents_ptr = m_field.get_ref_ents();
1079
1080 for (auto &fit : (*fields_ptr)) {
1081
1082 const EntityHandle meshset = fit->getMeshset();
1083 Range parent_ents;
1084 CHKERR moab.get_entities_by_handle(meshset, parent_ents, true);
1085
1086 if (verb >= VERY_VERBOSE) {
1088 MOFEM_LOG("BitRefSelf", Sev::noisy) << "Parnets:" << std::endl
1089 << parent_ents << std::endl;
1090 }
1091
1092 Range children_ents;
1093 CHKERR updateRangeByChildren(parent_ents, children_ents);
1094 CHKERR filterEntitiesByRefLevel(child_bit, BitRefLevel().set(),
1095 children_ents, verb);
1096
1097 CHKERR moab.add_entities(meshset, children_ents);
1098
1099 if (verb >= VERY_VERBOSE) {
1101 MOFEM_LOG("BitRefSelf", Sev::noisy) << "Children: " << std::endl
1102 << children_ents << std::endl;
1103 }
1104 }
1106}
virtual const Field_multiIndex * get_fields() const =0
Get the fields object.
MoFEMErrorCode updateRangeByChildren(const Range &parent, Range &child, MoFEMTypes bh=MF_ZERO)
Update range by childresn.

◆ updateFieldMeshsetByEntitiesChildren() [2/2]

MoFEMErrorCode MoFEM::BitRefManager::updateFieldMeshsetByEntitiesChildren ( const std::string  name,
const BitRefLevel child_bit,
int  verb = 0 
)

update field meshset by child entities

Definition at line 1108 of file BitRefManager.cpp.

1109 {
1110 MoFEM::Interface &m_field = cOre;
1111 moab::Interface &moab = m_field.get_moab();
1113
1114 EntityHandle field_meshset = m_field.get_field_structure(name)->getMeshset();
1115
1116 Range parent_ents;
1117 CHKERR moab.get_entities_by_handle(field_meshset, parent_ents, true);
1118
1119 if (verb >= VERBOSE) {
1121 MOFEM_LOG("BitRefSelf", Sev::noisy) << "Parnets:" << endl
1122 << parent_ents << std::endl;
1123 }
1124
1125 Range children_ents;
1126 CHKERR updateRangeByChildren(parent_ents, children_ents);
1127 CHKERR filterEntitiesByRefLevel(child_bit, BitRefLevel().set(), children_ents,
1128 verb);
1129
1130 CHKERR moab.add_entities(field_meshset, children_ents);
1131
1132 if (verb >= VERBOSE) {
1134 MOFEM_LOG("BitRefSelf", Sev::noisy) << "Children: " << endl
1135 << children_ents << std::endl;
1136 }
1137
1139}
virtual Field * get_field_structure(const std::string &name)=0
get field structure
EntityHandle getMeshset() const
Get field meshset.

◆ updateFiniteElementMeshsetByEntitiesChildren()

MoFEMErrorCode MoFEM::BitRefManager::updateFiniteElementMeshsetByEntitiesChildren ( const std::string  name,
const BitRefLevel child_bit,
const EntityType  fe_ent_type,
int  verb = 0 
)

update finite element meshset by child entities

Definition at line 1141 of file BitRefManager.cpp.

1143 {
1144 MoFEM::Interface &m_field = cOre;
1146 EntityHandle meshset = m_field.get_finite_element_meshset(name);
1147 CHKERR updateMeshsetByEntitiesChildren(meshset, child_bit, meshset,
1148 fe_ent_type, false, verb);
1150}
MoFEMErrorCode updateMeshsetByEntitiesChildren(const EntityHandle parent, const BitRefLevel &parent_bit, const BitRefLevel &parent_mask, const BitRefLevel &child_bit, const BitRefLevel &child_mask, const EntityHandle child, EntityType child_type, const bool recursive=false, int verb=0)
Get child entities form meshset containing parent entities.
virtual EntityHandle get_finite_element_meshset(const std::string &name) const =0

◆ updateMeshsetByEntitiesChildren() [1/2]

MoFEMErrorCode MoFEM::BitRefManager::updateMeshsetByEntitiesChildren ( const EntityHandle  parent,
const BitRefLevel child_bit,
const EntityHandle  child,
EntityType  child_type,
const bool  recursive = false,
int  verb = 0 
)

Get child entities form meshset containing parent entities.

Search for refined entities of given type whose parent are entities in the parent meshset. It can be used for example to transfer information about boundary conditions to refined mesh or split mesh by interface elements. It is used by function refineMeshset, to update MESHSET finite elements.

Parameters
parentmeshset
parent_bitrefinement level
maskof parent bit ref level
child_bitrefinement level
maskof child bit ref level
childmeshset where child entities are stored (if the child meshset is set to be the parent meshset, the parent would be updated with the refined entities)
child_typemeshset is update only by entities of specified type. if type is set to MBMAXTYPE all types are updated.
recursiveif true parent meshset is searched recursively

Definition at line 1061 of file BitRefManager.cpp.

1064 {
1067 parent, BitRefLevel().set(), BitRefLevel().set(), child_bit, child_bit,
1068 child, child_type, recursive, verb);
1070}

◆ updateMeshsetByEntitiesChildren() [2/2]

MoFEMErrorCode MoFEM::BitRefManager::updateMeshsetByEntitiesChildren ( const EntityHandle  parent,
const BitRefLevel parent_bit,
const BitRefLevel parent_mask,
const BitRefLevel child_bit,
const BitRefLevel child_mask,
const EntityHandle  child,
EntityType  child_type,
const bool  recursive = false,
int  verb = 0 
)

Get child entities form meshset containing parent entities.

Search for refined entities of given type whose parent are entities in the parent meshset. It can be used for example to transfer information about boundary conditions to refined mesh or split mesh by interface elements. It is used by function refineMeshset, to update MESHSET finite elements.

Parameters
parentmeshset
parent_bitrefinement level
maskof parent bit ref level
child_bitrefinement level
maskof child bit ref level
childmeshset where child entities are stored (if the child meshset is set to be the parent meshset, the parent would be updated with the refined entities)
child_typemeshset is update only by entities of specified type. if type is set to MBMAXTYPE all types are updated.
recursiveif true parent meshset is searched recursively
Examples
hanging_node_approx.cpp.

Definition at line 1036 of file BitRefManager.cpp.

1040 {
1041 MoFEM::Interface &m_field = cOre;
1042 moab::Interface &moab = m_field.get_moab();
1044 Range parent_ents;
1045 CHKERR moab.get_entities_by_handle(parent, parent_ents, recursive);
1046 CHKERR filterEntitiesByRefLevel(parent_bit, parent_mask, parent_ents, verb);
1047 if (verb >= VERY_VERBOSE) {
1049 MOFEM_LOG("BitRefSelf", Sev::noisy) << "Parnets: " << parent;
1050 }
1051 Range children_ents;
1052 CHKERR updateRangeByChildren(parent_ents, children_ents);
1053 if (child_type < MBMAXTYPE)
1054 children_ents = children_ents.subset_by_type(child_type);
1055 CHKERR filterEntitiesByRefLevel(child_bit, BitRefLevel().set(), children_ents,
1056 verb);
1057 CHKERR moab.add_entities(child, children_ents);
1059}

◆ updateRange()

DEPRECATED MoFEMErrorCode MoFEM::BitRefManager::updateRange ( const Range &  parent,
Range &  child,
MoFEMTypes  bh = MF_ZERO 
)
Deprecated:
use updateRangeByChildren

Definition at line 585 of file BitRefManager.hpp.

585 {
586 return updateRangeByChildren(parent, child, bh);
587 }

◆ updateRangeByChildren()

MoFEMErrorCode MoFEM::BitRefManager::updateRangeByChildren ( const Range &  parent,
Range &  child,
MoFEMTypes  bh = MF_ZERO 
)

Update range by childresn.

FIXME: NOT TESTED

Parameters
parentparent range
childchildren range
Returns
error code
Examples
hanging_node_approx.cpp.

Definition at line 1152 of file BitRefManager.cpp.

1153 {
1154 MoFEM::Interface &m_field = cOre;
1155 auto ref_ents_ptr = m_field.get_ref_ents();
1157 auto &ref_ents =
1158 const_cast<RefEntity_multiIndex *>(ref_ents_ptr)->get<Ent_Ent_mi_tag>();
1159 std::vector<EntityHandle> child_ents_vec;
1160 child_ents_vec.reserve(ref_ents.size());
1161 for (Range::const_pair_iterator pit = parent_ents.const_pair_begin();
1162 pit != parent_ents.const_pair_end(); pit++) {
1163 const auto f = pit->first;
1164 auto it = ref_ents.lower_bound(f);
1165 if (it != ref_ents.end()) {
1166 const auto s = pit->second;
1167 auto hi_it = ref_ents.upper_bound(s);
1168 if (bh == MF_EXIST) {
1169 if (std::distance(it, hi_it) != (s - f) + 1) {
1170 SETERRQ2(
1171 PETSC_COMM_SELF, MOFEM_NOT_FOUND,
1172 "Number of entities and entities parents is different %d != %d ",
1173 std::distance(it, hi_it), (s - f) + 1);
1174 }
1175 }
1176 for (; it != hi_it; it++) {
1177#ifndef NDEBUG
1178 if (it->get()->getEntType() == MBENTITYSET) {
1179 SETERRQ(PETSC_COMM_SELF, MOFEM_IMPOSIBLE_CASE,
1180 "This should not happen; Entity should not have part of the "
1181 "meshset. It has no children.");
1182 }
1183#endif
1184 child_ents_vec.emplace_back((*it)->getEnt());
1185 }
1186 } else if (bh == MF_EXIST) {
1187 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_FOUND, "Entities not found");
1188 }
1189 }
1190 child_ents.insert_list(child_ents_vec.begin(), child_ents_vec.end());
1192}
@ MF_EXIST
Definition: definitions.h:113
@ MOFEM_NOT_FOUND
Definition: definitions.h:46
@ MOFEM_IMPOSIBLE_CASE
Definition: definitions.h:48

◆ updateRangeByParent()

MoFEMErrorCode MoFEM::BitRefManager::updateRangeByParent ( const Range &  parent,
Range &  child,
MoFEMTypes  bh = MF_ZERO 
)

Update range by prents.

FIXME: NOT TESTED

Parameters
childparent range
parentchildren range
Returns
error code

Definition at line 1194 of file BitRefManager.cpp.

1196 {
1197 MoFEM::Interface &m_field = cOre;
1198 auto ref_ents_ptr = m_field.get_ref_ents();
1200 auto &ref_ents =
1201 const_cast<RefEntity_multiIndex *>(ref_ents_ptr)->get<Ent_mi_tag>();
1202 std::vector<EntityHandle> parent_ents_vec;
1203 parent_ents_vec.reserve(ref_ents.size());
1204 for (Range::const_pair_iterator pit = child_ents.const_pair_begin();
1205 pit != child_ents.const_pair_end(); pit++) {
1206 const auto f = pit->first;
1207 auto it = ref_ents.lower_bound(f);
1208 if (it != ref_ents.end()) {
1209 const auto s = pit->second;
1210 auto hi_it = ref_ents.upper_bound(s);
1211 if (bh == MF_EXIST) {
1212 if (std::distance(it, hi_it) != (s - f) + 1) {
1213 SETERRQ2(
1214 PETSC_COMM_SELF, MOFEM_NOT_FOUND,
1215 "Number of entities and enties parents is diffrent %d != %d ",
1216 std::distance(it, hi_it), (s - f) + 1);
1217 }
1218 }
1219 for (; it != hi_it; it++) {
1220#ifndef NDEBUG
1221 if (it->get()->getEntType() == MBENTITYSET) {
1222 SETERRQ(PETSC_COMM_SELF, MOFEM_IMPOSIBLE_CASE,
1223 "This should not happen; Entity should not have part of the "
1224 "meshset. It has no children.");
1225 }
1226#endif
1227 auto parent = (*it)->getParentEnt();
1228 if (parent)
1229 parent_ents_vec.emplace_back(parent);
1230 }
1231 } else if (bh == MF_EXIST) {
1232 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_FOUND, "Entities not found");
1233 }
1234 }
1235 parent_ents.insert_list(parent_ents_vec.begin(), parent_ents_vec.end());
1237}