v0.14.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
5
6namespace MoFEM {
7
8// moab base meshsets
11 CHKERR moab.tag_get_handle(DIRICHLET_SET_TAG_NAME, nsTag);
12 CHKERR moab.tag_get_handle(NEUMANN_SET_TAG_NAME, ssTag);
13 CHKERR moab.tag_get_handle(
14 (std::string(DIRICHLET_SET_TAG_NAME) + "__BC_DATA").c_str(), nsTag_data);
15 CHKERR moab.tag_get_handle(
16 (std::string(NEUMANN_SET_TAG_NAME) + "__BC_DATA").c_str(), ssTag_data);
17 CHKERR moab.tag_get_handle(MATERIAL_SET_TAG_NAME, bhTag);
18 CHKERR moab.tag_get_handle(BLOCK_HEADER, bhTag_header);
19 CHKERR moab.tag_get_handle(BLOCK_ATTRIBUTES, thBlockAttribs);
20 CHKERR moab.tag_get_handle(NAME_TAG_NAME, entityNameTag);
22}
23CubitMeshSets::CubitMeshSets(moab::Interface &moab, const EntityHandle meshset)
24 : meshset(meshset), cubitBcType(UNKNOWNSET), msId(nullptr),
25 tagBcData(nullptr), tagBcSize(0), tagBlockHeaderData(nullptr),
26 tagBlockAttributes(nullptr), tagBlockAttributesSize(0), tagName(nullptr),
27 meshsetsMask(NODESET | SIDESET | BLOCKSET) {
28
30 CHKERR moab.tag_get_tags_on_entity(meshset, tag_handles);
31 for (auto tit = tag_handles.begin(); tit != tag_handles.end(); tit++) {
32 if (*tit == nsTag || *tit == ssTag || *tit == bhTag) {
33 CHKERR moab.tag_get_by_ptr(*tit, &meshset, 1, (const void **)&msId);
34 }
35 if (*tit == nsTag) {
36 if (*msId != -1) {
38 }
39 }
40 if (*tit == ssTag) {
41 if (*msId != -1) {
43 }
44 }
45 if (*tit == bhTag) {
46 if (*msId != -1) {
48 }
49 }
50 if ((*tit == nsTag_data) || (*tit == ssTag_data)) {
51 CHKERR moab.tag_get_by_ptr(*tit, &meshset, 1, (const void **)&tagBcData,
52 &tagBcSize);
53
55 }
56 if (*tit == bhTag_header) {
57 int tag_length;
58 CHKERR moab.tag_get_length(*tit, tag_length);
59 CHKERR moab.tag_get_by_ptr(*tit, &meshset, 1,
60 (const void **)&tagBlockHeaderData);
61 if (tagBlockHeaderData[1] > 0)
63 }
64 if (*tit == thBlockAttribs) {
65 CHKERR moab.tag_get_by_ptr(*tit, &meshset, 1,
66 (const void **)&tagBlockAttributes,
68 }
69 if (*tit == entityNameTag) {
70 CHKERR moab.tag_get_by_ptr(entityNameTag, &meshset, 1,
71 (const void **)&tagName);
73 }
74 }
75
76 // If BC set has name, unset UNKNOWNNAME
77 if (cubitBcType.to_ulong() &
80 if ((cubitBcType & CubitBCType(UNKNOWNNAME)).any()) {
81 cubitBcType = cubitBcType & (~CubitBCType(UNKNOWNNAME));
82 }
83 }
84}
85CubitMeshSets::CubitMeshSets(moab::Interface &moab,
86 const CubitBCType cubit_bc_type, const int ms_id)
87 : cubitBcType(cubit_bc_type), msId(nullptr), tagBcData(nullptr),
88 tagBcSize(0), tagBlockHeaderData(nullptr), tagBlockAttributes(nullptr),
89 tagBlockAttributesSize(0), tagName(nullptr),
90 meshsetsMask(NODESET | SIDESET | BLOCKSET) {
91
93 CHKERR moab.create_meshset(MESHSET_SET, meshset);
94 switch (cubitBcType.to_ulong()) {
95 case NODESET:
96 CHKERR moab.tag_set_data(nsTag, &meshset, 1, &ms_id);
97 CHKERR moab.tag_get_by_ptr(nsTag, &meshset, 1, (const void **)&msId);
98 break;
99 case SIDESET:
100 CHKERR moab.tag_set_data(ssTag, &meshset, 1, &ms_id);
101 CHKERR moab.tag_get_by_ptr(ssTag, &meshset, 1, (const void **)&msId);
102 break;
103 case BLOCKSET:
104 CHKERR moab.tag_set_data(bhTag, &meshset, 1, &ms_id);
105 CHKERR moab.tag_get_by_ptr(bhTag, &meshset, 1, (const void **)&msId);
106 break;
107 default: {
108 PetscTraceBackErrorHandler(PETSC_COMM_WORLD, __LINE__, PETSC_FUNCTION_NAME,
109 __FILE__, MOFEM_DATA_INCONSISTENCY,
110 PETSC_ERROR_INITIAL, "Unknow meshset type",
111 PETSC_NULL);
112 PetscMPIAbortErrorHandler(PETSC_COMM_WORLD, __LINE__, PETSC_FUNCTION_NAME,
113 __FILE__, MOFEM_DATA_INCONSISTENCY,
114 PETSC_ERROR_INITIAL, "Unknown meshset type",
115 PETSC_NULL);
116 }
117 }
118}
120 moab::Interface &moab, const int dimension, Range &entities,
121 const bool recursive) const {
123 rval =
124 moab.get_entities_by_dimension(meshset, dimension, entities, recursive);
125 if (rval != MB_SUCCESS) {
126 std::ostringstream ss;
127 ss << "bc set " << *this << std::endl;
128 PetscPrintf(PETSC_COMM_WORLD, ss.str().c_str());
129 }
130 CHKERR rval;
132}
134 moab::Interface &moab, Range &entities, const bool recursive) const {
136 if ((cubitBcType & CubitBCType(BLOCKSET)).any()) {
137 if (tagBlockHeaderData != nullptr) {
139 entities, recursive);
140 } else {
141 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "dimension unknown");
142 }
143 }
144 if ((cubitBcType & CubitBCType(NODESET)).any()) {
145 return getMeshsetIdEntitiesByDimension(moab, 0, entities, recursive);
146 }
148}
150 moab::Interface &moab, const EntityType type, Range &entities,
151 const bool recursive) const {
153 rval = moab.get_entities_by_type(meshset, type, entities, recursive);
154 if (rval != MB_SUCCESS) {
155 std::ostringstream ss;
156 ss << "bc set " << *this << std::endl;
157 PetscPrintf(PETSC_COMM_WORLD, ss.str().c_str());
158 }
159 CHKERR rval;
161}
162
163MoFEMErrorCode CubitMeshSets::getBcData(std::vector<char> &bc_data) const {
165 bc_data.resize(tagBcSize);
166 copy(&tagBcData[0], &tagBcData[tagBcSize], bc_data.begin());
168}
169
171 std::vector<unsigned int> &material_data) const {
173 material_data.resize(3);
174 copy(&tagBlockHeaderData[0], &tagBlockHeaderData[3], material_data.begin());
176}
177
180 if (tagBlockHeaderData != nullptr) {
181 std::vector<unsigned int> material_data;
182 getBlockHeaderData(material_data);
183 os << "block_header_data = ";
184 std::vector<unsigned int>::iterator vit = material_data.begin();
185 for (; vit != material_data.end(); vit++) {
186 os << std::hex << (int)((unsigned int)*vit) << " ";
187 }
188 os << ": ";
189 vit = material_data.begin();
190 for (; vit != material_data.end(); vit++) {
191 os << *vit;
192 }
193 os << std::endl;
194 } else {
195 os << "no block header data" << std::endl;
196 }
198}
199
200std::string CubitMeshSets::getName() const {
201 if (tagName != nullptr) {
202 return std::string(tagName);
203 } else {
204 return "NoNameSet";
205 }
206}
207
210 std::string name = getName();
211 os << std::endl;
212 os << "Block name: " << name << std::endl;
214}
215
217CubitMeshSets::getTypeFromBcData(const std::vector<char> &bc_data,
218 CubitBCType &type) const {
220
221 // See CubitBCType in common.hpp
222 if (bc_data.size() == 0) {
224 }
225
226 if (strcmp(&bc_data[0], "Displacement") == 0)
227 type |= DISPLACEMENTSET;
228 else if (strcmp(&bc_data[0], "Force") == 0)
229 type |= FORCESET;
230 else if (strcmp(&bc_data[0], "Velocity") == 0)
231 type |= VELOCITYSET;
232 else if (strcmp(&bc_data[0], "Acceleration") == 0)
233 type |= ACCELERATIONSET;
234 else if (strcmp(&bc_data[0], "Temperature") == 0)
235 type |= TEMPERATURESET;
236 else if (strcmp(&bc_data[0], "Pressure") == 0)
237 type |= PRESSURESET;
238 else if (strcmp(&bc_data[0], "HeatFlux") == 0)
239 type |= HEATFLUXSET;
240 else if (strcmp(&bc_data[0], "cfd_bc") == 0)
241 type |= INTERFACESET;
242 else
243 type |= UNKNOWNNAME;
244
246}
249 std::vector<char> bc_data;
250 CHKERR getBcData(bc_data);
251 CHKERR getTypeFromBcData(bc_data, type);
253}
254
257 std::vector<char> bc_data;
258 CHKERR getBcData(bc_data);
259 os << "bc_data = ";
260 std::vector<char>::iterator vit = bc_data.begin();
261 for (; vit != bc_data.end(); vit++) {
262 os << std::hex << (int)((unsigned char)*vit) << " ";
263 }
264 os << ": ";
265 vit = bc_data.begin();
266 for (; vit != bc_data.end(); vit++) {
267 os << *vit;
268 }
269 os << std::endl;
271}
272
274CubitMeshSets::getAttributes(std::vector<double> &attributes) const {
276 attributes.resize(tagBlockAttributesSize);
277 if (tagBlockAttributesSize > 0) {
279 attributes.begin());
280 }
282}
283
285CubitMeshSets::setAttributes(moab::Interface &moab,
286 const std::vector<double> &attributes) {
287
289 int tag_size[] = {(int)attributes.size()};
290 void const *tag_data[] = {&*attributes.begin()};
291 CHKERR moab.tag_set_by_ptr(thBlockAttribs, &meshset, 1, tag_data, tag_size);
292 CHKERR moab.tag_get_by_ptr(thBlockAttribs, &meshset, 1,
293 (const void **)&tagBlockAttributes,
296}
297
300 std::vector<double> attributes;
301 CHKERR getAttributes(attributes);
302 os << std::endl;
303 os << "Block attributes" << std::endl;
304 os << "----------------" << std::endl;
305 for (unsigned int ii = 0; ii < attributes.size(); ii++) {
306 os << "attr. no: " << ii + 1 << " value: " << attributes[ii] << std::endl;
307 }
308 os << std::endl;
310}
311
313 CubitBCType &type) const {
315 // See CubitBCType in common.hpp
316 if (name.compare(0, 11, "MAT_ELASTIC") == 0) {
317 type |= MAT_ELASTICSET;
318 } else if (name.compare(0, 11, "MAT_THERMAL") == 0) {
319 type |= MAT_THERMALSET;
320 } else if (name.compare(0, 12, "MAT_MOISTURE") == 0) {
321 type |= MAT_MOISTURESET;
322 } else if (name.compare(0, 10, "MAT_INTERF") == 0) {
323 type |= MAT_INTERFSET;
324 } else if (name.compare(0, 11, "BODY_FORCES") == 0) {
325 type |= BODYFORCESSET;
326 }
327 // To be extended as appropriate
328 else {
329 type |= UNKNOWNNAME;
330 }
332}
333
336 std::string name = getName();
337 CHKERR getTypeFromName(name, type);
339}
340
341std::ostream &operator<<(std::ostream &os, const CubitMeshSets &e) {
342 // get name of cubit meshset
343 std::ostringstream ss;
344 unsigned jj = 0;
345 while (1 << jj != LASTSET_BC) {
346 const CubitBCType jj_bc_type = 1 << jj;
347 if ((e.cubitBcType & jj_bc_type).any()) {
348 string bc_type_name;
349 ss << " " << string(CubitBCNames[jj + 1]);
350 }
351 ++jj;
352 }
353
354 // push data to stream
355 os << "meshset " << e.meshset << " type" << ss.str();
356 if (e.msId != nullptr)
357 os << " msId " << *(e.msId);
358 if (e.tagName != nullptr) {
359 os << " name " << e.getName();
360 }
361 if (e.tagBlockHeaderData != nullptr) {
362 os << " block header: ";
363 os << " blockCol = " << e.tagBlockHeaderData[0];
364 os << " blockMat = " << e.tagBlockHeaderData[1];
365 os << " blockDimension = " << e.tagBlockHeaderData[2];
366 }
367 return os;
368}
369
371
372 switch (e.cubitBcType.to_ulong()) {
373 case BLOCKSET: {
374 nAme.resize(NAME_TAG_SIZE);
375 CHKERR mOab.tag_set_data(e.entityNameTag, &e.meshset, 1, nAme.c_str());
376 CHKERR mOab.tag_get_by_ptr(e.entityNameTag, &e.meshset, 1,
377 (const void **)&e.tagName);
378
379 CubitBCType type;
380 CHKERR e.getTypeFromName(type);
381 e.cubitBcType |= type;
382 }; break;
383 case NODESET:
384 case SIDESET: {
385 nAme.resize(NAME_TAG_SIZE);
386 CHKERR mOab.tag_set_data(e.entityNameTag, &e.meshset, 1, nAme.c_str());
387 CHKERR mOab.tag_get_by_ptr(e.entityNameTag, &e.meshset, 1,
388 (const void **)&e.tagName);
389 }; break;
390 default:
391 THROW_MESSAGE("not implemented for this CubitBC type");
392 }
393}
394
396 CubitMeshSets &e) {
397 e.cubitBcType |= bIt;
398}
399
402}
403
405 CubitMeshSets &e) {
406 // Need to run this to set tag size in number of doubles, don;t know nothing
407 // about structure
408 int tag_size[] = {(int)(aTtr.getSizeOfData() / sizeof(double))};
409 void const *tag_data[] = {aTtr.getDataPtr()};
410 CHKERR mOab.tag_set_by_ptr(e.thBlockAttribs, &e.meshset, 1, tag_data,
411 tag_size);
412 CHKERR mOab.tag_get_by_ptr(e.thBlockAttribs, &e.meshset, 1,
413 (const void **)&e.tagBlockAttributes,
415 // Here I know about structure
417}
418
420
421 // Need to run this to set tag size, don;t know nothing about structure
422 int tag_size[] = {(int)bcData.getSizeOfData()};
423 void const *tag_data[] = {bcData.getDataPtr()};
424 if ((e.cubitBcType & CubitBCType(NODESET)).any()) {
425 CHKERR mOab.tag_set_by_ptr(e.nsTag_data, &e.meshset, 1, tag_data, tag_size);
426 CHKERR mOab.tag_get_by_ptr(e.nsTag_data, &e.meshset, 1,
427 (const void **)&e.tagBcData, &e.tagBcSize);
428 } else if ((e.cubitBcType & CubitBCType(SIDESET)).any()) {
429 CHKERR mOab.tag_set_by_ptr(e.ssTag_data, &e.meshset, 1, tag_data, tag_size);
430 CHKERR mOab.tag_get_by_ptr(e.ssTag_data, &e.meshset, 1,
431 (const void **)&e.tagBcData, &e.tagBcSize);
432 } else {
433 THROW_MESSAGE("You have to have NODESET or SIDESET to apply BC data on it");
434 }
435 // Here I know about structure
437 // Get Type form BC data
439}
440
441} // namespace MoFEM
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:447
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
@ TEMPERATURESET
Definition: definitions.h:155
@ MATERIALSET
Definition: definitions.h:149
@ PRESSURESET
Definition: definitions.h:152
@ BODYFORCESSET
block name is "BODY_FORCES"
Definition: definitions.h:162
@ MAT_ELASTICSET
block name is "MAT_ELASTIC"
Definition: definitions.h:159
@ ACCELERATIONSET
Definition: definitions.h:154
@ FORCESET
Definition: definitions.h:151
@ HEATFLUXSET
Definition: definitions.h:156
@ NODESET
Definition: definitions.h:146
@ SIDESET
Definition: definitions.h:147
@ UNKNOWNNAME
Definition: definitions.h:158
@ VELOCITYSET
Definition: definitions.h:153
@ MAT_INTERFSET
Definition: definitions.h:160
@ DISPLACEMENTSET
Definition: definitions.h:150
@ LASTSET_BC
Definition: definitions.h:166
@ MAT_THERMALSET
block name is "MAT_THERMAL"
Definition: definitions.h:161
@ MAT_MOISTURESET
block name is "MAT_MOISTURE"
Definition: definitions.h:163
@ UNKNOWNSET
Definition: definitions.h:145
@ BLOCKSET
Definition: definitions.h:148
@ INTERFACESET
Definition: definitions.h:157
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:31
static const char *const CubitBCNames[]
Names of types of sets and boundary conditions.
Definition: definitions.h:175
#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
#define THROW_MESSAGE(msg)
Throw MoFEM exception.
Definition: definitions.h:561
static MoFEMErrorCodeGeneric< moab::ErrorCode > rval
Definition: Exceptions.hpp:74
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
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