v0.15.0
Loading...
Searching...
No Matches
BCMultiIndices.cpp
Go to the documentation of this file.
1/** \file BCMultiIndices.cpp
2 * \brief Structures for managing boundary conditions
3 */
4
5namespace MoFEM {
6
7// moab base meshsets
10 CHKERR moab.tag_get_handle(DIRICHLET_SET_TAG_NAME, nsTag);
11 CHKERR moab.tag_get_handle(NEUMANN_SET_TAG_NAME, ssTag);
12 CHKERR moab.tag_get_handle(
13 (std::string(DIRICHLET_SET_TAG_NAME) + "__BC_DATA").c_str(), nsTag_data);
14 CHKERR moab.tag_get_handle(
15 (std::string(NEUMANN_SET_TAG_NAME) + "__BC_DATA").c_str(), ssTag_data);
16 CHKERR moab.tag_get_handle(MATERIAL_SET_TAG_NAME, bhTag);
17 CHKERR moab.tag_get_handle(BLOCK_HEADER, bhTag_header);
18 CHKERR moab.tag_get_handle(BLOCK_ATTRIBUTES, thBlockAttribs);
19 CHKERR moab.tag_get_handle(NAME_TAG_NAME, entityNameTag);
21}
22CubitMeshSets::CubitMeshSets(moab::Interface &moab, const EntityHandle meshset)
23 : meshset(meshset), cubitBcType(UNKNOWNSET), msId(nullptr),
24 tagBcData(nullptr), tagBcSize(0), tagBlockHeaderData(nullptr),
25 tagBlockAttributes(nullptr), tagBlockAttributesSize(0), tagName(nullptr),
26 meshsetsMask(NODESET | SIDESET | BLOCKSET) {
27
28 auto construct = [&]() {
30
32 CHKERR moab.tag_get_tags_on_entity(meshset, tag_handles);
33 for (auto tit = tag_handles.begin(); tit != tag_handles.end(); tit++) {
34 if (*tit == nsTag || *tit == ssTag || *tit == bhTag) {
35 CHKERR moab.tag_get_by_ptr(*tit, &meshset, 1, (const void **)&msId);
36 }
37 if (*tit == nsTag) {
38 if (*msId != -1) {
40 }
41 }
42 if (*tit == ssTag) {
43 if (*msId != -1) {
45 }
46 }
47 if (*tit == bhTag) {
48 if (*msId != -1) {
50 }
51 }
52 }
53
54 for (auto tit = tag_handles.begin(); tit != tag_handles.end(); tit++) {
55 if ((cubitBcType & CubitBCType(NODESET | SIDESET)).any()) {
56 if ((*tit == nsTag_data) || (*tit == ssTag_data)) {
57 CHKERR moab.tag_get_by_ptr(*tit, &meshset, 1,
58 (const void **)&tagBcData, &tagBcSize);
60 }
61 }
62 if ((cubitBcType & CubitBCType(BLOCKSET)).any()) {
63 if (*tit == bhTag_header) {
64 int tag_length;
65 CHKERR moab.tag_get_length(*tit, tag_length);
66 CHKERR moab.tag_get_by_ptr(*tit, &meshset, 1,
67 (const void **)&tagBlockHeaderData);
68 if (tagBlockHeaderData[1] > 0)
70 }
71 if (*tit == thBlockAttribs) {
72 CHKERR moab.tag_get_by_ptr(*tit, &meshset, 1,
73 (const void **)&tagBlockAttributes,
75 }
76 if (*tit == entityNameTag) {
77 CHKERR moab.tag_get_by_ptr(entityNameTag, &meshset, 1,
78 (const void **)&tagName);
80 }
81 }
82 }
83
84 // If BC set has name, unset UNKNOWNNAME
85 if (cubitBcType.to_ulong() &
88 if ((cubitBcType & CubitBCType(UNKNOWNNAME)).any()) {
89 cubitBcType = cubitBcType & (~CubitBCType(UNKNOWNNAME));
90 }
91 }
92
94 };
95
96 CHK_THROW_MESSAGE(construct(), "construct CubitMeshSets");
97}
98CubitMeshSets::CubitMeshSets(moab::Interface &moab,
99 const CubitBCType cubit_bc_type, const int ms_id)
100 : cubitBcType(cubit_bc_type), msId(nullptr), tagBcData(nullptr),
101 tagBcSize(0), tagBlockHeaderData(nullptr), tagBlockAttributes(nullptr),
102 tagBlockAttributesSize(0), tagName(nullptr),
103 meshsetsMask(NODESET | SIDESET | BLOCKSET) {
104
106 CHKERR moab.create_meshset(MESHSET_SET, meshset);
107 switch (cubitBcType.to_ulong()) {
108 case NODESET:
109 CHKERR moab.tag_set_data(nsTag, &meshset, 1, &ms_id);
110 CHKERR moab.tag_get_by_ptr(nsTag, &meshset, 1, (const void **)&msId);
111 break;
112 case SIDESET:
113 CHKERR moab.tag_set_data(ssTag, &meshset, 1, &ms_id);
114 CHKERR moab.tag_get_by_ptr(ssTag, &meshset, 1, (const void **)&msId);
115 break;
116 case BLOCKSET:
117 CHKERR moab.tag_set_data(bhTag, &meshset, 1, &ms_id);
118 CHKERR moab.tag_get_by_ptr(bhTag, &meshset, 1, (const void **)&msId);
119 break;
120 default: {
121 PetscTraceBackErrorHandler(PETSC_COMM_WORLD, __LINE__, PETSC_FUNCTION_NAME,
122 __FILE__, MOFEM_DATA_INCONSISTENCY,
123 PETSC_ERROR_INITIAL, "Unknow meshset type",
124 PETSC_NULLPTR);
125 PetscMPIAbortErrorHandler(PETSC_COMM_WORLD, __LINE__, PETSC_FUNCTION_NAME,
126 __FILE__, MOFEM_DATA_INCONSISTENCY,
127 PETSC_ERROR_INITIAL, "Unknown meshset type",
128 PETSC_NULLPTR);
129 }
130 }
131}
133 moab::Interface &moab, const int dimension, Range &entities,
134 const bool recursive) const {
136 rval =
137 moab.get_entities_by_dimension(meshset, dimension, entities, recursive);
138 if (rval != MB_SUCCESS) {
139 MOFEM_LOG("SELF", Sev::error) << "bc set " << *this;
140 }
141 CHKERR rval;
143}
145 moab::Interface &moab, Range &entities, const bool recursive) const {
147 if ((cubitBcType & CubitBCType(BLOCKSET)).any()) {
148 if (tagBlockHeaderData != nullptr) {
150 entities, recursive);
151 } else {
152 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "dimension unknown");
153 }
154 }
155 if ((cubitBcType & CubitBCType(NODESET)).any()) {
156 return getMeshsetIdEntitiesByDimension(moab, 0, entities, recursive);
157 }
159}
161 moab::Interface &moab, const EntityType type, Range &entities,
162 const bool recursive) const {
164 rval = moab.get_entities_by_type(meshset, type, entities, recursive);
165 if (rval != MB_SUCCESS) {
166 MOFEM_LOG("SELF", Sev::error) << "bc set " << *this;
167 }
168 CHKERR rval;
170}
171
172MoFEMErrorCode CubitMeshSets::getBcData(std::vector<char> &bc_data) const {
174 bc_data.resize(tagBcSize);
175 copy(&tagBcData[0], &tagBcData[tagBcSize], bc_data.begin());
177}
178
180 std::vector<unsigned int> &material_data) const {
182 material_data.resize(3);
183 copy(&tagBlockHeaderData[0], &tagBlockHeaderData[3], material_data.begin());
185}
186
189 if (tagBlockHeaderData != nullptr) {
190 std::vector<unsigned int> material_data;
191 getBlockHeaderData(material_data);
192 os << "block_header_data = ";
193 std::vector<unsigned int>::iterator vit = material_data.begin();
194 for (; vit != material_data.end(); vit++) {
195 os << std::hex << (int)((unsigned int)*vit) << " ";
196 }
197 os << ": ";
198 vit = material_data.begin();
199 for (; vit != material_data.end(); vit++) {
200 os << *vit;
201 }
202 os << std::endl;
203 } else {
204 os << "no block header data" << std::endl;
205 }
207}
208
209std::string CubitMeshSets::getName() const {
210 if (tagName != nullptr) {
211 return std::string(tagName);
212 } else {
213 return "NoNameSet";
214 }
215}
216
219 std::string name = getName();
220 os << std::endl;
221 os << "Block name: " << name << std::endl;
223}
224
226CubitMeshSets::getTypeFromBcData(const std::vector<char> &bc_data,
227 CubitBCType &type) const {
229
230 // See CubitBCType in common.hpp
231 if (bc_data.size() == 0) {
233 }
234
238 if (strcmp(&bc_data[0], "Displacement") == 0)
239 type |= DISPLACEMENTSET;
240 else if (strcmp(&bc_data[0], "Force") == 0)
241 type |= FORCESET;
242 else if (strcmp(&bc_data[0], "Velocity") == 0)
243 type |= VELOCITYSET;
244 else if (strcmp(&bc_data[0], "Acceleration") == 0)
245 type |= ACCELERATIONSET;
246 else if (strcmp(&bc_data[0], "Temperature") == 0)
247 type |= TEMPERATURESET;
248 else if (strcmp(&bc_data[0], "Pressure") == 0)
249 type |= PRESSURESET;
250 else if (strcmp(&bc_data[0], "HeatFlux") == 0)
251 type |= HEATFLUXSET;
252 else if (strcmp(&bc_data[0], "cfd_bc") == 0)
253 type |= INTERFACESET;
254 else
255 type |= UNKNOWNNAME;
256
258}
261 std::vector<char> bc_data;
262 CHKERR getBcData(bc_data);
263 CHKERR getTypeFromBcData(bc_data, type);
265}
266
269 std::vector<char> bc_data;
270 CHKERR getBcData(bc_data);
271 os << "bc_data = ";
272 std::vector<char>::iterator vit = bc_data.begin();
273 for (; vit != bc_data.end(); vit++) {
274 os << std::hex << (int)((unsigned char)*vit) << " ";
275 }
276 os << ": ";
277 vit = bc_data.begin();
278 for (; vit != bc_data.end(); vit++) {
279 os << *vit;
280 }
281 os << std::endl;
283}
284
286CubitMeshSets::getAttributes(std::vector<double> &attributes) const {
288 attributes.resize(tagBlockAttributesSize);
289 if (tagBlockAttributesSize > 0) {
291 attributes.begin());
292 }
294}
295
297CubitMeshSets::setAttributes(moab::Interface &moab,
298 const std::vector<double> &attributes) {
299
301 int tag_size[] = {(int)attributes.size()};
302 void const *tag_data[] = {&*attributes.begin()};
303 CHKERR moab.tag_set_by_ptr(thBlockAttribs, &meshset, 1, tag_data, tag_size);
304 CHKERR moab.tag_get_by_ptr(thBlockAttribs, &meshset, 1,
305 (const void **)&tagBlockAttributes,
308}
309
312 std::vector<double> attributes;
313 CHKERR getAttributes(attributes);
314 os << std::endl;
315 os << "Block attributes" << std::endl;
316 os << "----------------" << std::endl;
317 for (unsigned int ii = 0; ii < attributes.size(); ii++) {
318 os << "attr. no: " << ii + 1 << " value: " << attributes[ii] << std::endl;
319 }
320 os << std::endl;
322}
323
325 CubitBCType &type) const {
327 // See CubitBCType in common.hpp
330
331 if (name.compare(0, 11, "MAT_ELASTIC") == 0) {
332 type |= MAT_ELASTICSET;
333 } else if (name.compare(0, 11, "MAT_THERMAL") == 0) {
334 type |= MAT_THERMALSET;
335 } else if (name.compare(0, 12, "MAT_MOISTURE") == 0) {
336 type |= MAT_MOISTURESET;
337 } else if (name.compare(0, 10, "MAT_INTERF") == 0) {
338 type |= MAT_INTERFSET;
339 } else if (name.compare(0, 11, "BODY_FORCES") == 0) {
340 type |= BODYFORCESSET;
341 }
342 // To be extended as appropriate
343 else {
344 type |= UNKNOWNNAME;
345 }
347}
348
351 std::string name = getName();
352 CHKERR getTypeFromName(name, type);
354}
355
356std::ostream &operator<<(std::ostream &os, const CubitMeshSets &e) {
357 // get name of cubit meshset
358 std::ostringstream ss;
359 unsigned jj = 0;
360 while (1 << jj != LASTSET_BC) {
361 const CubitBCType jj_bc_type = 1 << jj;
362 if ((e.cubitBcType & jj_bc_type).any()) {
363 string bc_type_name;
364 ss << " " << string(CubitBCNames[jj + 1]);
365 }
366 ++jj;
367 }
368
369 // push data to stream
370 os << "meshset " << e.meshset << " type" << ss.str();
371 if (e.msId != nullptr)
372 os << " msId " << *(e.msId);
373 if (e.tagName != nullptr) {
374 os << " name " << e.getName();
375 }
376 if (e.tagBlockHeaderData != nullptr) {
377 os << " block header: ";
378 os << " blockCol = " << e.tagBlockHeaderData[0];
379 os << " blockMat = " << e.tagBlockHeaderData[1];
380 os << " blockDimension = " << e.tagBlockHeaderData[2];
381 }
382 return os;
383}
384
386
387 switch (e.cubitBcType.to_ulong()) {
388 case BLOCKSET: {
389 nAme.resize(NAME_TAG_SIZE);
390 CHKERR mOab.tag_set_data(e.entityNameTag, &e.meshset, 1, nAme.c_str());
391 CHKERR mOab.tag_get_by_ptr(e.entityNameTag, &e.meshset, 1,
392 (const void **)&e.tagName);
393
394 CubitBCType type;
395 CHKERR e.getTypeFromName(type);
396 e.cubitBcType |= type;
397 }; break;
398 case NODESET:
399 case SIDESET: {
400 nAme.resize(NAME_TAG_SIZE);
401 CHKERR mOab.tag_set_data(e.entityNameTag, &e.meshset, 1, nAme.c_str());
402 CHKERR mOab.tag_get_by_ptr(e.entityNameTag, &e.meshset, 1,
403 (const void **)&e.tagName);
404 }; break;
405 default:
406 THROW_MESSAGE("not implemented for this CubitBC type");
407 }
408}
409
414
418
420 CubitMeshSets &e) {
421 // Need to run this to set tag size in number of doubles, don;t know nothing
422 // about structure
423 int tag_size[] = {(int)(aTtr.getSizeOfData() / sizeof(double))};
424 void const *tag_data[] = {aTtr.getDataPtr()};
425 CHKERR mOab.tag_set_by_ptr(e.thBlockAttribs, &e.meshset, 1, tag_data,
426 tag_size);
427 CHKERR mOab.tag_get_by_ptr(e.thBlockAttribs, &e.meshset, 1,
428 (const void **)&e.tagBlockAttributes,
430 // Here I know about structure
432}
433
435
436 // Need to run this to set tag size, don;t know nothing about structure
437 int tag_size[] = {(int)bcData.getSizeOfData()};
438 void const *tag_data[] = {bcData.getDataPtr()};
439 if ((e.cubitBcType & CubitBCType(NODESET)).any()) {
440 CHKERR mOab.tag_set_by_ptr(e.nsTag_data, &e.meshset, 1, tag_data, tag_size);
441 CHKERR mOab.tag_get_by_ptr(e.nsTag_data, &e.meshset, 1,
442 (const void **)&e.tagBcData, &e.tagBcSize);
443 } else if ((e.cubitBcType & CubitBCType(SIDESET)).any()) {
444 CHKERR mOab.tag_set_by_ptr(e.ssTag_data, &e.meshset, 1, tag_data, tag_size);
445 CHKERR mOab.tag_get_by_ptr(e.ssTag_data, &e.meshset, 1,
446 (const void **)&e.tagBcData, &e.tagBcSize);
447 } else {
448 THROW_MESSAGE("You have to have NODESET or SIDESET to apply BC data on it");
449 }
450 // Here I know about structure
452 // Get Type form BC data
454}
455
456} // namespace MoFEM
#define CHK_THROW_MESSAGE(err, msg)
Check and throw MoFEM exception.
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
@ TEMPERATURESET
@ MATERIALSET
@ PRESSURESET
@ BODYFORCESSET
block name is "BODY_FORCES"
@ MAT_ELASTICSET
block name is "MAT_ELASTIC"
@ ACCELERATIONSET
@ FORCESET
@ HEATFLUXSET
@ NODESET
@ SIDESET
@ UNKNOWNNAME
@ VELOCITYSET
@ MAT_INTERFSET
@ DISPLACEMENTSET
@ LASTSET_BC
@ MAT_THERMALSET
block name is "MAT_THERMAL"
@ MAT_MOISTURESET
block name is "MAT_MOISTURE"
@ UNKNOWNSET
@ BLOCKSET
@ INTERFACESET
@ MOFEM_DATA_INCONSISTENCY
Definition definitions.h:31
static const char *const CubitBCNames[]
Names of types of sets and boundary conditions.
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define CHKERR
Inline error check.
#define THROW_MESSAGE(msg)
Throw MoFEM exception.
#define MOFEM_LOG(channel, severity)
Log.
static MoFEMErrorCodeGeneric< moab::ErrorCode > rval
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
std::bitset< 32 > CubitBCType
Definition Types.hpp:52
implementation of Data Operators for Forces and Sources
Definition Common.hpp:10
std::ostream & operator<<(std::ostream &os, const EntitiesFieldData::EntData &e)
const std::vector< double > & aTtr
void operator()(CubitMeshSets &e)
this struct keeps basic methods for moab meshset about material and boundary conditions
MoFEMErrorCode getAttributes(std::vector< double > &attributes) const
get Cubit block attributes
MoFEMErrorCode getBcData(std::vector< char > &bc_data) const
get bc_data vector from MoFEM database
MoFEMErrorCode printName(std::ostream &os) const
print name of block, sideset etc. (this is set in Cubit setting properties)
MoFEMErrorCode printAttributes(std::ostream &os) const
print the attributes vector
unsigned int * tagBlockHeaderData
MoFEMErrorCode getTypeFromName(const std::string &name, CubitBCType &type) const
Function that returns the CubitBCType type of the block name, sideset name etc.
MoFEMErrorCode setBcDataStructure(CUBIT_BC_DATA_TYPE &data)
CubitMeshSets(Interface &moab, const EntityHandle meshset)
MoFEMErrorCode setAttributes(moab::Interface &moab, const std::vector< double > &attributes)
cet Cubit block attributes
MoFEMErrorCode getBlockHeaderData(std::vector< unsigned int > &material_data) const
get block_headers vector from MoFEM database
MoFEMErrorCode printBcData(std::ostream &os) const
print bc_data int stream given by os
MoFEMErrorCode getMeshsetIdEntitiesByDimension(Interface &moab, const int dimension, Range &entities, const bool recursive=false) const
get entities form meshset
std::string getName() const
get name of block, sideset etc. (this is set in Cubit block properties)
std::vector< Tag > tag_handles
vector of tag handles to types of data passed from cubit
MoFEMErrorCode setAttributeDataStructure(const ATTRIBUTE_TYPE &data)
fill meshset data with data on structure
MoFEMErrorCode getTagsHandlers(Interface &moab)
int * msId
cubit meshset ID
MoFEMErrorCode getMeshsetIdEntitiesByType(Interface &moab, const EntityType type, Range &entities, const bool recursive=false) const
get entities by type
MoFEMErrorCode getTypeFromBcData(const std::vector< char > &bc_data, CubitBCType &type) const
Function that returns the CubitBCType type of the contents of bc_data.
MoFEMErrorCode printBlockHeaderData(std::ostream &os) const
print material_data int stream given by os
virtual const void * getDataPtr() const =0
get pointer to data structure
virtual std::size_t getSizeOfData() const =0
get data structure size
virtual const void * getDataPtr() const =0
get pointer to data structure
virtual std::size_t getSizeOfData() const =0
get data structure size