v0.13.2
Loading...
Searching...
No Matches
Classes | Public Types | Public Member Functions | Static Public Member Functions | Private Attributes | List of all members
MoFEM::BcManager Struct Reference

Simple interface for fast problem set-up. More...

#include <src/interfaces/BcManager.hpp>

Inheritance diagram for MoFEM::BcManager:
[legend]
Collaboration diagram for MoFEM::BcManager:
[legend]

Classes

struct  BCs
 Data structure storing bc markers and atributes. More...
 

Public Types

using BcMapByBlockName = std::map< string, boost::shared_ptr< BCs > >
 
using BcMarkerPtr = boost::shared_ptr< std::vector< char unsigned > >
 

Public Member Functions

MoFEMErrorCode query_interface (boost::typeindex::type_index type_index, UnknownInterface **iface) const
 
 BcManager (const MoFEM::Core &core)
 
virtual ~BcManager ()=default
 
MoFEMErrorCode getOptions ()
 get options More...
 
MoFEMErrorCode removeBlockDOFsOnEntities (const std::string problem_name, const std::string block_name, const std::string field_name, int lo, int hi, bool get_low_dim_ents=true, bool is_distributed_mesh=true)
 Remove DOFs from problem. More...
 
MoFEMErrorCode pushMarkDOFsOnEntities (const std::string problem_name, const std::string block_name, const std::string field_name, int lo, int hi, bool get_low_dim_ents=true)
 Mark block DOFs. More...
 
template<typename T >
MoFEMErrorCode removeBlockDOFsOnEntities (const std::string problem_name, const std::string field_name, bool get_low_dim_ents=true, bool block_name_field_prefix=false, bool is_distributed_mesh=true)
 Mark block DOFs. More...
 
template<typename T >
MoFEMErrorCode pushMarkDOFsOnEntities (const std::string problem_name, const std::string field_name, bool get_low_dim_ents=true, bool block_name_field_prefix=false)
 Mark block DOFs. More...
 
template<typename T >
MoFEMErrorCode removeBlockDOFsOnEntities (const std::string problem_name, const std::string block_name, const std::string field_name, bool get_low_dim_ents=true, bool is_distributed_mesh=true)
 Mark block DOFs. More...
 
template<typename T >
MoFEMErrorCode pushMarkDOFsOnEntities (const std::string problem_name, const std::string block_name, const std::string field_name, bool get_low_dim_ents=true)
 Mark block DOFs. More...
 
boost::shared_ptr< BCspopMarkDOFsOnEntities (const std::string block_name)
 Get bc data and remove element. More...
 
auto getBcStructure (const std::string bc_id)
 Get the bc structure object. More...
 
BcMapByBlockNamegetBcMapByBlockName ()
 Get the bc map. More...
 
BcMarkerPtr getMergedBlocksMarker (std::vector< std::regex > bc_regex_vec)
 Get the Merged Boundary Marker object. More...
 
auto getMergedBlocksMarker (std::vector< string > bc_names)
 Get the Merged Boundary Marker object. More...
 
BcMarkerPtr getMergedBlocksMarker (const std::vector< BcMarkerPtr > &boundary_markers_ptr_vec)
 Get the Merged Blocks Marker object. More...
 
auto checkBlock (const std::pair< string, boost::shared_ptr< BCs > > &bc, std::regex reg)
 check if given boundary condition name is in the map bc element More...
 
auto checkBlock (const std::pair< std::string, boost::shared_ptr< BCs > > &bc, std::string name)
 check if given boundary condition name is in the map bc element More...
 
SmartPetscObj< IS > getBlockIS (const std::string block_prefix, const std::string block_name, const std::string field_name, const std::string problem_name, int lo, int hi, SmartPetscObj< IS > is_expand=SmartPetscObj< IS >())
 Get block IS. More...
 
SmartPetscObj< IS > getBlockIS (const std::string problem_name, const std::string block_name, const std::string field_name, int lo, int hi, SmartPetscObj< IS > is_expand=SmartPetscObj< IS >())
 Get block IS. More...
 
template<>
MoFEMErrorCode removeBlockDOFsOnEntities (const std::string problem_name, const std::string field_name, bool get_low_dim_ents, bool block_name_field_prefix, bool is_distributed_mesh)
 
template<>
MoFEMErrorCode removeBlockDOFsOnEntities (const std::string problem_name, const std::string field_name, bool get_low_dim_ents, bool block_name_field_prefix, bool is_distributed_mesh)
 
template<>
MoFEMErrorCode removeBlockDOFsOnEntities (const std::string problem_name, const std::string field_name, bool get_low_dim_ents, bool block_name_field_prefix, bool is_distributed_mesh)
 
template<>
MoFEMErrorCode removeBlockDOFsOnEntities (const std::string problem_name, const std::string field_name, bool get_low_dim_ents, bool block_name_field_prefix, bool is_distributed_mesh)
 
template<>
MoFEMErrorCode removeBlockDOFsOnEntities (const std::string problem_name, const std::string block_name, const std::string field_name, bool get_low_dim_ents, bool is_distributed_mesh)
 
template<>
MoFEMErrorCode pushMarkDOFsOnEntities (const std::string problem_name, const std::string field_name, bool get_low_dim_ents, bool block_name_field_prefix)
 
template<>
MoFEMErrorCode pushMarkDOFsOnEntities (const std::string problem_name, const std::string field_name, bool get_low_dim_ents, bool block_name_field_prefix)
 
template<>
MoFEMErrorCode pushMarkDOFsOnEntities (const std::string problem_name, const std::string field_name, bool get_low_dim_ents, bool block_name_field_prefix)
 
template<>
MoFEMErrorCode pushMarkDOFsOnEntities (const std::string problem_name, const std::string field_name, bool get_low_dim_ents, bool block_name_field_prefix)
 
template<>
MoFEMErrorCode pushMarkDOFsOnEntities (const std::string problem_name, const std::string field_name, const std::string block_name, bool get_low_dim_ents)
 
template<>
MoFEMErrorCode pushMarkDOFsOnEntities (const std::string problem_name, const std::string field_name, bool get_low_dim_ents, bool block_name_field_prefix)
 
template<>
MoFEMErrorCode removeBlockDOFsOnEntities (const std::string problem_name, const std::string field_name, bool get_low_dim_ents, bool block_name_field_prefix, bool is_distributed_mesh)
 
template<>
MoFEMErrorCode pushMarkDOFsOnEntities (const std::string problem_name, const std::string field_name, bool get_low_dim_ents, bool block_name_field_prefix)
 
template<>
MoFEMErrorCode removeBlockDOFsOnEntities (const std::string problem_name, const std::string field_name, bool get_low_dim_ents, bool block_name_field_prefix, bool is_distributed_mesh)
 
template<>
MoFEMErrorCode pushMarkDOFsOnEntities (const std::string problem_name, const std::string field_name, bool get_low_dim_ents, bool block_name_field_prefix)
 
template<>
MoFEMErrorCode removeBlockDOFsOnEntities (const std::string problem_name, const std::string field_name, bool get_low_dim_ents, bool block_name_field_prefix, bool is_distributed_mesh)
 
template<>
MoFEMErrorCode removeBlockDOFsOnEntities (const std::string problem_name, const std::string field_name, bool get_low_dim_ents, bool block_name_field_prefix, bool is_distributed_mesh)
 
template<>
MoFEMErrorCode removeBlockDOFsOnEntities (const std::string problem_name, const std::string field_name, bool get_low_dim_ents, bool block_name_field_prefix, bool is_distributed_mesh)
 
template<>
MoFEMErrorCode removeBlockDOFsOnEntities (const std::string problem_name, const std::string field_name, bool get_low_dim_ents, bool block_name_field_prefix, bool is_distributed_mesh)
 
template<>
MoFEMErrorCode removeBlockDOFsOnEntities (const std::string problem_name, const std::string field_name, bool get_low_dim_ents, bool block_name_field_prefix, bool is_distributed_mesh)
 
template<>
MoFEMErrorCode removeBlockDOFsOnEntities (const std::string problem_name, const std::string block_name, const std::string field_name, bool get_low_dim_ents, bool is_distributed_mesh)
 
template<>
MoFEMErrorCode pushMarkDOFsOnEntities (const std::string problem_name, const std::string field_name, bool get_low_dim_ents, bool block_name_field_prefix)
 
template<>
MoFEMErrorCode pushMarkDOFsOnEntities (const std::string problem_name, const std::string field_name, bool get_low_dim_ents, bool block_name_field_prefix)
 
template<>
MoFEMErrorCode pushMarkDOFsOnEntities (const std::string problem_name, const std::string field_name, bool get_low_dim_ents, bool block_name_field_prefix)
 
template<>
MoFEMErrorCode pushMarkDOFsOnEntities (const std::string problem_name, const std::string field_name, bool get_low_dim_ents, bool block_name_field_prefix)
 
template<>
MoFEMErrorCode pushMarkDOFsOnEntities (const std::string problem_name, const std::string block_name, const std::string field_name, bool get_low_dim_ents)
 
template<>
MoFEMErrorCode pushMarkDOFsOnEntities (const std::string problem_name, const std::string field_name, bool get_low_dim_ents, bool block_name_field_prefix)
 
template<>
MoFEMErrorCode removeBlockDOFsOnEntities (const std::string problem_name, const std::string field_name, bool get_low_dim_ents, bool block_name_field_prefix, bool is_distributed_mesh)
 
template<>
MoFEMErrorCode pushMarkDOFsOnEntities (const std::string problem_name, const std::string field_name, bool get_low_dim_ents, bool block_name_field_prefix)
 
template<>
MoFEMErrorCode removeBlockDOFsOnEntities (const std::string problem_name, const std::string field_name, bool get_low_dim_ents, bool block_name_field_prefix, bool is_distributed_mesh)
 
template<>
MoFEMErrorCode pushMarkDOFsOnEntities (const std::string problem_name, const std::string field_name, bool get_low_dim_ents, bool block_name_field_prefix)
 
template<>
MoFEMErrorCode removeBlockDOFsOnEntities (const std::string problem_name, const std::string field_name, bool get_low_dim_ents, bool block_name_field_prefix, bool is_distributed_mesh)
 
- Public Member Functions inherited from MoFEM::UnknownInterface
virtual MoFEMErrorCode query_interface (boost::typeindex::type_index type_index, UnknownInterface **iface) const =0
 
template<class IFACE >
MoFEMErrorCode registerInterface (bool error_if_registration_failed=true)
 Register interface. More...
 
template<class IFACE >
MoFEMErrorCode getInterface (IFACE *&iface) const
 Get interface refernce to pointer of interface. More...
 
template<class IFACE >
MoFEMErrorCode getInterface (IFACE **const iface) const
 Get interface pointer to pointer of interface. More...
 
template<class IFACE , typename boost::enable_if< boost::is_pointer< IFACE >, int >::type = 0>
IFACE getInterface () const
 Get interface pointer to pointer of interface. More...
 
template<class IFACE , typename boost::enable_if< boost::is_reference< IFACE >, int >::type = 0>
IFACE getInterface () const
 Get reference to interface. More...
 
template<class IFACE >
IFACE * getInterface () const
 Function returning pointer to interface. More...
 
virtual ~UnknownInterface ()=default
 

Static Public Member Functions

static std::pair< std::string, std::string > extractStringFromBlockId (const std::string block_id, const std::string prb_name)
 Extract block name and block name form block id. More...
 
- Static Public Member Functions inherited from MoFEM::UnknownInterface
static MoFEMErrorCode getLibVersion (Version &version)
 Get library version. More...
 
static MoFEMErrorCode getFileVersion (moab::Interface &moab, Version &version)
 Get database major version. More...
 
static MoFEMErrorCode setFileVersion (moab::Interface &moab, Version version=Version(MoFEM_VERSION_MAJOR, MoFEM_VERSION_MINOR, MoFEM_VERSION_BUILD))
 Get database major version. More...
 
static MoFEMErrorCode getInterfaceVersion (Version &version)
 Get database major version. More...
 

Private Attributes

MoFEM::CorecOre
 
BcMapByBlockName bcMapByBlockName
 

Detailed Description

Simple interface for fast problem set-up.

Examples
free_surface.cpp, heat_equation.cpp, nonlinear_elastic.cpp, photon_diffusion.cpp, plastic.cpp, poisson_2d_homogeneous.cpp, poisson_3d_homogeneous.cpp, shallow_wave.cpp, thermo_elastic.cpp, and wave_equation.cpp.

Definition at line 24 of file BcManager.hpp.

Member Typedef Documentation

◆ BcMapByBlockName

using MoFEM::BcManager::BcMapByBlockName = std::map<string, boost::shared_ptr<BCs> >
Todo:
Add markers for standard BCs from cubit on Nodests and Sidesets used bu cubit for displacements, forces, etc. Also add function to add blockset by id and type.

Definition at line 182 of file BcManager.hpp.

◆ BcMarkerPtr

using MoFEM::BcManager::BcMarkerPtr = boost::shared_ptr<std::vector<char unsigned> >

Definition at line 184 of file BcManager.hpp.

Constructor & Destructor Documentation

◆ BcManager()

MoFEM::BcManager::BcManager ( const MoFEM::Core core)

Definition at line 16 of file BcManager.cpp.

16 : cOre(const_cast<Core &>(core)) {
17
18 if (!LogManager::checkIfChannelExist("BcMngWorld")) {
19 auto core_log = logging::core::get();
20
21 core_log->add_sink(
23 core_log->add_sink(
25 core_log->add_sink(
27
28 LogManager::setLog("BcMngWorld");
29 LogManager::setLog("BcMngSync");
30 LogManager::setLog("BcMngSelf");
31
32 MOFEM_LOG_TAG("BcMngWorld", "BcMng");
33 MOFEM_LOG_TAG("BcMngSync", "BcMng");
34 MOFEM_LOG_TAG("BcMngSelf", "BcMng");
35 }
36
37 MOFEM_LOG("BcMngWorld", Sev::noisy) << "BC manager created";
38}
static LoggerType & setLog(const std::string channel)
Set ans resset chanel logger.
Definition: LogManager.cpp:391
#define MOFEM_LOG(channel, severity)
Log.
Definition: LogManager.hpp:308
#define MOFEM_LOG_TAG(channel, tag)
Tag channel.
Definition: LogManager.hpp:339
CoreTmp< 0 > Core
Definition: Core.hpp:1094
MoFEM::Core & cOre
Definition: BcManager.hpp:301
static boost::shared_ptr< SinkType > createSink(boost::shared_ptr< std::ostream > stream_ptr, std::string comm_filter)
Create a sink object.
Definition: LogManager.cpp:300
static boost::shared_ptr< std::ostream > getStrmWorld()
Get the strm world object.
Definition: LogManager.cpp:346
static boost::shared_ptr< std::ostream > getStrmSync()
Get the strm sync object.
Definition: LogManager.cpp:350
static bool checkIfChannelExist(const std::string channel)
Check if channel exist.
Definition: LogManager.cpp:401
static boost::shared_ptr< std::ostream > getStrmSelf()
Get the strm self object.
Definition: LogManager.cpp:342

◆ ~BcManager()

virtual MoFEM::BcManager::~BcManager ( )
virtualdefault

Member Function Documentation

◆ checkBlock() [1/2]

auto MoFEM::BcManager::checkBlock ( const std::pair< std::string, boost::shared_ptr< BCs > > &  bc,
std::string  name 
)
inline

check if given boundary condition name is in the map bc element

Parameters
bcelement of the map
namebc name
Returns
auto

Definition at line 250 of file BcManager.hpp.

251 {
252 auto full_name = std::string("(.*)_") + name + std::string("(.*)");
253 return checkBlock(bc, std::regex(full_name));
254 }
auto checkBlock(const std::pair< string, boost::shared_ptr< BCs > > &bc, std::regex reg)
check if given boundary condition name is in the map bc element
Definition: BcManager.hpp:238

◆ checkBlock() [2/2]

auto MoFEM::BcManager::checkBlock ( const std::pair< string, boost::shared_ptr< BCs > > &  bc,
std::regex  reg 
)
inline

check if given boundary condition name is in the map bc element

Parameters
bcelement of the map
regbc regex
Returns
auto

Definition at line 238 of file BcManager.hpp.

239 {
240 return std::regex_match(bc.first, reg);
241 }

◆ extractStringFromBlockId()

std::pair< std::string, std::string > MoFEM::BcManager::extractStringFromBlockId ( const std::string  block_id,
const std::string  prb_name 
)
static

Extract block name and block name form block id.

Parameters
block_id
prb_name
Returns
std::pair<std::string, std::string>

Definition at line 1045 of file BcManager.cpp.

1046 {
1047
1048 // Assumes that field name is consist with letters and numbers.
1049 // No special characters.
1050 auto field_rgx_str =
1051 (boost::format("%s_([a-zA-Z0-9]*)_(.*)") % prb_name).str();
1052 std::regex field_rgx(field_rgx_str);
1053 std::smatch match_field_name;
1054 std::string field_name;
1055 std::string block_name;
1056
1057 if (std::regex_search(block_id, match_field_name, field_rgx)) {
1058 field_name = match_field_name[1];
1059 block_name = match_field_name[2];
1060 } else {
1062 "Field name and block name can not be resolved");
1063 }
1064
1065 return std::make_pair(field_name, block_name);
1066}
#define CHK_THROW_MESSAGE(err, msg)
Check and throw MoFEM exception.
Definition: definitions.h:595
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:31
constexpr auto field_name

◆ getBcMapByBlockName()

BcMapByBlockName & MoFEM::BcManager::getBcMapByBlockName ( )
inline

Get the bc map.

Returns
auto
Examples
photon_diffusion.cpp.

Definition at line 200 of file BcManager.hpp.

200{ return bcMapByBlockName; }
BcMapByBlockName bcMapByBlockName
Definition: BcManager.hpp:303

◆ getBcStructure()

auto MoFEM::BcManager::getBcStructure ( const std::string  bc_id)
inline

Get the bc structure object.

Parameters
block_name
Returns
auto

Definition at line 191 of file BcManager.hpp.

191 {
192 return bcMapByBlockName.at(bc_id);
193 }

◆ getBlockIS() [1/2]

SmartPetscObj< IS > MoFEM::BcManager::getBlockIS ( const std::string  block_prefix,
const std::string  block_name,
const std::string  field_name,
const std::string  problem_name,
int  lo,
int  hi,
SmartPetscObj< IS >  is_expand = SmartPetscObj<IS>() 
)

Get block IS.

Parameters
block_prefixfor hashmap
block_namefor hash map
field_namefor hash map and IS
problem_namefor IS
lo
hi
is_expandis to extend
Returns
SmartPetscObj<IS>

Definition at line 218 of file BcManager.cpp.

222 {
223 Interface &m_field = cOre;
224
225 const std::string bc_id =
226 block_prefix + "_" + field_name + "_" + block_name + "(.*)";
227
228 Range bc_ents;
229 for (auto bc : getBcMapByBlockName()) {
230 if (std::regex_match(bc.first, std::regex(bc_id))) {
231 bc_ents.merge(*(bc.second->getBcEntsPtr()));
232 MOFEM_LOG("BcMngWorld", Sev::verbose)
233 << "Get entities from block and add to IS. Block name " << bc.first;
234 }
235 }
236
237 SmartPetscObj<IS> is_bc;
238 auto get_is = [&]() {
240 CHKERR m_field.getInterface<CommInterface>()->synchroniseEntities(bc_ents);
241 CHKERR m_field.getInterface<ISManager>()->isCreateProblemFieldAndRank(
242 problem_name, ROW, field_name, lo, hi, is_bc, &bc_ents);
243 if (is_expand) {
244 IS is_tmp;
245 CHKERR ISExpand(is_bc, is_expand, &is_tmp);
246 is_bc = SmartPetscObj<IS>(is_tmp);
247 }
248 CHKERR ISSort(is_bc);
250 };
251
252 if (get_is())
253 CHK_THROW_MESSAGE(MOFEM_DATA_INCONSISTENCY, "IS is not created");
254
255 return is_bc;
256}
@ ROW
Definition: definitions.h:123
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:416
#define CHKERR
Inline error check.
Definition: definitions.h:535
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1975
BcMapByBlockName & getBcMapByBlockName()
Get the bc map.
Definition: BcManager.hpp:200

◆ getBlockIS() [2/2]

SmartPetscObj< IS > MoFEM::BcManager::getBlockIS ( const std::string  problem_name,
const std::string  block_name,
const std::string  field_name,
int  lo,
int  hi,
SmartPetscObj< IS >  is_expand = SmartPetscObj<IS>() 
)

Get block IS.

Parameters
problem_name
block_name
field_name
lo
hi
is_expandis to extend
Returns
SmartPetscObj<IS>

Definition at line 258 of file BcManager.cpp.

261 {
262 return getBlockIS(problem_name, block_name, field_name, problem_name, lo, hi,
263 is_expand);
264}
SmartPetscObj< IS > getBlockIS(const std::string block_prefix, const std::string block_name, const std::string field_name, const std::string problem_name, int lo, int hi, SmartPetscObj< IS > is_expand=SmartPetscObj< IS >())
Get block IS.
Definition: BcManager.cpp:218

◆ getMergedBlocksMarker() [1/3]

BcMarkerPtr MoFEM::BcManager::getMergedBlocksMarker ( const std::vector< BcMarkerPtr > &  boundary_markers_ptr_vec)

Get the Merged Blocks Marker object.

Parameters
boundary_markers_ptr_vecvector of boundary markers to merge
Returns
BcMarkerPtr

◆ getMergedBlocksMarker() [2/3]

BcManager::BcMarkerPtr MoFEM::BcManager::getMergedBlocksMarker ( std::vector< std::regex >  bc_regex_vec)

Get the Merged Boundary Marker object.

Parameters
bc_regex_vecboundary name regex vector
Returns
boundaryMarker

Definition at line 187 of file BcManager.cpp.

187 {
188 BcManager::BcMarkerPtr boundary_marker_ptr;
189 if (bcMapByBlockName.size()) {
190 for (auto b : bcMapByBlockName) {
191 for (auto &reg_name : bc_regex_vec) {
192 if (std::regex_match(b.first, reg_name)) {
193 if (!boundary_marker_ptr)
194 boundary_marker_ptr =
195 boost::make_shared<std::vector<char unsigned>>();
196 boundary_marker_ptr->resize(b.second->bcMarkers.size(), 0);
197 for (int i = 0; i != b.second->bcMarkers.size(); ++i) {
198 (*boundary_marker_ptr)[i] |= b.second->bcMarkers[i];
199 }
200 }
201 }
202 }
203 }
204 return boundary_marker_ptr;
205}
FTensor::Index< 'i', SPACE_DIM > i
boost::shared_ptr< std::vector< char unsigned > > BcMarkerPtr
Definition: BcManager.hpp:184

◆ getMergedBlocksMarker() [3/3]

auto MoFEM::BcManager::getMergedBlocksMarker ( std::vector< string >  bc_names)
inline

Get the Merged Boundary Marker object.

Parameters
bc_namesvector of boundary names
Returns
boundaryMarker

Definition at line 215 of file BcManager.hpp.

215 {
216 std::vector<std::regex> reg_vec(bc_names.size());
217 for (int i = 0; i != bc_names.size(); ++i) {
218 auto full_name = std::string("(.*)_") + bc_names[i] + std::string("(.*)");
219 reg_vec[i] = std::regex(full_name);
220 }
221 return getMergedBlocksMarker(reg_vec);
222 }
BcMarkerPtr getMergedBlocksMarker(std::vector< std::regex > bc_regex_vec)
Get the Merged Boundary Marker object.
Definition: BcManager.cpp:187

◆ getOptions()

MoFEMErrorCode MoFEM::BcManager::getOptions ( )

get options

Returns
error code

Definition at line 40 of file BcManager.cpp.

40 {
41 PetscBool flg = PETSC_TRUE;
43 ierr = PetscOptionsBegin(PETSC_COMM_WORLD, "", "BcManager options", "none");
44 ierr = PetscOptionsEnd();
47}
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:447
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:483
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:440
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Exceptions.hpp:76

◆ popMarkDOFsOnEntities()

boost::shared_ptr< BcManager::BCs > MoFEM::BcManager::popMarkDOFsOnEntities ( const std::string  block_name)

Get bc data and remove element.

Parameters
block_name
Returns
boost::shared_ptr<BCs>

Definition at line 176 of file BcManager.cpp.

176 {
177 auto bc_it = bcMapByBlockName.find(block_name);
178 if (bc_it != bcMapByBlockName.end()) {
179 auto bc = bc_it->second;
180 bcMapByBlockName.erase(bc_it);
181 return bc;
182 }
183 return boost::shared_ptr<BCs>();
184}

◆ pushMarkDOFsOnEntities() [1/19]

template<>
MoFEMErrorCode MoFEM::BcManager::pushMarkDOFsOnEntities ( const std::string  problem_name,
const std::string  block_name,
const std::string  field_name,
bool  get_low_dim_ents 
)

Definition at line 915 of file BcManager.cpp.

917 {
918 Interface &m_field = cOre;
919 auto prb_mng = m_field.getInterface<ProblemsManager>();
921
922 CHKERR pushMarkDOFsOnEntities(problem_name, block_name, field_name, 0,
923 MAX_DOFS_ON_ENTITY, get_low_dim_ents);
924
925 auto regex_str =
926 (boost::format("%s_%s_%s(.*)") % problem_name % field_name % block_name)
927 .str();
928
929 for (auto &m : bcMapByBlockName) {
930
931 auto &bc_id = m.first;
932
933 if (std::regex_match(bc_id, std::regex(regex_str))) {
934
935 auto &bc = m.second;
936 bc->tempBcPtr = boost::make_shared<TemperatureCubitBcData>();
937 bc->tempBcPtr->data.flag1 = 1;
938 if (bc->bcAttributes.empty()) {
939 bc->tempBcPtr->data.value1 = 0;
940 MOFEM_LOG("BcMngWorld", Sev::warning)
941 << "Expected one atributes on block but have "
942 << bc->bcAttributes.size();
943 } else if (bc->bcAttributes.size() >= 1) {
944 bc->tempBcPtr->data.value1 = bc->bcAttributes[0];
945 }
946 }
947 }
948
950}
#define MAX_DOFS_ON_ENTITY
Maximal number of DOFs on entity.
Definition: definitions.h:236
FTensor::Index< 'm', SPACE_DIM > m
MoFEMErrorCode pushMarkDOFsOnEntities(const std::string problem_name, const std::string block_name, const std::string field_name, int lo, int hi, bool get_low_dim_ents=true)
Mark block DOFs.
Definition: BcManager.cpp:87
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface refernce to pointer of interface.

◆ pushMarkDOFsOnEntities() [2/19]

template<typename T >
MoFEMErrorCode MoFEM::BcManager::pushMarkDOFsOnEntities ( const std::string  problem_name,
const std::string  block_name,
const std::string  field_name,
bool  get_low_dim_ents = true 
)

Mark block DOFs.

Parameters
problem_name
field_name
block_name
get_low_dim_entsget lower dimension entities
Returns
MoFEMErrorCode

◆ pushMarkDOFsOnEntities() [3/19]

MoFEMErrorCode MoFEM::BcManager::pushMarkDOFsOnEntities ( const std::string  problem_name,
const std::string  block_name,
const std::string  field_name,
int  lo,
int  hi,
bool  get_low_dim_ents = true 
)

Mark block DOFs.

Parameters
problem_name
block_name
field_name
lolowest coefficient
hihighest coefficient
get_low_dim_entsget lower dimension entities field name
Returns
MoFEMErrorCode

Definition at line 87 of file BcManager.cpp.

91 {
92 Interface &m_field = cOre;
93 auto prb_mng = m_field.getInterface<ProblemsManager>();
95
96 auto get_dim = [&](const Range &ents) {
97 for (auto d : {3, 2, 1})
98 if (ents.num_of_dimension(d))
99 return d;
100 return 0;
101 };
102
103 auto get_adj_ents = [&](const Range &ents) {
104 Range verts;
105 CHKERR m_field.get_moab().get_connectivity(ents, verts, true);
106 const auto dim = get_dim(ents);
107 for (size_t d = 1; d < dim; ++d)
108 CHKERR m_field.get_moab().get_adjacencies(ents, d, false, verts,
109 moab::Interface::UNION);
110 verts.merge(ents);
111 CHKERR m_field.getInterface<CommInterface>()->synchroniseEntities(verts);
112 return verts;
113 };
114
115 auto fix_disp = [&]() {
117
118 auto mark_fix_dofs = [&](std::vector<unsigned char> &marked_field_dofs,
119 const auto lo, const auto hi) {
120 return prb_mng->modifyMarkDofs(problem_name, ROW, field_name, lo, hi,
122 marked_field_dofs);
123 };
124
125 auto iterate_meshsets = [&](auto &&meshset_vec_ptr) {
127 for (auto m : meshset_vec_ptr) {
128 auto bc = boost::make_shared<BCs>();
129 CHKERR m_field.get_moab().get_entities_by_handle(m->getMeshset(),
130 bc->bcEnts, true);
131 CHKERR m_field.getInterface<CommInterface>()->synchroniseEntities(
132 bc->bcEnts);
133 CHKERR m->getAttributes(bc->bcAttributes);
134
135 MOFEM_LOG("BcMngWorld", Sev::verbose)
136 << "Found block " << m->getName() << " number of entities "
137 << bc->bcEnts.size() << " number of attributes "
138 << bc->bcAttributes.size() << " highest dim of entities "
139 << get_dim(bc->bcEnts);
140 CHKERR mark_fix_dofs(bc->bcMarkers, lo, hi);
141 if (get_low_dim_ents) {
142 auto low_dim_ents = get_adj_ents(bc->bcEnts);
143 CHKERR prb_mng->markDofs(problem_name, ROW, ProblemsManager::AND,
144 low_dim_ents, bc->bcMarkers);
145 bc->bcEnts.swap(low_dim_ents);
146 } else
147 CHKERR prb_mng->markDofs(problem_name, ROW, ProblemsManager::AND,
148 bc->bcEnts, bc->bcMarkers);
149
150 const std::string bc_id =
151 problem_name + "_" + field_name + "_" + m->getName();
152 bcMapByBlockName[bc_id] = bc;
153 }
155 };
156
157 CHKERR iterate_meshsets(
158
159 m_field.getInterface<MeshsetsManager>()->getCubitMeshsetPtr(std::regex(
160
161 (boost::format("%s(.*)") % block_name).str()
162
163 ))
164
165 );
166
168 };
169
170 CHKERR fix_disp();
171
173}
const int dim
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

◆ pushMarkDOFsOnEntities() [4/19]

template<>
MoFEMErrorCode MoFEM::BcManager::pushMarkDOFsOnEntities ( const std::string  problem_name,
const std::string  field_name,
bool  get_low_dim_ents,
bool  block_name_field_prefix 
)

◆ pushMarkDOFsOnEntities() [5/19]

template<>
MoFEMErrorCode MoFEM::BcManager::pushMarkDOFsOnEntities ( const std::string  problem_name,
const std::string  field_name,
bool  get_low_dim_ents,
bool  block_name_field_prefix 
)

◆ pushMarkDOFsOnEntities() [6/19]

template<>
MoFEMErrorCode MoFEM::BcManager::pushMarkDOFsOnEntities ( const std::string  problem_name,
const std::string  field_name,
bool  get_low_dim_ents,
bool  block_name_field_prefix 
)

◆ pushMarkDOFsOnEntities() [7/19]

template<>
MoFEMErrorCode MoFEM::BcManager::pushMarkDOFsOnEntities ( const std::string  problem_name,
const std::string  field_name,
bool  get_low_dim_ents,
bool  block_name_field_prefix 
)

◆ pushMarkDOFsOnEntities() [8/19]

template<>
MoFEMErrorCode MoFEM::BcManager::pushMarkDOFsOnEntities ( const std::string  problem_name,
const std::string  field_name,
bool  get_low_dim_ents,
bool  block_name_field_prefix 
)

◆ pushMarkDOFsOnEntities() [9/19]

template<>
MoFEMErrorCode MoFEM::BcManager::pushMarkDOFsOnEntities ( const std::string  problem_name,
const std::string  field_name,
bool  get_low_dim_ents,
bool  block_name_field_prefix 
)

◆ pushMarkDOFsOnEntities() [10/19]

template<>
MoFEMErrorCode MoFEM::BcManager::pushMarkDOFsOnEntities ( const std::string  problem_name,
const std::string  field_name,
bool  get_low_dim_ents,
bool  block_name_field_prefix 
)

◆ pushMarkDOFsOnEntities() [11/19]

template<>
MoFEMErrorCode MoFEM::BcManager::pushMarkDOFsOnEntities ( const std::string  problem_name,
const std::string  field_name,
bool  get_low_dim_ents,
bool  block_name_field_prefix 
)

Definition at line 539 of file BcManager.cpp.

542 {
543 Interface &m_field = cOre;
544 auto prb_mng = m_field.getInterface<ProblemsManager>();
546
547 if (block_name_field_prefix)
548 MOFEM_LOG("BcMngWorld", Sev::warning)
549 << "Argument block_name_field_prefix=true has no effect";
550
551 auto get_dim = [&](const Range &ents) {
552 for (auto d : {3, 2, 1})
553 if (ents.num_of_dimension(d))
554 return d;
555 return 0;
556 };
557
558 auto get_adj_ents = [&](const Range &ents) {
559 Range verts;
560 CHKERR m_field.get_moab().get_connectivity(ents, verts, true);
561 const auto dim = get_dim(ents);
562 for (size_t d = 1; d < dim; ++d)
563 CHKERR m_field.get_moab().get_adjacencies(ents, d, false, verts,
564 moab::Interface::UNION);
565 verts.merge(ents);
566 CHKERR m_field.getInterface<CommInterface>()->synchroniseEntities(verts);
567 return verts;
568 };
569
570 auto fix_disp = [&]() {
572
573 auto iterate_meshsets = [&](auto &&meshset_vec_ptr) {
575 for (auto m : meshset_vec_ptr) {
576 auto bc = boost::make_shared<BCs>();
577 CHKERR m_field.get_moab().get_entities_by_handle(m->getMeshset(),
578 bc->bcEnts, true);
579
580 bc->dispBcPtr = boost::make_shared<DisplacementCubitBcData>();
581 CHKERR m->getBcDataStructure(*(bc->dispBcPtr));
582
583 MOFEM_LOG("BcMngWorld", Sev::verbose)
584 << "Found block DISPLACEMENTSET number of entities "
585 << bc->bcEnts.size() << " highest dim of entities "
586 << get_dim(bc->bcEnts);
587 MOFEM_LOG("BcMngWorld", Sev::verbose) << *bc->dispBcPtr;
588
589 if (bc->dispBcPtr->data.flag1)
590 CHKERR prb_mng->modifyMarkDofs(problem_name, ROW, field_name, 0, 0,
592 bc->bcMarkers);
593 if (bc->dispBcPtr->data.flag2)
594 CHKERR prb_mng->modifyMarkDofs(problem_name, ROW, field_name, 1, 1,
596 bc->bcMarkers);
597 if (bc->dispBcPtr->data.flag3)
598 CHKERR prb_mng->modifyMarkDofs(problem_name, ROW, field_name, 2, 2,
600 bc->bcMarkers);
601
602 if (get_low_dim_ents) {
603 auto low_dim_ents = get_adj_ents(bc->bcEnts);
604 CHKERR prb_mng->markDofs(problem_name, ROW, ProblemsManager::AND,
605 low_dim_ents, bc->bcMarkers);
606 bc->bcEnts.swap(low_dim_ents);
607 } else
608 CHKERR prb_mng->markDofs(problem_name, ROW, ProblemsManager::AND,
609 bc->bcEnts, bc->bcMarkers);
610
611 const std::string bc_id =
612 problem_name + "_" + field_name + "_DISPLACEMENTSET" +
613 boost::lexical_cast<std::string>(m->getMeshsetId());
614 bcMapByBlockName[bc_id] = bc;
615 }
617 };
618
619 CHKERR iterate_meshsets(
620
621 m_field.getInterface<MeshsetsManager>()->getCubitMeshsetPtr(
623
624 );
625
627 };
628
629 CHKERR fix_disp();
630
632}
@ NODESET
Definition: definitions.h:146
@ DISPLACEMENTSET
Definition: definitions.h:150

◆ pushMarkDOFsOnEntities() [12/19]

template<>
MoFEMErrorCode MoFEM::BcManager::pushMarkDOFsOnEntities ( const std::string  problem_name,
const std::string  field_name,
bool  get_low_dim_ents,
bool  block_name_field_prefix 
)

Definition at line 635 of file BcManager.cpp.

637 {
638 Interface &m_field = cOre;
639 auto prb_mng = m_field.getInterface<ProblemsManager>();
641
642 if (block_name_field_prefix)
643 MOFEM_LOG("BcMngWorld", Sev::warning)
644 << "Argument block_name_field_prefix=true has no effect";
645
646 auto get_dim = [&](const Range &ents) {
647 for (auto d : {3, 2, 1})
648 if (ents.num_of_dimension(d))
649 return d;
650 return 0;
651 };
652
653 auto get_adj_ents = [&](const Range &ents) {
654 Range verts;
655 CHKERR m_field.get_moab().get_connectivity(ents, verts, true);
656 const auto dim = get_dim(ents);
657 for (size_t d = 1; d < dim; ++d)
658 CHKERR m_field.get_moab().get_adjacencies(ents, d, false, verts,
659 moab::Interface::UNION);
660 verts.merge(ents);
661 CHKERR m_field.getInterface<CommInterface>()->synchroniseEntities(verts);
662 return verts;
663 };
664
665 auto fix_temp = [&]() {
667
668 auto iterate_meshsets = [&](auto &&meshset_vec_ptr) {
670 for (auto m : meshset_vec_ptr) {
671 auto bc = boost::make_shared<BCs>();
672 CHKERR m_field.get_moab().get_entities_by_handle(m->getMeshset(),
673 bc->bcEnts, true);
674 bc->tempBcPtr = boost::make_shared<TemperatureCubitBcData>();
675 CHKERR m->getBcDataStructure(*(bc->tempBcPtr));
676
677 MOFEM_LOG("BcMngWorld", Sev::verbose)
678 << "Found block TEMPERATURESET number of entities "
679 << bc->bcEnts.size() << " highest dim of entities "
680 << get_dim(bc->bcEnts);
681 MOFEM_LOG("BcMngWorld", Sev::verbose) << *bc->tempBcPtr;
682
683 CHKERR prb_mng->modifyMarkDofs(
684 problem_name, ROW, field_name, 0, MAX_DOFS_ON_ENTITY,
685 ProblemsManager::MarkOP::OR, 1, bc->bcMarkers);
686
687 if (get_low_dim_ents) {
688 auto low_dim_ents = get_adj_ents(bc->bcEnts);
689 CHKERR prb_mng->markDofs(problem_name, ROW, ProblemsManager::AND,
690 low_dim_ents, bc->bcMarkers);
691 bc->bcEnts.swap(low_dim_ents);
692 } else
693 CHKERR prb_mng->markDofs(problem_name, ROW, ProblemsManager::AND,
694 bc->bcEnts, bc->bcMarkers);
695
696 const std::string bc_id =
697 problem_name + "_" + field_name + "_TEMPERATURESET" +
698 boost::lexical_cast<std::string>(m->getMeshsetId());
699 bcMapByBlockName[bc_id] = bc;
700 }
702 };
703
704 CHKERR iterate_meshsets(
705
706 m_field.getInterface<MeshsetsManager>()->getCubitMeshsetPtr(
708
709 );
710
712 };
713
714 CHKERR fix_temp();
715
717}
@ TEMPERATURESET
Definition: definitions.h:155

◆ pushMarkDOFsOnEntities() [13/19]

template<>
MoFEMErrorCode MoFEM::BcManager::pushMarkDOFsOnEntities ( const std::string  problem_name,
const std::string  field_name,
bool  get_low_dim_ents,
bool  block_name_field_prefix 
)

Definition at line 720 of file BcManager.cpp.

722 {
723 Interface &m_field = cOre;
724 auto prb_mng = m_field.getInterface<ProblemsManager>();
726
727 if (block_name_field_prefix)
728 MOFEM_LOG("BcMngWorld", Sev::warning)
729 << "Argument block_name_field_prefix=true has no effect";
730
731 auto get_dim = [&](const Range &ents) {
732 for (auto d : {3, 2, 1})
733 if (ents.num_of_dimension(d))
734 return d;
735 return 0;
736 };
737
738 auto get_adj_ents = [&](const Range &ents) {
739 Range verts;
740 CHKERR m_field.get_moab().get_connectivity(ents, verts, true);
741 const auto dim = get_dim(ents);
742 for (size_t d = 1; d < dim; ++d)
743 CHKERR m_field.get_moab().get_adjacencies(ents, d, false, verts,
744 moab::Interface::UNION);
745 verts.merge(ents);
746 CHKERR m_field.getInterface<CommInterface>()->synchroniseEntities(verts);
747 return verts;
748 };
749
750 auto fix_disp = [&]() {
752
753 auto iterate_meshsets = [&](auto &&meshset_vec_ptr) {
755 for (auto m : meshset_vec_ptr) {
756 auto bc = boost::make_shared<BCs>();
757 CHKERR m_field.get_moab().get_entities_by_handle(m->getMeshset(),
758 bc->bcEnts, true);
759 bc->heatFluxBcPtr = boost::make_shared<HeatFluxCubitBcData>();
760 CHKERR m->getBcDataStructure(*(bc->heatFluxBcPtr));
761
762 MOFEM_LOG("BcMngWorld", Sev::verbose)
763 << "Found block HEATFLUX number of entities " << bc->bcEnts.size()
764 << " highest dim of entities " << get_dim(bc->bcEnts);
765 MOFEM_LOG("BcMngWorld", Sev::verbose) << *bc->heatFluxBcPtr;
766
767 CHKERR prb_mng->modifyMarkDofs(
768 problem_name, ROW, field_name, 0, MAX_DOFS_ON_ENTITY,
769 ProblemsManager::MarkOP::OR, 1, bc->bcMarkers);
770
771 if (get_low_dim_ents) {
772 auto low_dim_ents = get_adj_ents(bc->bcEnts);
773 CHKERR prb_mng->markDofs(problem_name, ROW, ProblemsManager::AND,
774 low_dim_ents, bc->bcMarkers);
775 bc->bcEnts.swap(low_dim_ents);
776 } else
777 CHKERR prb_mng->markDofs(problem_name, ROW, ProblemsManager::AND,
778 bc->bcEnts, bc->bcMarkers);
779
780 const std::string bc_id =
781 problem_name + "_" + field_name + "_HEATFLUXSET" +
782 boost::lexical_cast<std::string>(m->getMeshsetId());
783 bcMapByBlockName[bc_id] = bc;
784 }
786 };
787
788 CHKERR iterate_meshsets(
789
790 m_field.getInterface<MeshsetsManager>()->getCubitMeshsetPtr(SIDESET |
792
793 );
794
796 };
797
798 CHKERR fix_disp();
799
801}
@ HEATFLUXSET
Definition: definitions.h:156
@ SIDESET
Definition: definitions.h:147

◆ pushMarkDOFsOnEntities() [14/19]

template<>
MoFEMErrorCode MoFEM::BcManager::pushMarkDOFsOnEntities ( const std::string  problem_name,
const std::string  field_name,
bool  get_low_dim_ents,
bool  block_name_field_prefix 
)

Definition at line 804 of file BcManager.cpp.

806 {
807 Interface &m_field = cOre;
808 auto prb_mng = m_field.getInterface<ProblemsManager>();
810
811 if (block_name_field_prefix) {
813 problem_name, (boost::format("%s_FIX_X") % field_name).str(),
814 field_name, 0, 0, get_low_dim_ents);
816 problem_name, (boost::format("%s_FIX_Y") % field_name).str(),
817 field_name, 1, 1, get_low_dim_ents);
819 problem_name, (boost::format("%s_FIX_Z") % field_name).str(),
820 field_name, 2, 2, get_low_dim_ents);
822 problem_name, (boost::format("%s_FIX_ALL") % field_name).str(),
823 field_name, 0, MAX_DOFS_ON_ENTITY, get_low_dim_ents);
824 } else {
825 CHKERR pushMarkDOFsOnEntities(problem_name, "FIX_X", field_name, 0, 0,
826 get_low_dim_ents);
827 CHKERR pushMarkDOFsOnEntities(problem_name, "FIX_Y", field_name, 1, 1,
828 get_low_dim_ents);
829 CHKERR pushMarkDOFsOnEntities(problem_name, "FIX_Z", field_name, 2, 2,
830 get_low_dim_ents);
831 CHKERR pushMarkDOFsOnEntities(problem_name, "FIX_ALL", field_name, 0,
832 MAX_DOFS_ON_ENTITY, get_low_dim_ents);
833 }
834
835 std::string regex_str;
836 if (block_name_field_prefix) {
837 regex_str = (boost::format("%s_%s_%s_(.*)") % problem_name % field_name %
839 .str();
840 } else {
841 regex_str = (boost::format("%s_%s_(.*)") % problem_name % field_name).str();
842 }
843
844 for (auto &m : bcMapByBlockName) {
845 auto &bc_id = m.first;
846 if (std::regex_match(bc_id, std::regex(regex_str))) {
847 auto &bc = m.second;
848 if (std::regex_match(bc_id, std::regex("(.*)_FIX_X(.*)"))) {
849 bc->dispBcPtr = boost::make_shared<DisplacementCubitBcData>();
850 bc->dispBcPtr->data.flag1 = 1;
851 if (bc->bcAttributes.empty()) {
852 bc->dispBcPtr->data.value1 = 0;
853 MOFEM_LOG("BcMngWorld", Sev::warning)
854 << "Expected one atributes on block but have "
855 << bc->bcAttributes.size();
856 } else if (bc->bcAttributes.size() >= 1) {
857 bc->dispBcPtr->data.value1 = bc->bcAttributes[0];
858 }
859 MOFEM_LOG("BcMngWorld", Sev::inform) << "Add X " << bc_id;
860 MOFEM_LOG("BcMngWorld", Sev::inform) << *bc->dispBcPtr;
861 } else if (std::regex_match(bc_id, std::regex("(.*)_FIX_Y(.*)"))) {
862 bc->dispBcPtr = boost::make_shared<DisplacementCubitBcData>();
863 bc->dispBcPtr->data.flag2 = 1;
864 if (bc->bcAttributes.empty()) {
865 bc->dispBcPtr->data.value2 = 0;
866 MOFEM_LOG("BcMngWorld", Sev::warning)
867 << "Expected one atributes on block but have "
868 << bc->bcAttributes.size();
869 } else if (bc->bcAttributes.size() == 1) {
870 bc->dispBcPtr->data.value2 = bc->bcAttributes[0];
871 } else if (bc->bcAttributes.size() >= 2) {
872 bc->dispBcPtr->data.value2 = bc->bcAttributes[1];
873 }
874 MOFEM_LOG("BcMngWorld", Sev::inform) << "Add Y " << bc_id;
875 MOFEM_LOG("BcMngWorld", Sev::inform) << *(bc->dispBcPtr);
876 } else if (std::regex_match(bc_id, std::regex("(.*)_FIX_Z(.*)"))) {
877 bc->dispBcPtr = boost::make_shared<DisplacementCubitBcData>();
878 bc->dispBcPtr->data.flag3 = 1;
879 if (bc->bcAttributes.empty()) {
880 bc->dispBcPtr->data.value3 = 0;
881 MOFEM_LOG("BcMngWorld", Sev::warning)
882 << "Expected one atributes on block but have "
883 << bc->bcAttributes.size();
884 } else if (bc->bcAttributes.size() == 1) {
885 bc->dispBcPtr->data.value3 = bc->bcAttributes[0];
886 } else if (bc->bcAttributes.size() == 3) {
887 bc->dispBcPtr->data.value3 = bc->bcAttributes[2];
888 }
889 MOFEM_LOG("BcMngWorld", Sev::inform) << "Add Z " << bc_id;
890 MOFEM_LOG("BcMngWorld", Sev::inform) << *(bc->dispBcPtr);
891 } else if (std::regex_match(bc_id, std::regex("(.*)_FIX_ALL(.*)"))) {
892 bc->dispBcPtr = boost::make_shared<DisplacementCubitBcData>();
893 bc->dispBcPtr->data.flag1 = 1;
894 bc->dispBcPtr->data.flag2 = 1;
895 bc->dispBcPtr->data.flag3 = 1;
896 if (bc->bcAttributes.size() >= 1) {
897 bc->dispBcPtr->data.value1 = bc->bcAttributes[0];
898 }
899 if (bc->bcAttributes.size() >= 2) {
900 bc->dispBcPtr->data.value2 = bc->bcAttributes[1];
901 }
902 if (bc->bcAttributes.size() >= 3) {
903 bc->dispBcPtr->data.value3 = bc->bcAttributes[2];
904 }
905 MOFEM_LOG("BcMngWorld", Sev::inform) << "Add ALL " << bc_id;
906 MOFEM_LOG("BcMngWorld", Sev::inform) << *(bc->dispBcPtr);
907 }
908 }
909 }
910
912}

◆ pushMarkDOFsOnEntities() [15/19]

template<>
MoFEMErrorCode MoFEM::BcManager::pushMarkDOFsOnEntities ( const std::string  problem_name,
const std::string  field_name,
bool  get_low_dim_ents,
bool  block_name_field_prefix 
)

Definition at line 953 of file BcManager.cpp.

955 {
957 CHKERR pushMarkDOFsOnEntities<BcMeshsetType<DISPLACEMENTSET>>(
958 problem_name, field_name, get_low_dim_ents, block_name_field_prefix);
959 CHKERR pushMarkDOFsOnEntities<BcVectorMeshsetType<BLOCKSET>>(
960 problem_name, field_name, get_low_dim_ents, block_name_field_prefix);
962}

◆ pushMarkDOFsOnEntities() [16/19]

template<>
MoFEMErrorCode MoFEM::BcManager::pushMarkDOFsOnEntities ( const std::string  problem_name,
const std::string  field_name,
bool  get_low_dim_ents,
bool  block_name_field_prefix 
)

Definition at line 980 of file BcManager.cpp.

982 {
984 CHKERR pushMarkDOFsOnEntities<BcMeshsetType<TEMPERATURESET>>(
985 problem_name, field_name, get_low_dim_ents, block_name_field_prefix);
986
987 auto get_block_name = [&]() {
988 if (block_name_field_prefix)
989 return (boost::format("%s_FIX_SCALAR") % field_name).str();
990 else
991 return field_name;
992 };
993
994 CHKERR pushMarkDOFsOnEntities<BcScalarMeshsetType<BLOCKSET>>(
995 problem_name, get_block_name(), field_name, get_low_dim_ents);
997}

◆ pushMarkDOFsOnEntities() [17/19]

template<>
MoFEMErrorCode MoFEM::BcManager::pushMarkDOFsOnEntities ( const std::string  problem_name,
const std::string  field_name,
bool  get_low_dim_ents,
bool  block_name_field_prefix 
)

Definition at line 1023 of file BcManager.cpp.

1025 {
1027 CHKERR pushMarkDOFsOnEntities<BcMeshsetType<HEATFLUXSET>>(
1028 problem_name, field_name, get_low_dim_ents, block_name_field_prefix);
1030}

◆ pushMarkDOFsOnEntities() [18/19]

template<typename T >
MoFEMErrorCode MoFEM::BcManager::pushMarkDOFsOnEntities ( const std::string  problem_name,
const std::string  field_name,
bool  get_low_dim_ents = true,
bool  block_name_field_prefix = false 
)

Mark block DOFs.

Parameters
problem_name
field_name
get_low_dim_entsget lower dimension entities
block_name_field_prefix
Returns
MoFEMErrorCode

◆ pushMarkDOFsOnEntities() [19/19]

template<>
MoFEMErrorCode MoFEM::BcManager::pushMarkDOFsOnEntities ( const std::string  problem_name,
const std::string  field_name,
const std::string  block_name,
bool  get_low_dim_ents 
)

◆ query_interface()

MoFEMErrorCode MoFEM::BcManager::query_interface ( boost::typeindex::type_index  type_index,
UnknownInterface **  iface 
) const
virtual

Implements MoFEM::UnknownInterface.

Definition at line 9 of file BcManager.cpp.

10 {
12 *iface = const_cast<BcManager *>(this);
14}
BcManager(const MoFEM::Core &core)
Definition: BcManager.cpp:16

◆ removeBlockDOFsOnEntities() [1/19]

template<>
MoFEMErrorCode MoFEM::BcManager::removeBlockDOFsOnEntities ( const std::string  problem_name,
const std::string  block_name,
const std::string  field_name,
bool  get_low_dim_ents,
bool  is_distributed_mesh 
)

◆ removeBlockDOFsOnEntities() [2/19]

template<>
MoFEMErrorCode MoFEM::BcManager::removeBlockDOFsOnEntities ( const std::string  problem_name,
const std::string  block_name,
const std::string  field_name,
bool  get_low_dim_ents,
bool  is_distributed_mesh 
)

Definition at line 484 of file BcManager.cpp.

488 {
489 Interface &m_field = cOre;
490 auto prb_mng = m_field.getInterface<ProblemsManager>();
492
493 CHKERR pushMarkDOFsOnEntities<BcScalarMeshsetType<BLOCKSET>>(
494 problem_name, block_name, field_name, get_low_dim_ents);
495
496 Range ents_to_remove;
497
498 for (auto m :
499
500 m_field.getInterface<MeshsetsManager>()->getCubitMeshsetPtr(BLOCKSET)) {
501
502 std::string bc_id = problem_name + "_" + field_name + "_" + m->getName();
503
504 auto str = boost::format("%s_%s_%s(.*)")
505
506 % problem_name % field_name % block_name;
507
508 if (std::regex_match(bc_id, std::regex(str.str()))) {
509
510 auto bc = bcMapByBlockName.at(bc_id);
511
512 if (auto disp_bc = bc->tempBcPtr) {
513 if (disp_bc->data.flag1) {
514 ents_to_remove.merge(bc->bcEnts);
515 }
516 } else {
517 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
518 "BC type not implemented");
519 }
520 }
521 }
522
523 auto remove_dofs_on_ents = [&](const Range &ents, const int lo,
524 const int hi) {
525 if (is_distributed_mesh)
526 return prb_mng->removeDofsOnEntities(problem_name, field_name, ents, lo,
527 hi);
528 else
529 return prb_mng->removeDofsOnEntitiesNotDistributed(
530 problem_name, field_name, ents, lo, hi);
531 };
532
533 CHKERR remove_dofs_on_ents(ents_to_remove, 0, MAX_DOFS_ON_ENTITY);
534
536}
@ BLOCKSET
Definition: definitions.h:148

◆ removeBlockDOFsOnEntities() [3/19]

template<typename T >
MoFEMErrorCode MoFEM::BcManager::removeBlockDOFsOnEntities ( const std::string  problem_name,
const std::string  block_name,
const std::string  field_name,
bool  get_low_dim_ents = true,
bool  is_distributed_mesh = true 
)

Mark block DOFs.

Template Parameters
BCSET
Parameters
problem_name
field_name
block_name
get_low_dim_ents
is_distributed_mesh
Returns
MoFEMErrorCode

◆ removeBlockDOFsOnEntities() [4/19]

MoFEMErrorCode MoFEM::BcManager::removeBlockDOFsOnEntities ( const std::string  problem_name,
const std::string  block_name,
const std::string  field_name,
int  lo,
int  hi,
bool  get_low_dim_ents = true,
bool  is_distributed_mesh = true 
)

Remove DOFs from problem.

Parameters
problem_name
block_name
field_name
lolowest coefficient
hihighest coefficient
get_low_dim_entsget lower dimension entities
block_name_field_prefixblock name is expected to have prefix with field name
is_distributed_meshdistributed mesh
Returns
MoFEMErrorCode

Definition at line 49 of file BcManager.cpp.

52 {
53 Interface &m_field = cOre;
54 auto prb_mng = m_field.getInterface<ProblemsManager>();
56
57 CHKERR pushMarkDOFsOnEntities(problem_name, block_name, field_name, lo, hi,
58 get_low_dim_ents);
59
60 auto remove_dofs_on_ents = [&](const Range &ents, const int lo,
61 const int hi) {
62 if (is_distributed_mesh)
63 return prb_mng->removeDofsOnEntities(problem_name, field_name, ents, lo,
64 hi);
65 else
66 return prb_mng->removeDofsOnEntitiesNotDistributed(
67 problem_name, field_name, ents, lo, hi);
68 };
69
70 for (auto m :
71 m_field.getInterface<MeshsetsManager>()->getCubitMeshsetPtr(std::regex(
72
73 (boost::format("%s(.*)") % block_name).str()
74
75 ))
76
77 ) {
78 const std::string bc_id =
79 problem_name + "_" + field_name + "_" + m->getName();
80 CHKERR remove_dofs_on_ents(bcMapByBlockName.at(bc_id)->bcEnts, lo, hi);
81 bcMapByBlockName.at(bc_id)->bcMarkers = std::vector<unsigned char>();
82 }
83
85}

◆ removeBlockDOFsOnEntities() [5/19]

template<>
MoFEMErrorCode MoFEM::BcManager::removeBlockDOFsOnEntities ( const std::string  problem_name,
const std::string  field_name,
bool  get_low_dim_ents,
bool  block_name_field_prefix,
bool  is_distributed_mesh 
)

◆ removeBlockDOFsOnEntities() [6/19]

template<>
MoFEMErrorCode MoFEM::BcManager::removeBlockDOFsOnEntities ( const std::string  problem_name,
const std::string  field_name,
bool  get_low_dim_ents,
bool  block_name_field_prefix,
bool  is_distributed_mesh 
)

◆ removeBlockDOFsOnEntities() [7/19]

template<>
MoFEMErrorCode MoFEM::BcManager::removeBlockDOFsOnEntities ( const std::string  problem_name,
const std::string  field_name,
bool  get_low_dim_ents,
bool  block_name_field_prefix,
bool  is_distributed_mesh 
)

◆ removeBlockDOFsOnEntities() [8/19]

template<>
MoFEMErrorCode MoFEM::BcManager::removeBlockDOFsOnEntities ( const std::string  problem_name,
const std::string  field_name,
bool  get_low_dim_ents,
bool  block_name_field_prefix,
bool  is_distributed_mesh 
)

◆ removeBlockDOFsOnEntities() [9/19]

template<>
MoFEMErrorCode MoFEM::BcManager::removeBlockDOFsOnEntities ( const std::string  problem_name,
const std::string  field_name,
bool  get_low_dim_ents,
bool  block_name_field_prefix,
bool  is_distributed_mesh 
)

◆ removeBlockDOFsOnEntities() [10/19]

template<>
MoFEMErrorCode MoFEM::BcManager::removeBlockDOFsOnEntities ( const std::string  problem_name,
const std::string  field_name,
bool  get_low_dim_ents,
bool  block_name_field_prefix,
bool  is_distributed_mesh 
)

◆ removeBlockDOFsOnEntities() [11/19]

template<>
MoFEMErrorCode MoFEM::BcManager::removeBlockDOFsOnEntities ( const std::string  problem_name,
const std::string  field_name,
bool  get_low_dim_ents,
bool  block_name_field_prefix,
bool  is_distributed_mesh 
)

◆ removeBlockDOFsOnEntities() [12/19]

template<>
MoFEMErrorCode MoFEM::BcManager::removeBlockDOFsOnEntities ( const std::string  problem_name,
const std::string  field_name,
bool  get_low_dim_ents,
bool  block_name_field_prefix,
bool  is_distributed_mesh 
)

Definition at line 267 of file BcManager.cpp.

271 {
272 Interface &m_field = cOre;
273 auto prb_mng = m_field.getInterface<ProblemsManager>();
275
276 CHKERR pushMarkDOFsOnEntities<BcMeshsetType<DISPLACEMENTSET>>(
277 problem_name, field_name, get_low_dim_ents, block_name_field_prefix);
278
279 std::array<Range, 3> ents_to_remove;
280
281 for (auto m :
282
283 m_field.getInterface<MeshsetsManager>()->getCubitMeshsetPtr(
285
286 const std::string bc_id =
287 problem_name + "_" + field_name + "_DISPLACEMENTSET" +
288 boost::lexical_cast<std::string>(m->getMeshsetId());
289
290 auto bc = bcMapByBlockName.at(bc_id);
291
292 if (bc->dispBcPtr)
293 if (bc->dispBcPtr->data.flag1) {
294 ents_to_remove[0].merge(bc->bcEnts);
295 }
296 if (bc->dispBcPtr->data.flag2) {
297 ents_to_remove[1].merge(bc->bcEnts);
298 }
299 if (bc->dispBcPtr->data.flag3) {
300 ents_to_remove[2].merge(bc->bcEnts);
301 }
302
303 bc->bcMarkers = std::vector<unsigned char>();
304 }
305
306 auto remove_dofs_on_ents = [&](const Range &ents, const int lo,
307 const int hi) {
308 if (is_distributed_mesh)
309 return prb_mng->removeDofsOnEntities(problem_name, field_name, ents, lo,
310 hi);
311 else
312 return prb_mng->removeDofsOnEntitiesNotDistributed(
313 problem_name, field_name, ents, lo, hi);
314 };
315
316 CHKERR remove_dofs_on_ents(ents_to_remove[0], 0, 0);
317 CHKERR remove_dofs_on_ents(ents_to_remove[1], 1, 1);
318 CHKERR remove_dofs_on_ents(ents_to_remove[2], 2, 2);
319
321}

◆ removeBlockDOFsOnEntities() [13/19]

template<>
MoFEMErrorCode MoFEM::BcManager::removeBlockDOFsOnEntities ( const std::string  problem_name,
const std::string  field_name,
bool  get_low_dim_ents,
bool  block_name_field_prefix,
bool  is_distributed_mesh 
)

Definition at line 324 of file BcManager.cpp.

328 {
329 Interface &m_field = cOre;
330 auto prb_mng = m_field.getInterface<ProblemsManager>();
332
333 CHKERR pushMarkDOFsOnEntities<BcMeshsetType<TEMPERATURESET>>(
334 problem_name, field_name, get_low_dim_ents, block_name_field_prefix);
335
336 Range ents_to_remove;
337
338 for (auto m :
339
340 m_field.getInterface<MeshsetsManager>()->getCubitMeshsetPtr(
342
343 ) {
344 const std::string bc_id =
345 problem_name + "_" + field_name + "_TEMPERATURESET" +
346 boost::lexical_cast<std::string>(m->getMeshsetId());
347 auto bc = bcMapByBlockName.at(bc_id);
348 ents_to_remove.merge(bc->bcEnts);
349 bc->bcMarkers = std::vector<unsigned char>();
350 }
351
352 auto remove_dofs_on_ents = [&](const Range &ents, const int lo,
353 const int hi) {
354 if (is_distributed_mesh)
355 return prb_mng->removeDofsOnEntities(problem_name, field_name, ents, lo,
356 hi);
357 else
358 return prb_mng->removeDofsOnEntitiesNotDistributed(
359 problem_name, field_name, ents, lo, hi);
360 };
361
362 CHKERR remove_dofs_on_ents(ents_to_remove, 0, MAX_DOFS_ON_ENTITY);
363
365}

◆ removeBlockDOFsOnEntities() [14/19]

template<>
MoFEMErrorCode MoFEM::BcManager::removeBlockDOFsOnEntities ( const std::string  problem_name,
const std::string  field_name,
bool  get_low_dim_ents,
bool  block_name_field_prefix,
bool  is_distributed_mesh 
)

Definition at line 368 of file BcManager.cpp.

371 {
372 Interface &m_field = cOre;
373 auto prb_mng = m_field.getInterface<ProblemsManager>();
375
376 CHKERR pushMarkDOFsOnEntities<BcMeshsetType<HEATFLUXSET>>(
377 problem_name, field_name, get_low_dim_ents, block_name_field_prefix);
378
379 Range ents_to_remove;
380
381 for (auto m :
382
383 m_field.getInterface<MeshsetsManager>()->getCubitMeshsetPtr(NODESET |
385
386 ) {
387 const std::string bc_id =
388 problem_name + "_" + field_name + "_HEATFLUXSET" +
389 boost::lexical_cast<std::string>(m->getMeshsetId());
390 auto bc = bcMapByBlockName.at(bc_id);
391 ents_to_remove.merge(bc->bcEnts);
392 bc->bcMarkers = std::vector<unsigned char>();
393 }
394
395 auto remove_dofs_on_ents = [&](const Range &ents, const int lo,
396 const int hi) {
397 if (is_distributed_mesh)
398 return prb_mng->removeDofsOnEntities(problem_name, field_name, ents, lo,
399 hi);
400 else
401 return prb_mng->removeDofsOnEntitiesNotDistributed(
402 problem_name, field_name, ents, lo, hi);
403 };
404
405 CHKERR remove_dofs_on_ents(ents_to_remove, 0, MAX_DOFS_ON_ENTITY);
406
408}

◆ removeBlockDOFsOnEntities() [15/19]

template<>
MoFEMErrorCode MoFEM::BcManager::removeBlockDOFsOnEntities ( const std::string  problem_name,
const std::string  field_name,
bool  get_low_dim_ents,
bool  block_name_field_prefix,
bool  is_distributed_mesh 
)

Definition at line 411 of file BcManager.cpp.

415 {
416 Interface &m_field = cOre;
417 auto prb_mng = m_field.getInterface<ProblemsManager>();
419
420 CHKERR pushMarkDOFsOnEntities<BcVectorMeshsetType<BLOCKSET>>(
421 problem_name, field_name, get_low_dim_ents, block_name_field_prefix);
422
423 std::array<Range, 3> ents_to_remove;
424
425 for (auto m :
426
427 m_field.getInterface<MeshsetsManager>()->getCubitMeshsetPtr(BLOCKSET)) {
428
429 const auto block_name = m->getName();
430
431 std::string bc_id = problem_name + "_" + field_name + "_" + m->getName();
432 std::string regex_str;
433 if (block_name_field_prefix) {
434 regex_str = (boost::format("%s_%s_%s_((FIX_(ALL|X|Y|Z))|("
435 "DISPLACEMENT))(.*)") %
436 problem_name % field_name % field_name)
437 .str();
438 } else {
439 regex_str = (boost::format("%s_%s_((FIX_(ALL|X|Y|Z))|("
440 "DISPLACEMENT))(.*)") %
441 problem_name % field_name)
442 .str();
443 }
444
445 if (std::regex_match(bc_id, std::regex(regex_str))) {
446
447 auto bc = bcMapByBlockName.at(bc_id);
448
449 if (auto disp_bc = bc->dispBcPtr) {
450 if (disp_bc->data.flag1) {
451 ents_to_remove[0].merge(bc->bcEnts);
452 }
453 if (disp_bc->data.flag2) {
454 ents_to_remove[1].merge(bc->bcEnts);
455 }
456 if (disp_bc->data.flag3) {
457 ents_to_remove[2].merge(bc->bcEnts);
458 }
459 } else {
460 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
461 "BC type not implemented");
462 }
463 }
464 }
465
466 auto remove_dofs_on_ents = [&](const Range &ents, const int lo,
467 const int hi) {
468 if (is_distributed_mesh)
469 return prb_mng->removeDofsOnEntities(problem_name, field_name, ents, lo,
470 hi);
471 else
472 return prb_mng->removeDofsOnEntitiesNotDistributed(
473 problem_name, field_name, ents, lo, hi);
474 };
475
476 CHKERR remove_dofs_on_ents(ents_to_remove[0], 0, 0);
477 CHKERR remove_dofs_on_ents(ents_to_remove[1], 1, 1);
478 CHKERR remove_dofs_on_ents(ents_to_remove[2], 2, 2);
479
481}

◆ removeBlockDOFsOnEntities() [16/19]

template<>
MoFEMErrorCode MoFEM::BcManager::removeBlockDOFsOnEntities ( const std::string  problem_name,
const std::string  field_name,
bool  get_low_dim_ents,
bool  block_name_field_prefix,
bool  is_distributed_mesh 
)

Definition at line 965 of file BcManager.cpp.

968 {
970 CHKERR removeBlockDOFsOnEntities<BcMeshsetType<DISPLACEMENTSET>>(
971 problem_name, field_name, get_low_dim_ents, block_name_field_prefix,
972 is_distributed_mesh);
973 CHKERR removeBlockDOFsOnEntities<BcVectorMeshsetType<BLOCKSET>>(
974 problem_name, field_name, get_low_dim_ents, block_name_field_prefix,
975 is_distributed_mesh);
977}

◆ removeBlockDOFsOnEntities() [17/19]

template<>
MoFEMErrorCode MoFEM::BcManager::removeBlockDOFsOnEntities ( const std::string  problem_name,
const std::string  field_name,
bool  get_low_dim_ents,
bool  block_name_field_prefix,
bool  is_distributed_mesh 
)

Definition at line 1000 of file BcManager.cpp.

1003 {
1005 CHKERR removeBlockDOFsOnEntities<BcMeshsetType<TEMPERATURESET>>(
1006 problem_name, field_name, get_low_dim_ents, block_name_field_prefix,
1007 is_distributed_mesh);
1008
1009 auto get_block_name = [&]() {
1010 if (block_name_field_prefix)
1011 return (boost::format("%s_FIX_SCALAR") % field_name).str();
1012 else
1013 return field_name;
1014 };
1015
1016 CHKERR removeBlockDOFsOnEntities<BcScalarMeshsetType<BLOCKSET>>(
1017 problem_name, get_block_name(), field_name, get_low_dim_ents,
1018 is_distributed_mesh);
1020}

◆ removeBlockDOFsOnEntities() [18/19]

template<>
MoFEMErrorCode MoFEM::BcManager::removeBlockDOFsOnEntities ( const std::string  problem_name,
const std::string  field_name,
bool  get_low_dim_ents,
bool  block_name_field_prefix,
bool  is_distributed_mesh 
)

Definition at line 1033 of file BcManager.cpp.

1036 {
1038 CHKERR removeBlockDOFsOnEntities<BcMeshsetType<HEATFLUXSET>>(
1039 problem_name, field_name, get_low_dim_ents, block_name_field_prefix,
1040 is_distributed_mesh);
1042}

◆ removeBlockDOFsOnEntities() [19/19]

template<typename T >
MoFEMErrorCode MoFEM::BcManager::removeBlockDOFsOnEntities ( const std::string  problem_name,
const std::string  field_name,
bool  get_low_dim_ents = true,
bool  block_name_field_prefix = false,
bool  is_distributed_mesh = true 
)

Mark block DOFs.

Template Parameters
BCSET
Parameters
problem_name
field_name
get_low_dim_ents
is_distributed_mesh
block_name_field_prefix
Returns
MoFEMErrorCode

Member Data Documentation

◆ bcMapByBlockName

BcMapByBlockName MoFEM::BcManager::bcMapByBlockName
private

Definition at line 303 of file BcManager.hpp.

◆ cOre

MoFEM::Core& MoFEM::BcManager::cOre
private

Definition at line 301 of file BcManager.hpp.


The documentation for this struct was generated from the following files: