v0.15.0
Loading...
Searching...
No Matches
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
 
MoFEMErrorCode MoFEM::BitRefManager::setElementsBitRefLevel (const Range &ents, const BitRefLevel bit=BitRefLevel(), int verb=QUIET) const
 add entities to database and set bit ref level
 
MoFEMErrorCode MoFEM::BitRefManager::setEntitiesBitRefLevel (const Range &ents, const BitRefLevel bit=BitRefLevel(), int verb=QUIET) const
 add entities to database and set bit ref level
 
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
 
MoFEMErrorCode MoFEM::BitRefManager::setNthBitRefLevel (const int n, const bool b, int verb=QUIET) const
 Set nth bit ref level to all entities in database.
 
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
 
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
 
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.
 
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.
 
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.
 
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.
 
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.
 
MoFEMErrorCode MoFEM::BitRefManager::lambdaBitRefLevel (boost::function< void(EntityHandle ent, BitRefLevel &bit)> fun) const
 Process bit ref level by lambda function.
 
MoFEMErrorCode MoFEM::BitRefManager::lambdaBitRefLevel (const Range &ents, boost::function< void(EntityHandle ent, BitRefLevel &bit)> fun) const
 Process bit ref level by lambda function.
 
MoFEMErrorCode MoFEM::BitRefManager::addBitRefLevelByDim (const EntityHandle meshset, const int dim, const BitRefLevel bit, int verb=QUIET) const
 add bit ref level by dimension
 
MoFEMErrorCode MoFEM::BitRefManager::setNthBitRefLevel (const Range &ents, const int n, const bool b, int verb=QUIET) const
 Set nth bit ref level.
 

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
 
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
 
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
 
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
 
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
 
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
 

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.
 
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.
 
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.
 
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.
 

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.
 
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.
 
MoFEMErrorCode MoFEM::BitRefManager::updateFieldMeshsetByEntitiesChildren (const BitRefLevel &child_bit, int verb=0)
 update fields meshesets by child entities
 
MoFEMErrorCode MoFEM::BitRefManager::updateFieldMeshsetByEntitiesChildren (const std::string name, const BitRefLevel &child_bit, int verb=0)
 update field meshset by child entities
 
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
 
MoFEMErrorCode MoFEM::BitRefManager::updateRangeByChildren (const Range &parent, Range &child, MoFEMTypes bh=MF_ZERO)
 Update range by childrens.
 
MoFEMErrorCode MoFEM::BitRefManager::updateRangeByParent (const Range &child_ents, Range &parent_ents, MoFEMTypes bh=MF_ZERO)
 Update range by parents.
 
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

#include <src/interfaces/BitRefManager.hpp>

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 555 of file BitRefManager.cpp.

557 {
558 return lambdaBitRefLevel(
559 ents, [&](EntityHandle ent, BitRefLevel &ent_bit) { ent_bit |= bit; });
560}
MoFEMErrorCode lambdaBitRefLevel(boost::function< void(EntityHandle ent, BitRefLevel &bit)> fun) const
Process bit ref level by lambda function.
auto bit
set bit

◆ addBitRefLevelByDim()

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

#include <src/interfaces/BitRefManager.hpp>

add bit ref level by dimension

Parameters
meshset
dimdimension of entities
bitadded bit
verbverbosity level
Returns
MoFEMErrorCode

Definition at line 562 of file BitRefManager.cpp.

565 {
566 MoFEM::Interface &m_field = cOre;
567 moab::Interface &moab = m_field.get_moab();
568 Range ents, adj;
570 CHKERR moab.get_entities_by_dimension(meshset, dim, ents, true);
571 for (int dd = dim - 1; dd >= 0; dd--)
572 CHKERR moab.get_adjacencies(ents, dd, false, adj, moab::Interface::UNION);
573 ents.merge(adj);
574 if (verb == VERY_NOISY)
575 MOFEM_LOG("BitRefSelf", Sev::noisy) << "Add add bit ref level by dim";
576 CHKERR addBitRefLevel(ents, bit, verb);
578}
@ VERY_NOISY
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define CHKERR
Inline error check.
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.
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
virtual moab::Interface & get_moab()=0
Deprecated interface functions.

◆ addToDatabaseBitRefLevelByDim()

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

#include <src/interfaces/BitRefManager.hpp>

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 464 of file BitRefManager.cpp.

466 {
468 Range ents;
469 CHKERR getEntitiesByDimAndRefLevel(bit, mask, dim, ents);
470 CHKERR setBitRefLevel(ents, BitRefLevel(), false, verb);
472}
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
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition Types.hpp:40

◆ addToDatabaseBitRefLevelByType()

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

#include <src/interfaces/BitRefManager.hpp>

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 454 of file BitRefManager.cpp.

456 {
458 Range ents;
459 CHKERR getEntitiesByTypeAndRefLevel(bit, mask, type, ents);
460 CHKERR setBitRefLevel(ents, BitRefLevel(), false, verb);
462}
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

#include <src/interfaces/BitRefManager.hpp>

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 989 of file BitRefManager.cpp.

992 {
993 MoFEM::Interface &m_field = cOre;
994 moab::Interface &moab(m_field.get_moab());
996
997 if (verb > VERBOSE) {
999 MOFEM_LOG("BitRef", Sev::noisy) << "from: " << bit;
1000 }
1001
1002 CHKERR moab.get_adjacencies(from_entities, num_entities, to_dimension, false,
1003 adj_entities, operation_type);
1004 std::vector<BitRefLevel> bit_levels(adj_entities.size());
1005 CHKERR moab.tag_get_data(cOre.get_th_RefBitLevel(), adj_entities,
1006 &*bit_levels.begin());
1007 std::vector<BitRefLevel>::iterator b_it = bit_levels.begin();
1008 for (Range::iterator eit = adj_entities.begin(); eit != adj_entities.end();
1009 b_it++) {
1010
1011 if (verb > VERBOSE) {
1012 RefEntity adj_entity(m_field.get_basic_entity_data_ptr(), *eit);
1014 MOFEM_LOG("BitRef", Sev::noisy)
1015 << "to: " << adj_entity.getBitRefLevel() << " : " << adj_entity;
1016 }
1017
1018 if (!((*b_it) & bit).any()) {
1019 eit = adj_entities.erase(eit);
1020 } else {
1021 eit++;
1022 }
1023 }
1024 if (b_it != bit_levels.end()) {
1025 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "Data inconsistency");
1026 }
1028}
@ VERBOSE
@ MOFEM_DATA_INCONSISTENCY
Definition definitions.h:31
#define MOFEM_LOG_FUNCTION()
Set scope.
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:198

◆ 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

#include <src/interfaces/BitRefManager.hpp>

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 978 of file BitRefManager.cpp.

981 {
983 BitRefLevel bit = problem_ptr->getBitRefLevel();
984 CHKERR getAdjacencies(bit, from_entities, num_entities, to_dimension,
985 adj_entities, operation_type);
987}
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

#include <src/interfaces/BitRefManager.hpp>

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 952 of file BitRefManager.cpp.

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

◆ getAdjacenciesEquality()

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

#include <src/interfaces/BitRefManager.hpp>

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 926 of file BitRefManager.cpp.

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

◆ getEntitiesByDimAndRefLevel() [1/2]

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

#include <src/interfaces/BitRefManager.hpp>

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
free_surface.cpp, and hanging_node_approx.cpp.

Definition at line 823 of file BitRefManager.cpp.

825 {
826 MoFEM::Interface &m_field = cOre;
827 moab::Interface &moab(m_field.get_moab());
829 Range ents;
830 CHKERR getEntitiesByDimAndRefLevel(bit, mask, dim, ents, verb);
831 CHKERR moab.add_entities(meshset, ents);
833}

◆ getEntitiesByDimAndRefLevel() [2/2]

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

#include <src/interfaces/BitRefManager.hpp>

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 835 of file BitRefManager.cpp.

837 {
838 MoFEM::Interface &m_field = cOre;
839 moab::Interface &moab(m_field.get_moab());
841 CHKERR moab.get_entities_by_dimension(0, dim, ents, false);
842 CHKERR filterEntitiesByRefLevel(bit, mask, ents, verb);
844}
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

#include <src/interfaces/BitRefManager.hpp>

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 846 of file BitRefManager.cpp.

849 {
850 MoFEM::Interface &m_field = cOre;
851 moab::Interface &moab(m_field.get_moab());
853 Range ents;
854 CHKERR getEntitiesByRefLevel(bit, mask, ents, verb);
855 CHKERR moab.add_entities(meshset, ents);
857}
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

#include <src/interfaces/BitRefManager.hpp>

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 859 of file BitRefManager.cpp.

862 {
863 MoFEM::Interface &m_field = cOre;
864 moab::Interface &moab(m_field.get_moab());
866 Range meshset_ents;
867 CHKERR moab.get_entities_by_type(0, MBENTITYSET, meshset_ents, false);
868 CHKERR moab.get_entities_by_handle(0, ents, false);
869 ents.merge(meshset_ents);
870 CHKERR filterEntitiesByRefLevel(bit, mask, ents, verb);
872}

◆ getEntitiesByTypeAndRefLevel() [1/2]

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

#include <src/interfaces/BitRefManager.hpp>

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 735 of file BitRefManager.cpp.

737 {
738 MoFEM::Interface &m_field = cOre;
739 moab::Interface &moab(m_field.get_moab());
741 Range ents;
742 CHKERR getEntitiesByTypeAndRefLevel(bit, mask, type, ents, verb);
743 CHKERR moab.add_entities(meshset, ents);
745}

◆ getEntitiesByTypeAndRefLevel() [2/2]

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

#include <src/interfaces/BitRefManager.hpp>

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 812 of file BitRefManager.cpp.

814 {
815 MoFEM::Interface &m_field = cOre;
816 moab::Interface &moab(m_field.get_moab());
818 CHKERR moab.get_entities_by_type(0, type, ents, false);
819 CHKERR filterEntitiesByRefLevel(bit, mask, ents, verb);
821}

◆ lambdaBitRefLevel() [1/2]

MoFEMErrorCode MoFEM::BitRefManager::lambdaBitRefLevel ( boost::function< void(EntityHandle ent, BitRefLevel &bit)> fun) const

#include <src/interfaces/BitRefManager.hpp>

Process bit ref level by lambda function.

Note
That not apply to type of MBENTITY. To avoid problems with problem meshsets.
Parameters
fun
Returns
MoFEMErrorCode

Definition at line 525 of file BitRefManager.cpp.

526 {
527 MoFEM::Interface &m_field = cOre;
529 auto get_ents = [&]() {
530 Range ents;
531 CHKERR m_field.get_moab().get_entities_by_handle(
532 m_field.get_moab().get_root_set(), ents, true);
533 ents = subtract(ents, ents.subset_by_type(MBENTITYSET));
534 return ents;
535 };
536 CHKERR lambdaBitRefLevel(get_ents(), fun);
538}
auto fun
Function to approximate.

◆ lambdaBitRefLevel() [2/2]

MoFEMErrorCode MoFEM::BitRefManager::lambdaBitRefLevel ( const Range & ents,
boost::function< void(EntityHandle ent, BitRefLevel &bit)> fun ) const

#include <src/interfaces/BitRefManager.hpp>

Process bit ref level by lambda function.

Parameters
fun
Returns
MoFEMErrorCode

Definition at line 540 of file BitRefManager.cpp.

542 {
543 MoFEM::Interface &m_field = cOre;
545 std::vector<const BitRefLevel *> ents_bits_vec;
546 CHKERR RefEntity::getBitRefLevel(m_field.get_moab(), ents, ents_bits_vec);
547 auto eit = ents.begin();
548 for (auto &it : ents_bits_vec) {
549 fun(*eit, const_cast<BitRefLevel &>(*it));
550 ++eit;
551 }
553};
const BitRefLevel & getBitRefLevel() const
Get entity ref bit refinement signature.

◆ setBitLevelToMeshset()

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

#include <src/interfaces/BitRefManager.hpp>

brief add meshset and set bit ref level

Parameters
EntityHandleMeshSet
BitRefLevelbitLevel

Definition at line 474 of file BitRefManager.cpp.

476 {
477 MoFEM::Interface &m_field = cOre;
478 auto ref_ents_ptr = m_field.get_ref_ents();
479 auto ref_fe_ptr = m_field.get_ref_finite_elements();
481 // Add ref entity
482 std::pair<RefEntity_multiIndex::iterator, bool> p_ent =
483 const_cast<RefEntity_multiIndex *>(ref_ents_ptr)
484 ->insert(boost::shared_ptr<RefEntity>(
485 new RefEntity(m_field.get_basic_entity_data_ptr(), meshset)));
486 *(const_cast<RefEntity *>(p_ent.first->get())->getBitRefLevelPtr()) |= bit;
487 // Add ref element
488 boost::shared_ptr<RefElement> fe_ptr =
489 boost::shared_ptr<RefElement>(new RefElement_MESHSET(*p_ent.first));
490 std::pair<RefElement_multiIndex::iterator, bool> p_fe =
491 const_cast<RefElement_multiIndex *>(ref_fe_ptr)->insert(fe_ptr);
492
494 MOFEM_LOG("BitRefSelf", Sev::noisy)
495 << "Add meshset as ref_ent " << **p_fe.first;
496
498}
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

#include <src/interfaces/BitRefManager.hpp>

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 281 of file BitRefManager.cpp.

284 {
285 MoFEM::Interface &m_field = cOre;
286 auto ref_ents_ptr = m_field.get_ref_ents();
287 auto ref_fe_ptr = m_field.get_ref_finite_elements();
289
291 MOFEM_LOG_C("BitRefSelf", Sev::noisy, "Number of entities to add %d",
292 ents.size());
293
294 CHKERR setElementsBitRefLevel(ents, bit, verb);
295
296 if (!ents.empty()) {
297 for (int d = 3; d >= 1; --d) {
298 Range dim_ents;
299 if (only_tets && d == 3) {
300 dim_ents = ents.subset_by_type(MBTET);
301 } else {
302 dim_ents = ents.subset_by_dimension(d);
303 }
304
306 MOFEM_LOG_C("BitRefSelf", Sev::noisy,
307 " Number of dim %d entities to add %d", d,
308 dim_ents.size());
309
310 if (!dim_ents.empty()) {
311 for (int dd = 0; dd < d; ++dd) {
312 Range adj_ents;
313
314 if (dd == 0) {
315 rval = m_field.get_moab().get_connectivity(ents, adj_ents, true);
316
317 } else {
318 if (adj_ents_ptr) {
319 if (dd == 1) {
320 adj_ents = adj_ents_ptr->subset_by_dimension(MBEDGE);
321 } else if (dd == 2) {
322 adj_ents = adj_ents_ptr->subset_by_dimension(MBTRI);
323 }
324 } else {
325 rval = m_field.get_moab().get_adjacencies(
326 dim_ents, dd, true, adj_ents, moab::Interface::UNION);
327 }
328 }
329
330
332 MOFEM_LOG_C("BitRefSelf", Sev::noisy,
333 " Number of dim %d adj entities for dim %d to add %d",
334 d, dd, adj_ents.size());
335
336 if (rval == MB_MULTIPLE_ENTITIES_FOUND) {
337 auto log_message = [&](const auto sev) {
340 MOFEM_LOG("BitRefSelf", sev)
341 << "When get adjacencies moab return MB_MULTIPLE_ENTITIES_ "
342 "FOUND for dim = "
343 << dd << " and dim of entities " << d;
344 MOFEM_LOG_CHANNEL("BitRefSelf"); // reset channel
345 };
346
347 if (verb <= QUIET)
348 log_message(Sev::noisy);
349 else
350 log_message(Sev::warning);
351
352 rval = MB_SUCCESS;
353 }
355 for (Range::pair_iterator pit = adj_ents.pair_begin();
356 pit != adj_ents.pair_end(); ++pit) {
357 Range seed_ents_range;
358 // get first and last element of range
359 EntityHandle f = pit->first;
360 EntityHandle s = pit->second;
361 CHKERR SetBitRefLevelTool(m_field, bit, ref_ents_ptr, ref_fe_ptr)
362 .findEntsToAdd(f, s, seed_ents_range);
363 if (!seed_ents_range.empty())
364 CHKERR SetBitRefLevelTool(m_field, bit, ref_ents_ptr,
365 ref_fe_ptr)
366 .addEntsToDatabase(seed_ents_range);
367 }
368 }
369 }
370 }
371 }
372
374}
#define MOFEM_LOG_C(channel, severity, format,...)
@ QUIET
#define MOAB_THROW(err)
Check error code of MoAB function and throw MoFEM exception.
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.
#define MOFEM_LOG_ATTRIBUTES(channel, bit)
Add attributes to channel.
static MoFEMErrorCodeGeneric< moab::ErrorCode > rval

◆ setBitRefLevelByDim()

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

#include <src/interfaces/BitRefManager.hpp>

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 500 of file BitRefManager.cpp.

503 {
504 MoFEM::Interface &m_field = cOre;
506 Range ents;
507 CHKERR m_field.get_moab().get_entities_by_dimension(meshset, dim, ents,
508 false);
509 CHKERR setBitRefLevel(ents, bit, false, verb);
511}

◆ setBitRefLevelByType()

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

#include <src/interfaces/BitRefManager.hpp>

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 513 of file BitRefManager.cpp.

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

◆ setElementsBitRefLevel()

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

#include <src/interfaces/BitRefManager.hpp>

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 376 of file BitRefManager.cpp.

378 {
379 MoFEM::Interface &m_field = cOre;
380 auto ref_ents_ptr = m_field.get_ref_ents();
381 auto ref_fe_ptr = m_field.get_ref_finite_elements();
382
384
385 for (Range::const_pair_iterator pit = ents.pair_begin();
386 pit != ents.pair_end(); pit++) {
387 // get first and last element of range
388 EntityHandle f = pit->first;
389 EntityHandle s = pit->second;
390 Range seed_ents_range; // entities seeded not in database
391 // find ents to add
392 CHKERR SetBitRefLevelTool(m_field, bit, ref_ents_ptr, ref_fe_ptr)
393 .findEntsToAdd(f, s, seed_ents_range);
394 // add elements
395 if (!seed_ents_range.empty())
396 CHKERR SetBitRefLevelTool(m_field, bit, ref_ents_ptr, ref_fe_ptr)
397 .addEntsToDatabase(seed_ents_range);
398
399 Range seed_fe_range;
400 CHKERR SetBitRefLevelTool(m_field, bit, ref_ents_ptr, ref_fe_ptr)
401 .findElementsToAdd(f, s, seed_fe_range);
402 if (!seed_fe_range.empty()) {
403 CHKERR SetBitRefLevelTool(m_field, bit, ref_ents_ptr, ref_fe_ptr)
404 .addElementsToDatabase(seed_fe_range);
405 }
406 }
407
409 MOFEM_LOG("BitRefSelf", Sev::noisy)
410 << "Number of entities in database " << ref_ents_ptr->size();
411 MOFEM_LOG("BitRefSelf", Sev::noisy)
412 << "Number of finite element entities in database " << ref_fe_ptr->size();
413
415}

◆ setEntitiesBitRefLevel()

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

#include <src/interfaces/BitRefManager.hpp>

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 417 of file BitRefManager.cpp.

419 {
420 MoFEM::Interface &m_field = cOre;
421 auto ref_ents_ptr = m_field.get_ref_ents();
422 auto ref_fe_ptr = m_field.get_ref_finite_elements();
424
425 for (Range::const_pair_iterator pit = ents.pair_begin();
426 pit != ents.pair_end(); pit++) {
427 // get first and last element of range
428 EntityHandle f = pit->first;
429 EntityHandle s = pit->second;
430 Range seed_ents_range; // entities seeded not in database
431 // find ents to add
432 CHKERR SetBitRefLevelTool(m_field, bit, ref_ents_ptr, ref_fe_ptr)
433 .findEntsToAdd(f, s, seed_ents_range);
434 // add elements
435 if (!seed_ents_range.empty())
436 CHKERR SetBitRefLevelTool(m_field, bit, ref_ents_ptr, ref_fe_ptr)
437 .addEntsToDatabase(seed_ents_range);
438 }
440}

◆ setFieldEntitiesBitRefLevel()

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

#include <src/interfaces/BitRefManager.hpp>

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

Parameters
field_name
bit
verb
Returns
MoFEMErrorCode
Examples
dm_partitioned_no_field.cpp.

Definition at line 442 of file BitRefManager.cpp.

443 {
444 MoFEM::Interface &m_field = cOre;
446 EntityHandle field_meshset = m_field.get_field_meshset(field_name);
447 Range field_ents;
448 CHKERR m_field.get_moab().get_entities_by_handle(field_meshset, field_ents,
449 true);
450 CHKERR setEntitiesBitRefLevel(field_ents, bit, verb);
452}
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

#include <src/interfaces/BitRefManager.hpp>

Set nth bit ref level to all entities in database.

Parameters
nnth bit
bvalue to set
Returns
error code

Definition at line 588 of file BitRefManager.cpp.

589 {
590 if (verb == VERY_NOISY)
591 MOFEM_LOG("BitRefSelf", Sev::noisy) << "Set bit to all entities";
592 return lambdaBitRefLevel(
593 [&](EntityHandle ent, BitRefLevel &ent_bit) { ent_bit[n] = b; });
594}
const double n
refractive index of diffusive medium

◆ setNthBitRefLevel() [2/2]

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

#include <src/interfaces/BitRefManager.hpp>

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 580 of file BitRefManager.cpp.

581 {
582 if (verb == VERY_NOISY)
583 MOFEM_LOG("BitRefSelf", Sev::noisy) << "Set bit to " << ents;
584 return lambdaBitRefLevel(
585 ents, [&](EntityHandle ent, BitRefLevel &ent_bit) { ent_bit[n] = b; });
586}

◆ shiftLeftBitRef()

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

#include <src/interfaces/BitRefManager.hpp>

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

Note
Not implemented

Definition at line 596 of file BitRefManager.cpp.

598 {
600 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
602}
@ MOFEM_NOT_IMPLEMENTED
Definition definitions.h:32

◆ shiftRightBitRef()

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

#include <src/interfaces/BitRefManager.hpp>

right shift bit ref level

Definition at line 604 of file BitRefManager.cpp.

606 {
607 MoFEM::Interface &m_field = cOre;
608 auto ref_ents_ptr = m_field.get_ref_ents();
610 RefEntity_change_right_shift right_shift(1, mask);
611 for (int ii = 0; ii < shift; ii++) {
612 // delete bits on the right which are shifted to zero
613 BitRefLevel delete_bits = BitRefLevel().set(0) & mask;
614 if (delete_bits.any()) {
615 CHKERR m_field.delete_ents_by_bit_ref(delete_bits, delete_bits, true,
616 verb, mf);
617 }
618 for (RefEntity_multiIndex::iterator ent_it = ref_ents_ptr->begin();
619 ent_it != ref_ents_ptr->end(); ent_it++) {
620 if (verb >= NOISY) {
622 MOFEM_LOG("BitRefSelf", Sev::noisy)
623 << (*ent_it)->getBitRefLevel() << " : ";
624 }
625 right_shift(const_cast<boost::shared_ptr<RefEntity> &>(*ent_it));
626 if (verb == VERY_NOISY) {
628 MOFEM_LOG("BitRefSelf", Sev::noisy) << (*ent_it)->getBitRefLevel();
629 }
630 }
631 }
633}
@ NOISY
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 )

#include <src/interfaces/BitRefManager.hpp>

update fields meshesets by child entities

Note
This calls updateMeshsetByEntitiesChildren for all entity types.

Definition at line 1066 of file BitRefManager.cpp.

1067 {
1068 MoFEM::Interface &m_field = cOre;
1069 moab::Interface &moab = m_field.get_moab();
1070 auto fields_ptr = m_field.get_fields();
1072
1073 for (auto &fit : (*fields_ptr)) {
1074
1075 const EntityHandle meshset = fit->getMeshset();
1076 Range parent_ents;
1077 CHKERR moab.get_entities_by_handle(meshset, parent_ents, true);
1078
1079 if (verb >= VERY_VERBOSE) {
1081 MOFEM_LOG("BitRefSelf", Sev::noisy) << "Parnets:" << std::endl
1082 << parent_ents << std::endl;
1083 }
1084
1085 Range children_ents;
1086 CHKERR updateRangeByChildren(parent_ents, children_ents);
1087 CHKERR filterEntitiesByRefLevel(child_bit, BitRefLevel().set(),
1088 children_ents, verb);
1089
1090 CHKERR moab.add_entities(meshset, children_ents);
1091
1092 if (verb >= VERY_VERBOSE) {
1094 MOFEM_LOG("BitRefSelf", Sev::noisy) << "Children: " << std::endl
1095 << children_ents << std::endl;
1096 }
1097 }
1099}
@ VERY_VERBOSE
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 childrens.

◆ updateFieldMeshsetByEntitiesChildren() [2/2]

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

#include <src/interfaces/BitRefManager.hpp>

update field meshset by child entities

Definition at line 1101 of file BitRefManager.cpp.

1102 {
1103 MoFEM::Interface &m_field = cOre;
1104 moab::Interface &moab = m_field.get_moab();
1106
1107 EntityHandle field_meshset = m_field.get_field_structure(name)->getMeshset();
1108
1109 Range parent_ents;
1110 CHKERR moab.get_entities_by_handle(field_meshset, parent_ents, true);
1111
1112 if (verb >= VERBOSE) {
1114 MOFEM_LOG("BitRefSelf", Sev::noisy) << "Parnets:" << endl
1115 << parent_ents << std::endl;
1116 }
1117
1118 Range children_ents;
1119 CHKERR updateRangeByChildren(parent_ents, children_ents);
1120 CHKERR filterEntitiesByRefLevel(child_bit, BitRefLevel().set(), children_ents,
1121 verb);
1122
1123 CHKERR moab.add_entities(field_meshset, children_ents);
1124
1125 if (verb >= VERBOSE) {
1127 MOFEM_LOG("BitRefSelf", Sev::noisy) << "Children: " << endl
1128 << children_ents << std::endl;
1129 }
1130
1132}
virtual const Field * get_field_structure(const std::string &name, enum MoFEMTypes bh=MF_EXIST) const =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 )

#include <src/interfaces/BitRefManager.hpp>

update finite element meshset by child entities

Definition at line 1134 of file BitRefManager.cpp.

1136 {
1137 MoFEM::Interface &m_field = cOre;
1139 EntityHandle meshset = m_field.get_finite_element_meshset(name);
1140 CHKERR updateMeshsetByEntitiesChildren(meshset, child_bit, meshset,
1141 fe_ent_type, false, verb);
1143}
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 )

#include <src/interfaces/BitRefManager.hpp>

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 1055 of file BitRefManager.cpp.

1058 {
1061 parent, BitRefLevel().set(), BitRefLevel().set(), child_bit, child_bit,
1062 child, child_type, recursive, verb);
1064}

◆ 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 )

#include <src/interfaces/BitRefManager.hpp>

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 1030 of file BitRefManager.cpp.

1034 {
1035 MoFEM::Interface &m_field = cOre;
1036 moab::Interface &moab = m_field.get_moab();
1038 Range parent_ents;
1039 CHKERR moab.get_entities_by_handle(parent, parent_ents, recursive);
1040 CHKERR filterEntitiesByRefLevel(parent_bit, parent_mask, parent_ents, verb);
1041 if (verb >= VERY_VERBOSE) {
1043 MOFEM_LOG("BitRefSelf", Sev::noisy) << "Parents: " << parent;
1044 }
1045 Range children_ents;
1046 CHKERR updateRangeByChildren(parent_ents, children_ents);
1047 if (child_type < MBMAXTYPE)
1048 children_ents = children_ents.subset_by_type(child_type);
1049 CHKERR filterEntitiesByRefLevel(child_bit, BitRefLevel().set(), children_ents,
1050 verb);
1051 CHKERR moab.add_entities(child, children_ents);
1053}

◆ updateRange()

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

#include <src/interfaces/BitRefManager.hpp>

Deprecated
use updateRangeByChildren

Definition at line 598 of file BitRefManager.hpp.

598 {
599 return updateRangeByChildren(parent, child, bh);
600 }

◆ updateRangeByChildren()

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

#include <src/interfaces/BitRefManager.hpp>

Update range by childrens.

FIXME: NOT TESTED

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

Definition at line 1145 of file BitRefManager.cpp.

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

◆ updateRangeByParent()

MoFEMErrorCode MoFEM::BitRefManager::updateRangeByParent ( const Range & child_ents,
Range & parent_ents,
MoFEMTypes bh = MF_ZERO )

#include <src/interfaces/BitRefManager.hpp>

Update range by parents.

FIXME: NOT TESTED

Parameters
childparent range
parentchildren range
Returns
error code

Definition at line 1188 of file BitRefManager.cpp.

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