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
29 auto construct = [&]() {
31
33 CHKERR moab.tag_get_tags_on_entity(meshset, tag_handles);
34 for (auto tit = tag_handles.begin(); tit != tag_handles.end(); tit++) {
35 if (*tit == nsTag || *tit == ssTag || *tit == bhTag) {
36 CHKERR moab.tag_get_by_ptr(*tit, &meshset, 1, (const void **)&msId);
37 }
38 if (*tit == nsTag) {
39 if (*msId != -1) {
41 }
42 }
43 if (*tit == ssTag) {
44 if (*msId != -1) {
46 }
47 }
48 if (*tit == bhTag) {
49 if (*msId != -1) {
51 }
52 }
53 if ((*tit == nsTag_data) || (*tit == ssTag_data)) {
54 CHKERR moab.tag_get_by_ptr(*tit, &meshset, 1, (const void **)&tagBcData,
55 &tagBcSize);
56
58 }
59 if (*tit == bhTag_header) {
60 int tag_length;
61 CHKERR moab.tag_get_length(*tit, tag_length);
62 CHKERR moab.tag_get_by_ptr(*tit, &meshset, 1,
63 (const void **)&tagBlockHeaderData);
64 if (tagBlockHeaderData[1] > 0)
66 }
67 if (*tit == thBlockAttribs) {
68 CHKERR moab.tag_get_by_ptr(*tit, &meshset, 1,
69 (const void **)&tagBlockAttributes,
71 }
72 if (*tit == entityNameTag) {
73 CHKERR moab.tag_get_by_ptr(entityNameTag, &meshset, 1,
74 (const void **)&tagName);
76 }
77 }
78
79 // If BC set has name, unset UNKNOWNNAME
80 if (cubitBcType.to_ulong() &
83 if ((cubitBcType & CubitBCType(UNKNOWNNAME)).any()) {
84 cubitBcType = cubitBcType & (~CubitBCType(UNKNOWNNAME));
85 }
86 }
87
89 };
90
91 CHK_THROW_MESSAGE(construct(), "construct CubitMeshSets");
92}
93CubitMeshSets::CubitMeshSets(moab::Interface &moab,
94 const CubitBCType cubit_bc_type, const int ms_id)
95 : cubitBcType(cubit_bc_type), msId(nullptr), tagBcData(nullptr),
96 tagBcSize(0), tagBlockHeaderData(nullptr), tagBlockAttributes(nullptr),
97 tagBlockAttributesSize(0), tagName(nullptr),
98 meshsetsMask(NODESET | SIDESET | BLOCKSET) {
99
101 CHKERR moab.create_meshset(MESHSET_SET, meshset);
102 switch (cubitBcType.to_ulong()) {
103 case NODESET:
104 CHKERR moab.tag_set_data(nsTag, &meshset, 1, &ms_id);
105 CHKERR moab.tag_get_by_ptr(nsTag, &meshset, 1, (const void **)&msId);
106 break;
107 case SIDESET:
108 CHKERR moab.tag_set_data(ssTag, &meshset, 1, &ms_id);
109 CHKERR moab.tag_get_by_ptr(ssTag, &meshset, 1, (const void **)&msId);
110 break;
111 case BLOCKSET:
112 CHKERR moab.tag_set_data(bhTag, &meshset, 1, &ms_id);
113 CHKERR moab.tag_get_by_ptr(bhTag, &meshset, 1, (const void **)&msId);
114 break;
115 default: {
116 PetscTraceBackErrorHandler(PETSC_COMM_WORLD, __LINE__, PETSC_FUNCTION_NAME,
117 __FILE__, MOFEM_DATA_INCONSISTENCY,
118 PETSC_ERROR_INITIAL, "Unknow meshset type",
119 PETSC_NULL);
120 PetscMPIAbortErrorHandler(PETSC_COMM_WORLD, __LINE__, PETSC_FUNCTION_NAME,
121 __FILE__, MOFEM_DATA_INCONSISTENCY,
122 PETSC_ERROR_INITIAL, "Unknown meshset type",
123 PETSC_NULL);
124 }
125 }
126}
128 moab::Interface &moab, const int dimension, Range &entities,
129 const bool recursive) const {
131 rval =
132 moab.get_entities_by_dimension(meshset, dimension, entities, recursive);
133 if (rval != MB_SUCCESS) {
134 std::ostringstream ss;
135 ss << "bc set " << *this << std::endl;
136 PetscPrintf(PETSC_COMM_WORLD, ss.str().c_str());
137 }
138 CHKERR rval;
140}
142 moab::Interface &moab, Range &entities, const bool recursive) const {
144 if ((cubitBcType & CubitBCType(BLOCKSET)).any()) {
145 if (tagBlockHeaderData != nullptr) {
147 entities, recursive);
148 } else {
149 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "dimension unknown");
150 }
151 }
152 if ((cubitBcType & CubitBCType(NODESET)).any()) {
153 return getMeshsetIdEntitiesByDimension(moab, 0, entities, recursive);
154 }
156}
158 moab::Interface &moab, const EntityType type, Range &entities,
159 const bool recursive) const {
161 rval = moab.get_entities_by_type(meshset, type, entities, recursive);
162 if (rval != MB_SUCCESS) {
163 std::ostringstream ss;
164 ss << "bc set " << *this << std::endl;
165 PetscPrintf(PETSC_COMM_WORLD, ss.str().c_str());
166 }
167 CHKERR rval;
169}
170
171MoFEMErrorCode CubitMeshSets::getBcData(std::vector<char> &bc_data) const {
173 bc_data.resize(tagBcSize);
174 copy(&tagBcData[0], &tagBcData[tagBcSize], bc_data.begin());
176}
177
179 std::vector<unsigned int> &material_data) const {
181 material_data.resize(3);
182 copy(&tagBlockHeaderData[0], &tagBlockHeaderData[3], material_data.begin());
184}
185
188 if (tagBlockHeaderData != nullptr) {
189 std::vector<unsigned int> material_data;
190 getBlockHeaderData(material_data);
191 os << "block_header_data = ";
192 std::vector<unsigned int>::iterator vit = material_data.begin();
193 for (; vit != material_data.end(); vit++) {
194 os << std::hex << (int)((unsigned int)*vit) << " ";
195 }
196 os << ": ";
197 vit = material_data.begin();
198 for (; vit != material_data.end(); vit++) {
199 os << *vit;
200 }
201 os << std::endl;
202 } else {
203 os << "no block header data" << std::endl;
204 }
206}
207
208std::string CubitMeshSets::getName() const {
209 if (tagName != nullptr) {
210 return std::string(tagName);
211 } else {
212 return "NoNameSet";
213 }
214}
215
218 std::string name = getName();
219 os << std::endl;
220 os << "Block name: " << name << std::endl;
222}
223
225CubitMeshSets::getTypeFromBcData(const std::vector<char> &bc_data,
226 CubitBCType &type) const {
228
229 // See CubitBCType in common.hpp
230 if (bc_data.size() == 0) {
232 }
233
234 if (strcmp(&bc_data[0], "Displacement") == 0)
235 type |= DISPLACEMENTSET;
236 else if (strcmp(&bc_data[0], "Force") == 0)
237 type |= FORCESET;
238 else if (strcmp(&bc_data[0], "Velocity") == 0)
239 type |= VELOCITYSET;
240 else if (strcmp(&bc_data[0], "Acceleration") == 0)
241 type |= ACCELERATIONSET;
242 else if (strcmp(&bc_data[0], "Temperature") == 0)
243 type |= TEMPERATURESET;
244 else if (strcmp(&bc_data[0], "Pressure") == 0)
245 type |= PRESSURESET;
246 else if (strcmp(&bc_data[0], "HeatFlux") == 0)
247 type |= HEATFLUXSET;
248 else if (strcmp(&bc_data[0], "cfd_bc") == 0)
249 type |= INTERFACESET;
250 else
251 type |= UNKNOWNNAME;
252
254}
257 std::vector<char> bc_data;
258 CHKERR getBcData(bc_data);
259 CHKERR getTypeFromBcData(bc_data, type);
261}
262
265 std::vector<char> bc_data;
266 CHKERR getBcData(bc_data);
267 os << "bc_data = ";
268 std::vector<char>::iterator vit = bc_data.begin();
269 for (; vit != bc_data.end(); vit++) {
270 os << std::hex << (int)((unsigned char)*vit) << " ";
271 }
272 os << ": ";
273 vit = bc_data.begin();
274 for (; vit != bc_data.end(); vit++) {
275 os << *vit;
276 }
277 os << std::endl;
279}
280
282CubitMeshSets::getAttributes(std::vector<double> &attributes) const {
284 attributes.resize(tagBlockAttributesSize);
285 if (tagBlockAttributesSize > 0) {
287 attributes.begin());
288 }
290}
291
293CubitMeshSets::setAttributes(moab::Interface &moab,
294 const std::vector<double> &attributes) {
295
297 int tag_size[] = {(int)attributes.size()};
298 void const *tag_data[] = {&*attributes.begin()};
299 CHKERR moab.tag_set_by_ptr(thBlockAttribs, &meshset, 1, tag_data, tag_size);
300 CHKERR moab.tag_get_by_ptr(thBlockAttribs, &meshset, 1,
301 (const void **)&tagBlockAttributes,
304}
305
308 std::vector<double> attributes;
309 CHKERR getAttributes(attributes);
310 os << std::endl;
311 os << "Block attributes" << std::endl;
312 os << "----------------" << std::endl;
313 for (unsigned int ii = 0; ii < attributes.size(); ii++) {
314 os << "attr. no: " << ii + 1 << " value: " << attributes[ii] << std::endl;
315 }
316 os << std::endl;
318}
319
321 CubitBCType &type) const {
323 // See CubitBCType in common.hpp
324 if (name.compare(0, 11, "MAT_ELASTIC") == 0) {
325 type |= MAT_ELASTICSET;
326 } else if (name.compare(0, 11, "MAT_THERMAL") == 0) {
327 type |= MAT_THERMALSET;
328 } else if (name.compare(0, 12, "MAT_MOISTURE") == 0) {
329 type |= MAT_MOISTURESET;
330 } else if (name.compare(0, 10, "MAT_INTERF") == 0) {
331 type |= MAT_INTERFSET;
332 } else if (name.compare(0, 11, "BODY_FORCES") == 0) {
333 type |= BODYFORCESSET;
334 }
335 // To be extended as appropriate
336 else {
337 type |= UNKNOWNNAME;
338 }
340}
341
344 std::string name = getName();
345 CHKERR getTypeFromName(name, type);
347}
348
349std::ostream &operator<<(std::ostream &os, const CubitMeshSets &e) {
350 // get name of cubit meshset
351 std::ostringstream ss;
352 unsigned jj = 0;
353 while (1 << jj != LASTSET_BC) {
354 const CubitBCType jj_bc_type = 1 << jj;
355 if ((e.cubitBcType & jj_bc_type).any()) {
356 string bc_type_name;
357 ss << " " << string(CubitBCNames[jj + 1]);
358 }
359 ++jj;
360 }
361
362 // push data to stream
363 os << "meshset " << e.meshset << " type" << ss.str();
364 if (e.msId != nullptr)
365 os << " msId " << *(e.msId);
366 if (e.tagName != nullptr) {
367 os << " name " << e.getName();
368 }
369 if (e.tagBlockHeaderData != nullptr) {
370 os << " block header: ";
371 os << " blockCol = " << e.tagBlockHeaderData[0];
372 os << " blockMat = " << e.tagBlockHeaderData[1];
373 os << " blockDimension = " << e.tagBlockHeaderData[2];
374 }
375 return os;
376}
377
379
380 switch (e.cubitBcType.to_ulong()) {
381 case BLOCKSET: {
382 nAme.resize(NAME_TAG_SIZE);
383 CHKERR mOab.tag_set_data(e.entityNameTag, &e.meshset, 1, nAme.c_str());
384 CHKERR mOab.tag_get_by_ptr(e.entityNameTag, &e.meshset, 1,
385 (const void **)&e.tagName);
386
387 CubitBCType type;
388 CHKERR e.getTypeFromName(type);
389 e.cubitBcType |= type;
390 }; break;
391 case NODESET:
392 case SIDESET: {
393 nAme.resize(NAME_TAG_SIZE);
394 CHKERR mOab.tag_set_data(e.entityNameTag, &e.meshset, 1, nAme.c_str());
395 CHKERR mOab.tag_get_by_ptr(e.entityNameTag, &e.meshset, 1,
396 (const void **)&e.tagName);
397 }; break;
398 default:
399 THROW_MESSAGE("not implemented for this CubitBC type");
400 }
401}
402
407
411
413 CubitMeshSets &e) {
414 // Need to run this to set tag size in number of doubles, don;t know nothing
415 // about structure
416 int tag_size[] = {(int)(aTtr.getSizeOfData() / sizeof(double))};
417 void const *tag_data[] = {aTtr.getDataPtr()};
418 CHKERR mOab.tag_set_by_ptr(e.thBlockAttribs, &e.meshset, 1, tag_data,
419 tag_size);
420 CHKERR mOab.tag_get_by_ptr(e.thBlockAttribs, &e.meshset, 1,
421 (const void **)&e.tagBlockAttributes,
423 // Here I know about structure
425}
426
428
429 // Need to run this to set tag size, don;t know nothing about structure
430 int tag_size[] = {(int)bcData.getSizeOfData()};
431 void const *tag_data[] = {bcData.getDataPtr()};
432 if ((e.cubitBcType & CubitBCType(NODESET)).any()) {
433 CHKERR mOab.tag_set_by_ptr(e.nsTag_data, &e.meshset, 1, tag_data, tag_size);
434 CHKERR mOab.tag_get_by_ptr(e.nsTag_data, &e.meshset, 1,
435 (const void **)&e.tagBcData, &e.tagBcSize);
436 } else if ((e.cubitBcType & CubitBCType(SIDESET)).any()) {
437 CHKERR mOab.tag_set_by_ptr(e.ssTag_data, &e.meshset, 1, tag_data, tag_size);
438 CHKERR mOab.tag_get_by_ptr(e.ssTag_data, &e.meshset, 1,
439 (const void **)&e.tagBcData, &e.tagBcSize);
440 } else {
441 THROW_MESSAGE("You have to have NODESET or SIDESET to apply BC data on it");
442 }
443 // Here I know about structure
445 // Get Type form BC data
447}
448
449} // 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.
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