v0.9.1
BCMultiIndices.cpp
Go to the documentation of this file.
1 /** \file BCMultiIndices.cpp
2  * \brief Structures for managing boundary conditions
3  */
4 
5 /* MoFEM is free software: you can redistribute it and/or modify it under
6  * the terms of the GNU Lesser General Public License as published by the
7  * Free Software Foundation, either version 3 of the License, or (at your
8  * option) any later version.
9  *
10  * MoFEM is distributed in the hope that it will be useful, but WITHOUT
11  * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
12  * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
13  * License for more details.
14  *
15  * You should have received a copy of the GNU Lesser General Public
16  * License along with MoFEM. If not, see <http://www.gnu.org/licenses/>
17  */
18 
19 namespace MoFEM {
20 
21 // moab base meshsets
24  CHKERR moab.tag_get_handle(DIRICHLET_SET_TAG_NAME, nsTag);
25  CHKERR moab.tag_get_handle(NEUMANN_SET_TAG_NAME, ssTag);
26  CHKERR moab.tag_get_handle(
27  (std::string(DIRICHLET_SET_TAG_NAME) + "__BC_DATA").c_str(), nsTag_data);
28  CHKERR moab.tag_get_handle(
29  (std::string(NEUMANN_SET_TAG_NAME) + "__BC_DATA").c_str(), ssTag_data);
30  CHKERR moab.tag_get_handle(MATERIAL_SET_TAG_NAME, bhTag);
31  CHKERR moab.tag_get_handle(BLOCK_HEADER, bhTag_header);
32  CHKERR moab.tag_get_handle(BLOCK_ATTRIBUTES, thBlockAttribs);
33  CHKERR moab.tag_get_handle(NAME_TAG_NAME, entityNameTag);
35 }
37  : meshset(_meshset), cubitBcType(UNKNOWNSET), msId(nullptr),
38  tag_bc_data(nullptr), tag_bc_size(0), tag_block_header_data(nullptr),
39  tag_block_attributes(nullptr), tag_block_attributes_size(0),
40  tagName(nullptr), meshsets_mask(NODESET | SIDESET | BLOCKSET) {
41 
42  CHKERR getTagsHandlers(moab);
43  CHKERR moab.tag_get_tags_on_entity(meshset, tag_handles);
44  std::vector<Tag>::iterator tit = tag_handles.begin();
45  for (; tit != tag_handles.end(); tit++) {
46  if (*tit == nsTag || *tit == ssTag || *tit == bhTag) {
47  CHKERR moab.tag_get_by_ptr(*tit, &meshset, 1, (const void **)&msId);
48  }
49  if (*tit == nsTag) {
50  if (*msId != -1) {
52  }
53  }
54  if (*tit == ssTag) {
55  if (*msId != -1) {
57  }
58  }
59  if (*tit == bhTag) {
60  if (*msId != -1) {
62  }
63  }
64  if ((*tit == nsTag_data) || (*tit == ssTag_data)) {
65  CHKERR moab.tag_get_by_ptr(*tit, &meshset, 1, (const void **)&tag_bc_data,
66  &tag_bc_size);
67 
69  }
70  if (*tit == bhTag_header) {
71  int tag_length;
72  CHKERR moab.tag_get_length(*tit, tag_length);
73  CHKERR moab.tag_get_by_ptr(*tit, &meshset, 1,
74  (const void **)&tag_block_header_data);
75  if (tag_block_header_data[1] > 0)
77  }
78  if (*tit == thBlockAttribs) {
79  CHKERR moab.tag_get_by_ptr(*tit, &meshset, 1,
80  (const void **)&tag_block_attributes,
82  }
83  if (*tit == entityNameTag) {
84  CHKERR moab.tag_get_by_ptr(entityNameTag, &meshset, 1,
85  (const void **)&tagName);
87  }
88  }
89 
90  // If BC set has name, unset UNKNOWNNAME
91  if (cubitBcType.to_ulong() &
94  if ((cubitBcType & CubitBCType(UNKNOWNNAME)).any()) {
96  }
97  }
98 }
100  const CubitBCType _cubit_bc_type, const int _ms_id)
101  : cubitBcType(_cubit_bc_type), msId(nullptr), tag_bc_data(nullptr),
102  tag_bc_size(0), tag_block_header_data(nullptr),
103  tag_block_attributes(nullptr), tag_block_attributes_size(0),
104  tagName(nullptr), meshsets_mask(NODESET | SIDESET | BLOCKSET) {
105 
106  CHKERR getTagsHandlers(moab);
107  CHKERR moab.create_meshset(MESHSET_SET | MESHSET_TRACK_OWNER, meshset);
108  switch (cubitBcType.to_ulong()) {
109  case NODESET:
110  CHKERR moab.tag_set_data(nsTag, &meshset, 1, &_ms_id);
111  CHKERR moab.tag_get_by_ptr(nsTag, &meshset, 1, (const void **)&msId);
112  break;
113  case SIDESET:
114  CHKERR moab.tag_set_data(ssTag, &meshset, 1, &_ms_id);
115  CHKERR moab.tag_get_by_ptr(ssTag, &meshset, 1, (const void **)&msId);
116  break;
117  case BLOCKSET:
118  CHKERR moab.tag_set_data(bhTag, &meshset, 1, &_ms_id);
119  CHKERR moab.tag_get_by_ptr(bhTag, &meshset, 1, (const void **)&msId);
120  break;
121  default: {
122  PetscTraceBackErrorHandler(PETSC_COMM_WORLD, __LINE__, PETSC_FUNCTION_NAME,
123  __FILE__, MOFEM_DATA_INCONSISTENCY,
124  PETSC_ERROR_INITIAL, "Unknow meshset type",
125  PETSC_NULL);
126  PetscMPIAbortErrorHandler(PETSC_COMM_WORLD, __LINE__, PETSC_FUNCTION_NAME,
127  __FILE__, MOFEM_DATA_INCONSISTENCY,
128  PETSC_ERROR_INITIAL, "Unknown meshset type",
129  PETSC_NULL);
130  }
131  }
132 }
134  moab::Interface &moab, const int dimension, Range &entities,
135  const bool recursive) const {
137  rval =
138  moab.get_entities_by_dimension(meshset, dimension, entities, recursive);
139  if (rval != MB_SUCCESS) {
140  std::ostringstream ss;
141  ss << "bc set " << *this << std::endl;
142  PetscPrintf(PETSC_COMM_WORLD, ss.str().c_str());
143  }
144  CHKERR rval;
146 }
148  moab::Interface &moab, Range &entities, const bool recursive) const {
150  if ((cubitBcType & CubitBCType(BLOCKSET)).any()) {
151  if (tag_block_header_data != nullptr) {
153  entities, recursive);
154  } else {
155  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "dimension unknown");
156  }
157  }
158  if ((cubitBcType & CubitBCType(NODESET)).any()) {
159  return getMeshsetIdEntitiesByDimension(moab, 0, entities, recursive);
160  }
162 }
164  moab::Interface &moab, const EntityType type, Range &entities,
165  const bool recursive) const {
167  rval = moab.get_entities_by_type(meshset, type, entities, recursive);
168  if (rval != MB_SUCCESS) {
169  std::ostringstream ss;
170  ss << "bc set " << *this << std::endl;
171  PetscPrintf(PETSC_COMM_WORLD, ss.str().c_str());
172  }
173  CHKERR rval;
175 }
176 
177 MoFEMErrorCode CubitMeshSets::getBcData(std::vector<char> &bc_data) const {
179  bc_data.resize(tag_bc_size);
180  copy(&tag_bc_data[0], &tag_bc_data[tag_bc_size], bc_data.begin());
182 }
183 
185  std::vector<unsigned int> &material_data) const {
187  material_data.resize(3);
189  material_data.begin());
191 }
192 
195  if (tag_block_header_data != nullptr) {
196  std::vector<unsigned int> material_data;
197  getBlockHeaderData(material_data);
198  os << "block_header_data = ";
199  std::vector<unsigned int>::iterator vit = material_data.begin();
200  for (; vit != material_data.end(); vit++) {
201  os << std::hex << (int)((unsigned int)*vit) << " ";
202  }
203  os << ": ";
204  vit = material_data.begin();
205  for (; vit != material_data.end(); vit++) {
206  os << *vit;
207  }
208  os << std::endl;
209  } else {
210  os << "no block header data" << std::endl;
211  }
213 }
214 
215 std::string CubitMeshSets::getName() const {
216  if (tagName != nullptr) {
217  return std::string(tagName);
218  } else {
219  return "NoNameSet";
220  }
221 }
222 
223 MoFEMErrorCode CubitMeshSets::printName(std::ostream &os) const {
225  std::string name = getName();
226  os << std::endl;
227  os << "Block name: " << name << std::endl;
229 }
230 
232 CubitMeshSets::getTypeFromBcData(const std::vector<char> &bc_data,
233  CubitBCType &type) const {
235 
236  // See CubitBCType in common.hpp
237  if (bc_data.size() == 0) {
239  }
240 
241  if (strcmp(&bc_data[0], "Displacement") == 0)
243  else if (strcmp(&bc_data[0], "Force") == 0)
244  type |= FORCESET;
245  else if (strcmp(&bc_data[0], "Velocity") == 0)
246  type |= VELOCITYSET;
247  else if (strcmp(&bc_data[0], "Acceleration") == 0)
249  else if (strcmp(&bc_data[0], "Temperature") == 0)
250  type |= TEMPERATURESET;
251  else if (strcmp(&bc_data[0], "Pressure") == 0)
252  type |= PRESSURESET;
253  else if (strcmp(&bc_data[0], "HeatFlux") == 0)
254  type |= HEATFLUXSET;
255  else if (strcmp(&bc_data[0], "cfd_bc") == 0)
256  type |= INTERFACESET;
257  else
258  type |= UNKNOWNNAME;
259 
261 }
264  std::vector<char> bc_data;
265  CHKERR getBcData(bc_data);
266  CHKERR getTypeFromBcData(bc_data, type);
268 }
269 
272  std::vector<char> bc_data;
273  CHKERR getBcData(bc_data);
274  os << "bc_data = ";
275  std::vector<char>::iterator vit = bc_data.begin();
276  for (; vit != bc_data.end(); vit++) {
277  os << std::hex << (int)((unsigned char)*vit) << " ";
278  }
279  os << ": ";
280  vit = bc_data.begin();
281  for (; vit != bc_data.end(); vit++) {
282  os << *vit;
283  }
284  os << std::endl;
286 }
287 
289 CubitMeshSets::getAttributes(std::vector<double> &attributes) const {
291  attributes.resize(tag_block_attributes_size);
292  if (tag_block_attributes_size > 0) {
293  copy(&tag_block_attributes[0],
294  &tag_block_attributes[tag_block_attributes_size], attributes.begin());
295  }
297 }
298 
301  const std::vector<double> &attributes) {
302 
304  int tag_size[] = {(int)attributes.size()};
305  void const *tag_data[] = {&*attributes.begin()};
306  CHKERR moab.tag_set_by_ptr(thBlockAttribs, &meshset, 1, tag_data, tag_size);
307  CHKERR moab.tag_get_by_ptr(thBlockAttribs, &meshset, 1,
308  (const void **)&tag_block_attributes,
311 }
312 
315  std::vector<double> attributes;
316  CHKERR getAttributes(attributes);
317  os << std::endl;
318  os << "Block attributes" << std::endl;
319  os << "----------------" << std::endl;
320  for (unsigned int ii = 0; ii < attributes.size(); ii++) {
321  os << "attr. no: " << ii + 1 << " value: " << attributes[ii] << std::endl;
322  }
323  os << std::endl;
325 }
326 
328  CubitBCType &type) const {
330  // See CubitBCType in common.hpp
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) {
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 
356 std::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.tag_block_header_data != nullptr) {
377  os << " block header: ";
378  os << " blockCol = " << e.tag_block_header_data[0];
379  os << " blockMat = " << e.tag_block_header_data[1];
380  os << " blockDimension = " << e.tag_block_header_data[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 
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 
412  e.cubitBcType |= bIt;
413 }
414 
417 }
418 
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.tag_block_attributes,
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.tag_bc_data, &e.tag_bc_size);
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.tag_bc_data, &e.tag_bc_size);
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
MoFEMErrorCode printBcData(std::ostream &os) const
print bc_data int stream given by os
block name is "MAT_MOISTURE"
Definition: definitions.h:231
MoFEMErrorCode printAttributes(std::ostream &os) const
print the attributes vector
CubitMeshSets(Interface &moab, const EntityHandle _meshset)
MoFEMErrorCode printBlockHeaderData(std::ostream &os) const
print material_data int stream given by os
MoFEMErrorCode getMeshsetIdEntitiesByDimension(Interface &moab, const int dimension, Range &entities, const bool recursive=false) const
get entities form meshset
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
this struct keeps basic methods for moab meshset about material and boundary conditions
virtual std::size_t getSizeOfData() const =0
get data structure size
const std::vector< double > & aTtr
#define THROW_MESSAGE(a)
Throw MoFEM exception.
Definition: definitions.h:625
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:513
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
std::vector< Tag > tag_handles
vector of tag handles to types of data passed from cubit
static MoFEMErrorCodeGeneric< moab::ErrorCode > rval
Definition: Exceptions.hpp:84
std::ostream & operator<<(std::ostream &os, const DataForcesAndSourcesCore::EntData &e)
virtual const void * getDataPtr() const =0
get pointer to data structure
MoFEMErrorCode getBlockHeaderData(std::vector< unsigned int > &material_data) const
get block_headers vector from MoFEM database
MoFEMErrorCode setBcDataStructure(CUBIT_BC_DATA_TYPE &data)
MoFEMErrorCode getTagsHandlers(Interface &moab)
MoFEMErrorCode getTypeFromName(const std::string &name, CubitBCType &type) const
Function that returns the CubitBCType type of the block name, sideset name etc.
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:66
std::string getName() const
get name of block, sideset etc. (this is set in Cubit block properties)
block name is "MAT_THERMAL"
Definition: definitions.h:229
block name is "BODY_FORCES"
Definition: definitions.h:230
MoFEMErrorCode getAttributes(std::vector< double > &attributes) const
get Cubit block attributes
#define CHKERR
Inline error check.
Definition: definitions.h:601
block name is "MAT_ELASTIC"
Definition: definitions.h:227
MoFEMErrorCode printName(std::ostream &os) const
print name of block, sideset etc. (this is set in Cubit setting properties)
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.
std::bitset< 32 > CubitBCType
Definition: Types.hpp:63
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1791
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412
void operator()(CubitMeshSets &e)
virtual std::size_t getSizeOfData() const =0
get data structure size
unsigned int * tag_block_header_data
static const char *const CubitBCNames[]
Names of types of sets and boundary conditions.
Definition: definitions.h:243
int * msId
cubit meshset ID
virtual const void * getDataPtr() const =0
get pointer to data structure
MoFEMErrorCode setAttributes(moab::Interface &moab, const std::vector< double > &attributes)
cet Cubit block attributes
MoFEMErrorCode setAttributeDataStructure(const ATTRIBUTE_TYPE &data)
fill meshset data with data on structure
MoFEMErrorCode getBcData(std::vector< char > &bc_data) const
get bc_data vector from MoFEM database