v0.14.0
BcManager.cpp
Go to the documentation of this file.
1 /** \file BcManager.cpp
2  * \brief Manages boundary conditions
3  * \ingroup bc_manager
4  */
5 
6 namespace MoFEM {
7 
8 namespace BcManagerImplTools {
9 
10 auto get_dim(const Range &ents) {
11  for (auto d : {3, 2, 1})
12  if (ents.num_of_dimension(d))
13  return d;
14  return 0;
15 };
16 
17 auto get_adj_ents(moab::Interface &moab, const Range &ents) {
18  Range verts;
19  CHK_MOAB_THROW(moab.get_connectivity(ents, verts, true), "get verts");
20  const auto dim = get_dim(ents);
21  for (size_t d = 1; d < dim; ++d) {
22  for (auto dd = d + 1; dd <= dim; ++dd) {
23  CHK_MOAB_THROW(moab.get_adjacencies(ents.subset_by_dimension(dd), d,
24  false, verts, moab::Interface::UNION),
25  "get adj");
26  }
27  }
28  verts.merge(ents);
29  return verts;
30 }
31 } // namespace BcManagerImplTools
32 
33 MoFEMErrorCode
34 BcManager::query_interface(boost::typeindex::type_index type_index,
35  UnknownInterface **iface) const {
36  MoFEMFunctionBeginHot;
37  *iface = const_cast<BcManager *>(this);
38  MoFEMFunctionReturnHot(0);
39 }
40 
41 BcManager::BcManager(const Core &core) : cOre(const_cast<Core &>(core)) {
42 
43  if (!LogManager::checkIfChannelExist("BcMngWorld")) {
44  auto core_log = logging::core::get();
45 
46  core_log->add_sink(
47  LogManager::createSink(LogManager::getStrmWorld(), "BcMngWorld"));
48  core_log->add_sink(
49  LogManager::createSink(LogManager::getStrmSync(), "BcMngSync"));
50  core_log->add_sink(
51  LogManager::createSink(LogManager::getStrmSelf(), "BcMngSelf"));
52 
53  LogManager::setLog("BcMngWorld");
54  LogManager::setLog("BcMngSync");
55  LogManager::setLog("BcMngSelf");
56 
57  MOFEM_LOG_TAG("BcMngWorld", "BcMng");
58  MOFEM_LOG_TAG("BcMngSync", "BcMng");
59  MOFEM_LOG_TAG("BcMngSelf", "BcMng");
60  }
61 
62  MOFEM_LOG("BcMngWorld", Sev::noisy) << "BC manager created";
63 }
64 
65 MoFEMErrorCode BcManager::getOptions() {
66  MoFEMFunctionBeginHot;
67  ierr = PetscOptionsBegin(PETSC_COMM_WORLD, "", "BcManager options", "none");
68  ierr = PetscOptionsEnd();
69  CHKERRG(ierr);
70  MoFEMFunctionReturnHot(0);
71 }
72 
73 MoFEMErrorCode BcManager::removeBlockDOFsOnEntities(
74  const std::string problem_name, const std::string block_name,
75  const std::string field_name, int lo, int hi, bool get_low_dim_ents,
76  bool is_distributed_mesh) {
77  Interface &m_field = cOre;
78  auto prb_mng = m_field.getInterface<ProblemsManager>();
79  MoFEMFunctionBegin;
80 
81  CHKERR pushMarkDOFsOnEntities(problem_name, block_name, field_name, lo, hi,
82  get_low_dim_ents);
83 
84  auto remove_dofs_on_ents = [&](const Range &ents, const int lo,
85  const int hi) {
86  if (is_distributed_mesh)
87  return prb_mng->removeDofsOnEntities(problem_name, field_name, ents, lo,
88  hi);
89  else
90  return prb_mng->removeDofsOnEntitiesNotDistributed(
91  problem_name, field_name, ents, lo, hi);
92  };
93 
94  for (auto m :
95  m_field.getInterface<MeshsetsManager>()->getCubitMeshsetPtr(std::regex(
96 
97  (boost::format("%s(.*)") % block_name).str()
98 
99  ))
100 
101  ) {
102  const std::string bc_id =
103  problem_name + "_" + field_name + "_" + m->getName();
104  CHKERR remove_dofs_on_ents(bcMapByBlockName.at(bc_id)->bcEnts, lo, hi);
105  bcMapByBlockName.at(bc_id)->bcMarkers = std::vector<unsigned char>();
106  }
107 
108  MoFEMFunctionReturn(0);
109 }
110 
111 MoFEMErrorCode BcManager::pushMarkDOFsOnEntities(const std::string problem_name,
112  const std::string block_name,
113  const std::string field_name,
114  int lo, int hi,
115  bool get_low_dim_ents) {
116  Interface &m_field = cOre;
117  auto prb_mng = m_field.getInterface<ProblemsManager>();
118  MoFEMFunctionBegin;
119 
120  // if(problem_name.size())
121  // MOFEM_LOG("BcMngWorld", Sev::warning)
122  // << "Argument problem_name has no effect";
123 
124  auto fix_disp = [&]() {
125  MoFEMFunctionBegin;
126 
127  auto mark_fix_dofs = [&](std::vector<unsigned char> &marked_field_dofs,
128  const auto lo, const auto hi) {
129  return prb_mng->modifyMarkDofs(problem_name, ROW, field_name, lo, hi,
130  ProblemsManager::MarkOP::OR, 1,
131  marked_field_dofs);
132  };
133 
134  auto iterate_meshsets = [&](auto &&meshset_vec_ptr) {
135  MoFEMFunctionBegin;
136  for (auto m : meshset_vec_ptr) {
137  auto bc = boost::make_shared<BCs>();
138  CHKERR m_field.get_moab().get_entities_by_handle(m->getMeshset(),
139  bc->bcEnts, true);
140  CHKERR m->getAttributes(bc->bcAttributes);
141 
142  if (problem_name.size())
143  CHKERR mark_fix_dofs(bc->bcMarkers, lo, hi);
144  MOFEM_LOG("BcMngWorld", Sev::verbose)
145  << "Found block " << m->getName() << " number of attributes "
146  << bc->bcAttributes.size();
147 
148  if (get_low_dim_ents) {
149  auto low_dim_ents =
150  BcManagerImplTools::get_adj_ents(m_field.get_moab(), bc->bcEnts);
151  bc->bcEnts.swap(low_dim_ents);
152  }
153 
154  CHKERR m_field.getInterface<CommInterface>()->synchroniseEntities(
155  bc->bcEnts);
156  if (problem_name.size())
157  CHKERR prb_mng->markDofs(problem_name, ROW, ProblemsManager::AND,
158  bc->bcEnts, bc->bcMarkers);
159 
160  const std::string bc_id =
161  problem_name + "_" + field_name + "_" + m->getName();
162  bcMapByBlockName[bc_id] = bc;
163  }
164  MoFEMFunctionReturn(0);
165  };
166 
167  CHKERR iterate_meshsets(
168 
169  m_field.getInterface<MeshsetsManager>()->getCubitMeshsetPtr(std::regex(
170 
171  (boost::format("%s(.*)") % block_name).str()
172 
173  ))
174 
175  );
176 
177  MoFEMFunctionReturn(0);
178  };
179 
180  CHKERR fix_disp();
181 
182  MoFEMFunctionReturn(0);
183 }
184 
185 MoFEMErrorCode BcManager::addBlockDOFsToMPCs(const std::string problem_name,
186  const std::string field_name,
187  bool get_low_dim_ents,
188  bool block_name_field_prefix,
189  bool is_distributed_mesh) {
190  Interface &m_field = cOre;
191  MoFEMFunctionBegin;
192 
193  if (block_name_field_prefix)
194  MOFEM_LOG("BcMngWorld", Sev::warning)
195  << "Argument block_name_field_prefix=true has no effect";
196  if (is_distributed_mesh)
197  MOFEM_LOG("BcMngWorld", Sev::warning)
198  << "Argument is_distributed_mesh=true has no effect";
199  if (get_low_dim_ents)
200  MOFEM_LOG("BcMngWorld", Sev::warning)
201  << "Argument get_low_dim_ents=true has no effect";
202 
203  auto get_dim = [&](const Range &ents) {
204  for (auto d : {3, 2, 1})
205  if (ents.num_of_dimension(d))
206  return d;
207  return 0;
208  };
209 
210  // auto mark_fix_dofs = [&](std::vector<unsigned char> &marked_field_dofs,
211  // const auto lo, const auto hi) {
212  // return prb_mng->modifyMarkDofs(problem_name, ROW, field_name, lo, hi,
213  // ProblemsManager::MarkOP::OR, 1,
214  // marked_field_dofs);
215  // };
216 
217  auto iterate_mpc_meshsets = [&]() {
218  MoFEMFunctionBegin;
219 
220  auto mpc_meshset_ptr =
221  m_field.getInterface<MeshsetsManager>()->getCubitMeshsetPtr(
222  std::regex((boost::format("%s(.*)") % "MPC_(.*)").str()));
223 
224  for (auto m : mpc_meshset_ptr) {
225 
226  if (std::regex_match(m->getName(),
227  std::regex("(.*)COUPLING_LINKS(.*)"))) {
228 
229  auto bc = boost::make_shared<BCs>();
230  bc->mpcPtr = boost::make_shared<MPCsType>();
231  bc->mpcPtr->mpcType = MPC::COUPLING;
232 
233  std::string const corresponding_master_ms =
234  std::regex_replace(m->getName(), std::regex("LINKS"), "MASTER");
235 
236  Range links_ents;
237  CHKERR m_field.get_moab().get_entities_by_handle(m->getMeshset(),
238  links_ents, true);
239 
240  Range master_nodes;
241  if (m_field.getInterface<MeshsetsManager>()->checkMeshset(
242  corresponding_master_ms)) {
243  const CubitMeshSets *l;
244  CHKERR m_field.getInterface<MeshsetsManager>()->getCubitMeshsetPtr(
245  corresponding_master_ms, &l);
246  bc->mpcPtr->isReprocitical = false;
247 
248  // std::cout << "master meshset is: " << master_nodes << "\n";
249  CHKERR m_field.get_moab().get_entities_by_handle(l->getMeshset(),
250  master_nodes, true);
251  // std::cout << "master meshset is: " << master_nodes << "\n";
252  // if (master_nodes.subset_by_dimension(0).size() <
253  // links_ents.subset_by_dimension(1))
254  {
255  auto low_dim_ents = BcManagerImplTools::get_adj_ents(
256  m_field.get_moab(), master_nodes);
257  // std::cout << "lower_dim meshset is: " << low_dim_ents << "\n";
258  low_dim_ents = low_dim_ents.subset_by_dimension(0);
259  master_nodes.swap(low_dim_ents);
260  // std::cout << "master meshset 2 is: " << master_nodes << "\n";
261  }
262 
263  MOFEM_LOG("BcMngWorld", Sev::verbose)
264  << "Found block MASTER LINKS block: " << l->getName()
265  << " Entities size: " << master_nodes.size();
266 
267  } else {
268  MOFEM_LOG("BcMngWorld", Sev::warning)
269  << "MASTER LINKS block not found: " << corresponding_master_ms
270  << " setting reprocitical constraint. ("
271  << bc->mpcPtr->isReprocitical << ").";
272  }
273 
274  // if (std::regex_match(bc_id, std::regex("(.*)TIE(.*)")))
275  // bc->mpcPtr->mpcType = MPC::TIE;
276  // if (std::regex_match(bc_id, std::regex("(.*)RIGID_BODY(.*)")))
277  // bc->mpcPtr->mpcType = MPC::RIGID_BODY;
278 
279  // std::cout << "links_ents range is: " << links_ents << "\n";
280  for (auto &link : links_ents.subset_by_dimension(1)) {
281  Range verts;
282  CHKERR m_field.get_moab().get_connectivity(&link, 1, verts, true);
283  // std::cout << "verts range is: " << verts << "\n";
284  if (bc->mpcPtr->isReprocitical) {
285  bc->mpcPtr->mpcMasterEnts.insert(verts[0]);
286  bc->mpcPtr->mpcSlaveEnts.insert(verts[1]);
287  } else {
288  for (auto &m_node : verts)
289  if (master_nodes.find(m_node) != master_nodes.end()) {
290  // std::cout << "found ent: " << m_node << "\n";
291  bc->mpcPtr->mpcMasterEnts.insert(m_node);
292  bc->mpcPtr->mpcSlaveEnts.merge(
293  subtract(verts, Range(m_node, m_node)));
294  break;
295  }
296  }
297  }
298 
299  MOFEM_LOG("BcMngWorld", Sev::verbose)
300  << "Found block MPC LINKS block: " << m->getName()
301  << " Entities size (edges): " << links_ents.size()
302  << " Entities size (nodes): "
303  << bc->mpcPtr->mpcMasterEnts.size() +
304  bc->mpcPtr->mpcSlaveEnts.size()
305  << " (" << bc->mpcPtr->mpcMasterEnts.size() << " "
306  << bc->mpcPtr->mpcSlaveEnts.size() << ")";
307 
308  MOFEM_LOG("BcMngSync", Sev::noisy) << *bc->mpcPtr;
309  MOFEM_LOG_SEVERITY_SYNC(m_field.get_comm(), Sev::noisy);
310 
311  // CHKERR m_field.getInterface<CommInterface>()->synchroniseEntities(
312  // bc->mpcPtr->mpcMasterEnts);
313  // CHKERR m_field.getInterface<CommInterface>()->synchroniseEntities(
314  // bc->mpcPtr->mpcSlaveEnts);
315  vector<double> mAttributes;
316  CHKERR m->getAttributes(mAttributes);
317 
318  auto setFlags = [&](const auto &flags) {
319  auto &d = bc->mpcPtr->data;
320  if (flags.empty()) {
321  d.flag1 = d.flag2 = d.flag3 = d.flag4 = d.flag5 = d.flag6 = true;
322  return;
323  }
324  for (size_t i = 0; i < std::min(flags.size(), size_t(6)); ++i)
325  (&d.flag1)[i] = flags[i] > 0.0;
326  };
327 
328  setFlags(mAttributes);
329 
330  MOFEM_LOG("BcMngWorld", Sev::verbose)
331  << "Found block " << m->getName() << " number of entities "
332  << bc->mpcPtr->mpcMasterEnts.size() << " number of attributes "
333  << mAttributes.size() << " highest dim of entities "
334  << get_dim(bc->mpcPtr->mpcMasterEnts);
335 
336  // NOTE: we are not using markers at the moment
337  // for (int i = 0; i != mAttributes.size(); ++i) {
338  // if (mAttributes[i] > 0.0)
339  // CHKERR mark_fix_dofs(bc->mpcPtr->bcMasterMarkers, i, i);
340  // }
341 
342  // if (mAttributes.empty())
343  // CHKERR mark_fix_dofs(bc->mpcPtr->bcMasterMarkers, 0,
344  // MAX_DOFS_ON_ENTITY);
345 
346  const std::string bc_id =
347  problem_name + "_" + field_name + "_" + m->getName();
348  bcMapByBlockName[bc_id] = bc;
349  }
350  }
351  MoFEMFunctionReturn(0);
352  };
353 
354  CHKERR iterate_mpc_meshsets();
355 
356  MoFEMFunctionReturn(0);
357 }
358 
359 boost::shared_ptr<BcManager::BCs>
360 BcManager::popMarkDOFsOnEntities(const std::string block_name) {
361  auto bc_it = bcMapByBlockName.find(block_name);
362  if (bc_it != bcMapByBlockName.end()) {
363  auto bc = bc_it->second;
364  bcMapByBlockName.erase(bc_it);
365  return bc;
366  }
367  return boost::shared_ptr<BCs>();
368 }
369 
370 Range BcManager::getMergedBlocksRange(std::vector<std::regex> bc_regex_vec) {
371  Range ents;
372  if (bcMapByBlockName.size()) {
373  for (auto b : bcMapByBlockName) {
374  for (auto &reg_name : bc_regex_vec) {
375  if (std::regex_match(b.first, reg_name)) {
376  ents.merge(b.second->bcEnts);
377  }
378  }
379  }
380  }
381  return ents;
382 }
383 
384 BcManager::BcMarkerPtr
385 BcManager::getMergedBlocksMarker(std::vector<std::regex> bc_regex_vec) {
386  BcManager::BcMarkerPtr boundary_marker_ptr;
387  if (bcMapByBlockName.size()) {
388  for (auto b : bcMapByBlockName) {
389  for (auto &reg_name : bc_regex_vec) {
390  if (std::regex_match(b.first, reg_name)) {
391  if (!boundary_marker_ptr)
392  boundary_marker_ptr =
393  boost::make_shared<std::vector<char unsigned>>();
394  boundary_marker_ptr->resize(b.second->bcMarkers.size(), 0);
395  for (int i = 0; i != b.second->bcMarkers.size(); ++i) {
396  (*boundary_marker_ptr)[i] |= b.second->bcMarkers[i];
397  }
398  }
399  }
400  }
401  }
402  return boundary_marker_ptr;
403 }
404 
405 BcManager::BcMarkerPtr BcManager::getMergedBlocksMarker(
406  const std::vector<BcManager::BcMarkerPtr> &boundary_markers_ptr_vec) {
407  auto boundary_marker_ptr = boost::make_shared<std::vector<char unsigned>>();
408  for (auto &bcm : boundary_markers_ptr_vec) {
409  boundary_marker_ptr->resize(bcm->size(), 0);
410  for (int i = 0; i != bcm->size(); ++i)
411  (*boundary_marker_ptr)[i] |= (*bcm)[i];
412  }
413  return boundary_marker_ptr;
414 }
415 
416 SmartPetscObj<IS> BcManager::getBlockIS(const std::string block_prefix,
417  const std::string block_name,
418  const std::string field_name,
419  const std::string problem_name, int lo,
420  int hi, SmartPetscObj<IS> is_expand) {
421  Interface &m_field = cOre;
422 
423  const std::string bc_id =
424  block_prefix + "_" + field_name + "_" + block_name + "(.*)";
425 
426  Range bc_ents;
427  for (auto bc : getBcMapByBlockName()) {
428  if (std::regex_match(bc.first, std::regex(bc_id))) {
429  bc_ents.merge(*(bc.second->getBcEntsPtr()));
430  MOFEM_LOG("BcMngWorld", Sev::verbose)
431  << "Get entities from block and add to IS. Block name " << bc.first;
432  }
433  }
434 
435  SmartPetscObj<IS> is_bc;
436  auto get_is = [&]() {
437  MoFEMFunctionBegin;
438  CHKERR m_field.getInterface<CommInterface>()->synchroniseEntities(bc_ents);
439  CHKERR m_field.getInterface<ISManager>()->isCreateProblemFieldAndRank(
440  problem_name, ROW, field_name, lo, hi, is_bc, &bc_ents);
441  if (is_expand) {
442  IS is_tmp;
443  CHKERR ISExpand(is_bc, is_expand, &is_tmp);
444  is_bc = SmartPetscObj<IS>(is_tmp);
445  }
446  CHKERR ISSort(is_bc);
447  MoFEMFunctionReturn(0);
448  };
449 
450  if (get_is())
451  CHK_THROW_MESSAGE(MOFEM_DATA_INCONSISTENCY, "IS is not created");
452 
453  return is_bc;
454 }
455 
456 SmartPetscObj<IS> BcManager::getBlockIS(const std::string problem_name,
457  const std::string block_name,
458  const std::string field_name, int lo,
459  int hi, SmartPetscObj<IS> is_expand) {
460  return getBlockIS(problem_name, block_name, field_name, problem_name, lo, hi,
461  is_expand);
462 }
463 
464 template <>
465 MoFEMErrorCode
466 BcManager::removeBlockDOFsOnEntities<BcMeshsetType<DISPLACEMENTSET>>(
467  const std::string problem_name, const std::string field_name,
468  bool get_low_dim_ents, bool block_name_field_prefix,
469  bool is_distributed_mesh) {
470  Interface &m_field = cOre;
471  auto prb_mng = m_field.getInterface<ProblemsManager>();
472  MoFEMFunctionBegin;
473 
474  CHKERR pushMarkDOFsOnEntities<BcMeshsetType<DISPLACEMENTSET>>(
475  problem_name, field_name, get_low_dim_ents, block_name_field_prefix);
476 
477  std::array<Range, 3> ents_to_remove;
478 
479  for (auto m :
480 
481  m_field.getInterface<MeshsetsManager>()->getCubitMeshsetPtr(
482  NODESET | DISPLACEMENTSET)) {
483 
484  const std::string bc_id =
485  problem_name + "_" + field_name + "_DISPLACEMENTSET" +
486  boost::lexical_cast<std::string>(m->getMeshsetId());
487 
488  auto bc = bcMapByBlockName.at(bc_id);
489 
490  if (bc->dispBcPtr) {
491  if (bc->dispBcPtr->data.flag1) {
492  ents_to_remove[0].merge(bc->bcEnts);
493  }
494  if (bc->dispBcPtr->data.flag2) {
495  ents_to_remove[1].merge(bc->bcEnts);
496  }
497  if (bc->dispBcPtr->data.flag3) {
498  ents_to_remove[2].merge(bc->bcEnts);
499  }
500  if (bc->dispBcPtr->data.flag4) {
501  ents_to_remove[1].merge(bc->bcEnts);
502  ents_to_remove[2].merge(bc->bcEnts);
503  }
504  if (bc->dispBcPtr->data.flag5) {
505  ents_to_remove[0].merge(bc->bcEnts);
506  ents_to_remove[2].merge(bc->bcEnts);
507  }
508  if (bc->dispBcPtr->data.flag6) {
509  ents_to_remove[0].merge(bc->bcEnts);
510  ents_to_remove[1].merge(bc->bcEnts);
511  }
512  }
513  bc->bcMarkers = std::vector<unsigned char>();
514  }
515 
516  auto remove_dofs_on_ents = [&](const Range &ents, const int lo,
517  const int hi) {
518  if (is_distributed_mesh)
519  return prb_mng->removeDofsOnEntities(problem_name, field_name, ents, lo,
520  hi);
521  else
522  return prb_mng->removeDofsOnEntitiesNotDistributed(
523  problem_name, field_name, ents, lo, hi);
524  };
525 
526  CHKERR remove_dofs_on_ents(ents_to_remove[0], 0, 0);
527  CHKERR remove_dofs_on_ents(ents_to_remove[1], 1, 1);
528  CHKERR remove_dofs_on_ents(ents_to_remove[2], 2, 2);
529 
530  MoFEMFunctionReturn(0);
531 }
532 
533 template <>
534 MoFEMErrorCode
535 BcManager::removeBlockDOFsOnEntities<BcMeshsetType<TEMPERATURESET>>(
536  const std::string problem_name, const std::string field_name,
537  bool get_low_dim_ents, bool block_name_field_prefix,
538  bool is_distributed_mesh) {
539  Interface &m_field = cOre;
540  auto prb_mng = m_field.getInterface<ProblemsManager>();
541  MoFEMFunctionBegin;
542 
543  CHKERR pushMarkDOFsOnEntities<BcMeshsetType<TEMPERATURESET>>(
544  problem_name, field_name, get_low_dim_ents, block_name_field_prefix);
545 
546  Range ents_to_remove;
547 
548  for (auto m :
549 
550  m_field.getInterface<MeshsetsManager>()->getCubitMeshsetPtr(
551  NODESET | TEMPERATURESET)
552 
553  ) {
554  const std::string bc_id =
555  problem_name + "_" + field_name + "_TEMPERATURESET" +
556  boost::lexical_cast<std::string>(m->getMeshsetId());
557  auto bc = bcMapByBlockName.at(bc_id);
558  ents_to_remove.merge(bc->bcEnts);
559  bc->bcMarkers = std::vector<unsigned char>();
560  }
561 
562  auto remove_dofs_on_ents = [&](const Range &ents, const int lo,
563  const int hi) {
564  if (is_distributed_mesh)
565  return prb_mng->removeDofsOnEntities(problem_name, field_name, ents, lo,
566  hi);
567  else
568  return prb_mng->removeDofsOnEntitiesNotDistributed(
569  problem_name, field_name, ents, lo, hi);
570  };
571 
572  CHKERR remove_dofs_on_ents(ents_to_remove, 0, MAX_DOFS_ON_ENTITY);
573 
574  MoFEMFunctionReturn(0);
575 }
576 
577 template <>
578 MoFEMErrorCode BcManager::removeBlockDOFsOnEntities<BcMeshsetType<HEATFLUXSET>>(
579  const std::string problem_name, const std::string field_name,
580  bool get_low_dim_ents, bool block_name_field_prefix,
581  bool is_distributed_mesh) {
582  Interface &m_field = cOre;
583  auto prb_mng = m_field.getInterface<ProblemsManager>();
584  MoFEMFunctionBegin;
585 
586  CHKERR pushMarkDOFsOnEntities<BcMeshsetType<HEATFLUXSET>>(
587  problem_name, field_name, get_low_dim_ents, block_name_field_prefix);
588 
589  Range ents_to_remove;
590 
591  for (auto m :
592 
593  m_field.getInterface<MeshsetsManager>()->getCubitMeshsetPtr(NODESET |
594  HEATFLUXSET)
595 
596  ) {
597  const std::string bc_id =
598  problem_name + "_" + field_name + "_HEATFLUXSET" +
599  boost::lexical_cast<std::string>(m->getMeshsetId());
600  auto bc = bcMapByBlockName.at(bc_id);
601  ents_to_remove.merge(bc->bcEnts);
602  bc->bcMarkers = std::vector<unsigned char>();
603  }
604 
605  auto remove_dofs_on_ents = [&](const Range &ents, const int lo,
606  const int hi) {
607  if (is_distributed_mesh)
608  return prb_mng->removeDofsOnEntities(problem_name, field_name, ents, lo,
609  hi);
610  else
611  return prb_mng->removeDofsOnEntitiesNotDistributed(
612  problem_name, field_name, ents, lo, hi);
613  };
614 
615  CHKERR remove_dofs_on_ents(ents_to_remove, 0, MAX_DOFS_ON_ENTITY);
616 
617  MoFEMFunctionReturn(0);
618 }
619 
620 template <>
621 MoFEMErrorCode
622 BcManager::removeBlockDOFsOnEntities<BcDisplacementMeshsetType<BLOCKSET>>(
623  const std::string problem_name, const std::string field_name,
624  bool get_low_dim_ents, bool block_name_field_prefix,
625  bool is_distributed_mesh) {
626  Interface &m_field = cOre;
627  auto prb_mng = m_field.getInterface<ProblemsManager>();
628  MoFEMFunctionBegin;
629 
630  CHKERR pushMarkDOFsOnEntities<BcDisplacementMeshsetType<BLOCKSET>>(
631  problem_name, field_name, get_low_dim_ents, block_name_field_prefix);
632 
633  std::array<Range, 3> ents_to_remove;
634 
635  for (auto m :
636 
637  m_field.getInterface<MeshsetsManager>()->getCubitMeshsetPtr(
638  BLOCKSET | UNKNOWNNAME)) {
639 
640  const auto block_name = m->getName();
641 
642  std::string bc_id = problem_name + "_" + field_name + "_" + block_name;
643  std::string regex_str;
644  if (block_name_field_prefix) {
645  regex_str = (boost::format("%s_%s_%s_((FIX_(ALL|X|Y|Z))|("
646  "DISPLACEMENT|ROTATE))(.*)") %
647  problem_name % field_name % field_name)
648  .str();
649  } else {
650  regex_str = (boost::format("%s_%s_((FIX_(ALL|X|Y|Z))|("
651  "DISPLACEMENT|ROTATE))(.*)") %
652  problem_name % field_name)
653  .str();
654  }
655 
656  if (std::regex_match(bc_id, std::regex(regex_str))) {
657 
658  auto bc = bcMapByBlockName.at(bc_id);
659 
660  if (auto disp_bc = bc->dispBcPtr) {
661  if (disp_bc->data.flag1) {
662  ents_to_remove[0].merge(bc->bcEnts);
663  }
664  if (disp_bc->data.flag2) {
665  ents_to_remove[1].merge(bc->bcEnts);
666  }
667  if (disp_bc->data.flag3) {
668  ents_to_remove[2].merge(bc->bcEnts);
669  }
670  if (disp_bc->data.flag4) {
671  ents_to_remove[1].merge(bc->bcEnts);
672  ents_to_remove[2].merge(bc->bcEnts);
673  }
674  if (disp_bc->data.flag5) {
675  ents_to_remove[0].merge(bc->bcEnts);
676  ents_to_remove[2].merge(bc->bcEnts);
677  }
678  if (disp_bc->data.flag6) {
679  ents_to_remove[0].merge(bc->bcEnts);
680  ents_to_remove[1].merge(bc->bcEnts);
681  }
682  } else {
683  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
684  "BC type not implemented");
685  }
686  }
687  }
688 
689  auto remove_dofs_on_ents = [&](const Range &ents, const int lo,
690  const int hi) {
691  if (is_distributed_mesh)
692  return prb_mng->removeDofsOnEntities(problem_name, field_name, ents, lo,
693  hi);
694  else
695  return prb_mng->removeDofsOnEntitiesNotDistributed(
696  problem_name, field_name, ents, lo, hi);
697  };
698 
699  CHKERR remove_dofs_on_ents(ents_to_remove[0], 0, 0);
700  CHKERR remove_dofs_on_ents(ents_to_remove[1], 1, 1);
701  CHKERR remove_dofs_on_ents(ents_to_remove[2], 2, 2);
702 
703  MoFEMFunctionReturn(0);
704 }
705 
706 template <>
707 MoFEMErrorCode
708 BcManager::removeBlockDOFsOnEntities<BcScalarMeshsetType<BLOCKSET>>(
709  const std::string problem_name, const std::string block_name,
710  const std::string field_name, bool get_low_dim_ents,
711  bool is_distributed_mesh) {
712  Interface &m_field = cOre;
713  auto prb_mng = m_field.getInterface<ProblemsManager>();
714  MoFEMFunctionBegin;
715 
716  CHKERR pushMarkDOFsOnEntities<BcScalarMeshsetType<BLOCKSET>>(
717  problem_name, block_name, field_name, get_low_dim_ents);
718 
719  Range ents_to_remove;
720 
721  for (auto m :
722 
723  m_field.getInterface<MeshsetsManager>()->getCubitMeshsetPtr(
724  BLOCKSET | UNKNOWNNAME)) {
725 
726  std::string bc_id = problem_name + "_" + field_name + "_" + m->getName();
727 
728  auto str = boost::format("%s_%s_%s(.*)")
729 
730  % problem_name % field_name % block_name;
731 
732  if (std::regex_match(bc_id, std::regex(str.str()))) {
733 
734  auto bc = bcMapByBlockName.at(bc_id);
735 
736  if (auto temp_bc = bc->tempBcPtr) {
737  if (temp_bc->data.flag1) {
738  ents_to_remove.merge(bc->bcEnts);
739  }
740  } else {
741  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
742  "BC type not implemented");
743  }
744  }
745  }
746 
747  auto remove_dofs_on_ents = [&](const Range &ents, const int lo,
748  const int hi) {
749  if (is_distributed_mesh)
750  return prb_mng->removeDofsOnEntities(problem_name, field_name, ents, lo,
751  hi);
752  else
753  return prb_mng->removeDofsOnEntitiesNotDistributed(
754  problem_name, field_name, ents, lo, hi);
755  };
756 
757  CHKERR remove_dofs_on_ents(ents_to_remove, 0, MAX_DOFS_ON_ENTITY);
758 
759  MoFEMFunctionReturn(0);
760 }
761 
762 template <>
763 MoFEMErrorCode
764 BcManager::pushMarkDOFsOnEntities<BcMeshsetType<DISPLACEMENTSET>>(
765  const std::string problem_name, const std::string field_name,
766  bool get_low_dim_ents, bool block_name_field_prefix) {
767  Interface &m_field = cOre;
768  auto prb_mng = m_field.getInterface<ProblemsManager>();
769  MoFEMFunctionBegin;
770 
771  // if(problem_name.size()==0)
772  // MOFEM_LOG("BcMngWorld", Sev::warning)
773  // << "Argument problem_name has no effect";
774 
775  if (block_name_field_prefix)
776  MOFEM_LOG("BcMngWorld", Sev::warning)
777  << "Argument block_name_field_prefix=true has no effect";
778 
779  auto fix_disp = [&]() {
780  MoFEMFunctionBegin;
781 
782  auto iterate_meshsets = [&](auto &&meshset_vec_ptr) {
783  MoFEMFunctionBegin;
784  for (auto m : meshset_vec_ptr) {
785  auto bc = boost::make_shared<BCs>();
786  CHKERR m_field.get_moab().get_entities_by_handle(m->getMeshset(),
787  bc->bcEnts, true);
788  bc->dispBcPtr = boost::make_shared<DisplacementCubitBcData>();
789  CHKERR m->getBcDataStructure(*(bc->dispBcPtr));
790 
791  MOFEM_LOG("BcMngWorld", Sev::verbose)
792  << "Found block DISPLACEMENTSET id = " << m->getMeshsetId();
793  MOFEM_LOG("BcMngWorld", Sev::verbose) << *bc->dispBcPtr;
794 
795  MOFEM_LOG("BcMngSync", Sev::noisy)
796  << "Found block DISPLACEMENTSET id = " << m->getMeshsetId()
797  << " nb. of entities " << bc->bcEnts.size()
798  << " highest dim of entities "
799  << BcManagerImplTools::get_dim(bc->bcEnts);
800  MOFEM_LOG("BcMngSync", Sev::noisy) << *bc->dispBcPtr;
801  MOFEM_LOG_SEVERITY_SYNC(m_field.get_comm(), Sev::noisy);
802 
803  if (problem_name.size()) {
804 
805  if (bc->dispBcPtr->data.flag1)
806  CHKERR prb_mng->modifyMarkDofs(problem_name, ROW, field_name, 0, 0,
807  ProblemsManager::MarkOP::OR, 1,
808  bc->bcMarkers);
809  if (bc->dispBcPtr->data.flag2)
810  CHKERR prb_mng->modifyMarkDofs(problem_name, ROW, field_name, 1, 1,
811  ProblemsManager::MarkOP::OR, 1,
812  bc->bcMarkers);
813  if (bc->dispBcPtr->data.flag3)
814  CHKERR prb_mng->modifyMarkDofs(problem_name, ROW, field_name, 2, 2,
815  ProblemsManager::MarkOP::OR, 1,
816  bc->bcMarkers);
817  if (bc->dispBcPtr->data.flag4) {
818 
819  CHKERR prb_mng->modifyMarkDofs(problem_name, ROW, field_name, 1, 1,
820  ProblemsManager::MarkOP::OR, 1,
821  bc->bcMarkers);
822  CHKERR prb_mng->modifyMarkDofs(problem_name, ROW, field_name, 2, 2,
823  ProblemsManager::MarkOP::OR, 1,
824  bc->bcMarkers);
825  }
826  if (bc->dispBcPtr->data.flag5) {
827 
828  CHKERR prb_mng->modifyMarkDofs(problem_name, ROW, field_name, 0, 0,
829  ProblemsManager::MarkOP::OR, 1,
830  bc->bcMarkers);
831  CHKERR prb_mng->modifyMarkDofs(problem_name, ROW, field_name, 2, 2,
832  ProblemsManager::MarkOP::OR, 1,
833  bc->bcMarkers);
834  }
835  if (bc->dispBcPtr->data.flag6) {
836 
837  CHKERR prb_mng->modifyMarkDofs(problem_name, ROW, field_name, 0, 0,
838  ProblemsManager::MarkOP::OR, 1,
839  bc->bcMarkers);
840  CHKERR prb_mng->modifyMarkDofs(problem_name, ROW, field_name, 1, 1,
841  ProblemsManager::MarkOP::OR, 1,
842  bc->bcMarkers);
843  }
844  }
845 
846  if (get_low_dim_ents) {
847  auto low_dim_ents =
848  BcManagerImplTools::get_adj_ents(m_field.get_moab(), bc->bcEnts);
849  bc->bcEnts.swap(low_dim_ents);
850  }
851 
852  CHKERR m_field.getInterface<CommInterface>()->synchroniseEntities(
853  bc->bcEnts);
854  if (problem_name.size())
855  CHKERR prb_mng->markDofs(problem_name, ROW, ProblemsManager::AND,
856  bc->bcEnts, bc->bcMarkers);
857 
858  const std::string bc_id =
859  problem_name + "_" + field_name + "_DISPLACEMENTSET" +
860  boost::lexical_cast<std::string>(m->getMeshsetId());
861  bcMapByBlockName[bc_id] = bc;
862  }
863  MoFEMFunctionReturn(0);
864  };
865 
866  CHKERR iterate_meshsets(
867 
868  m_field.getInterface<MeshsetsManager>()->getCubitMeshsetPtr(
869  NODESET | DISPLACEMENTSET)
870 
871  );
872 
873  MoFEMFunctionReturn(0);
874  };
875 
876  CHKERR fix_disp();
877 
878  MoFEMFunctionReturn(0);
879 }
880 
881 template <>
882 MoFEMErrorCode BcManager::pushMarkDOFsOnEntities<BcMeshsetType<TEMPERATURESET>>(
883  const std::string problem_name, const std::string field_name,
884  bool get_low_dim_ents, bool block_name_field_prefix) {
885  Interface &m_field = cOre;
886  auto prb_mng = m_field.getInterface<ProblemsManager>();
887  MoFEMFunctionBegin;
888 
889  if (block_name_field_prefix)
890  MOFEM_LOG("BcMngWorld", Sev::warning)
891  << "Argument block_name_field_prefix=true has no effect";
892 
893  auto fix_temp = [&]() {
894  MoFEMFunctionBegin;
895 
896  auto iterate_meshsets = [&](auto &&meshset_vec_ptr) {
897  MoFEMFunctionBegin;
898  for (auto m : meshset_vec_ptr) {
899  auto bc = boost::make_shared<BCs>();
900  CHKERR m_field.get_moab().get_entities_by_handle(m->getMeshset(),
901  bc->bcEnts, true);
902  bc->tempBcPtr = boost::make_shared<TemperatureCubitBcData>();
903  CHKERR m->getBcDataStructure(*(bc->tempBcPtr));
904 
905  MOFEM_LOG("BcMngWorld", Sev::verbose)
906  << "Found block TEMPERATURESET id = " << m->getMeshsetId();
907  MOFEM_LOG("BcMngWorld", Sev::verbose) << *bc->tempBcPtr;
908 
909  CHKERR prb_mng->modifyMarkDofs(
910  problem_name, ROW, field_name, 0, MAX_DOFS_ON_ENTITY,
911  ProblemsManager::MarkOP::OR, 1, bc->bcMarkers);
912 
913  if (get_low_dim_ents) {
914  auto low_dim_ents =
915  BcManagerImplTools::get_adj_ents(m_field.get_moab(), bc->bcEnts);
916  bc->bcEnts.swap(low_dim_ents);
917  }
918 
919  CHKERR m_field.getInterface<CommInterface>()->synchroniseEntities(
920  bc->bcEnts);
921  CHKERR prb_mng->markDofs(problem_name, ROW, ProblemsManager::AND,
922  bc->bcEnts, bc->bcMarkers);
923 
924  const std::string bc_id =
925  problem_name + "_" + field_name + "_TEMPERATURESET" +
926  boost::lexical_cast<std::string>(m->getMeshsetId());
927  bcMapByBlockName[bc_id] = bc;
928  }
929  MoFEMFunctionReturn(0);
930  };
931 
932  CHKERR iterate_meshsets(
933 
934  m_field.getInterface<MeshsetsManager>()->getCubitMeshsetPtr(
935  NODESET | TEMPERATURESET)
936 
937  );
938 
939  MoFEMFunctionReturn(0);
940  };
941 
942  CHKERR fix_temp();
943 
944  MoFEMFunctionReturn(0);
945 }
946 
947 template <>
948 MoFEMErrorCode BcManager::pushMarkDOFsOnEntities<BcMeshsetType<HEATFLUXSET>>(
949  const std::string problem_name, const std::string field_name,
950  bool get_low_dim_ents, bool block_name_field_prefix) {
951  Interface &m_field = cOre;
952  auto prb_mng = m_field.getInterface<ProblemsManager>();
953  MoFEMFunctionBegin;
954 
955  if (block_name_field_prefix)
956  MOFEM_LOG("BcMngWorld", Sev::warning)
957  << "Argument block_name_field_prefix=true has no effect";
958 
959  auto fix_disp = [&]() {
960  MoFEMFunctionBegin;
961 
962  auto iterate_meshsets = [&](auto &&meshset_vec_ptr) {
963  MoFEMFunctionBegin;
964  for (auto m : meshset_vec_ptr) {
965  auto bc = boost::make_shared<BCs>();
966  CHKERR m_field.get_moab().get_entities_by_handle(m->getMeshset(),
967  bc->bcEnts, true);
968  bc->heatFluxBcPtr = boost::make_shared<HeatFluxCubitBcData>();
969  CHKERR m->getBcDataStructure(*(bc->heatFluxBcPtr));
970 
971  CHKERR prb_mng->modifyMarkDofs(
972  problem_name, ROW, field_name, 0, MAX_DOFS_ON_ENTITY,
973  ProblemsManager::MarkOP::OR, 1, bc->bcMarkers);
974 
975  MOFEM_LOG("BcMngWorld", Sev::verbose)
976  << "Found block HEATFLUX id = " << m->getMeshsetId();
977  MOFEM_LOG("BcMngWorld", Sev::verbose) << *bc->heatFluxBcPtr;
978 
979  if (get_low_dim_ents) {
980  auto low_dim_ents =
981  BcManagerImplTools::get_adj_ents(m_field.get_moab(), bc->bcEnts);
982  bc->bcEnts.swap(low_dim_ents);
983  }
984 
985  CHKERR m_field.getInterface<CommInterface>()->synchroniseEntities(
986  bc->bcEnts);
987  CHKERR prb_mng->markDofs(problem_name, ROW, ProblemsManager::AND,
988  bc->bcEnts, bc->bcMarkers);
989 
990  const std::string bc_id =
991  problem_name + "_" + field_name + "_HEATFLUXSET" +
992  boost::lexical_cast<std::string>(m->getMeshsetId());
993  bcMapByBlockName[bc_id] = bc;
994  }
995  MoFEMFunctionReturn(0);
996  };
997 
998  CHKERR iterate_meshsets(
999 
1000  m_field.getInterface<MeshsetsManager>()->getCubitMeshsetPtr(SIDESET |
1001  HEATFLUXSET)
1002 
1003  );
1004 
1005  MoFEMFunctionReturn(0);
1006  };
1007 
1008  CHKERR fix_disp();
1009 
1010  MoFEMFunctionReturn(0);
1011 }
1012 
1013 template <>
1014 MoFEMErrorCode
1015 BcManager::pushMarkDOFsOnEntities<BcDisplacementMeshsetType<BLOCKSET>>(
1016  const std::string problem_name, const std::string field_name,
1017  bool get_low_dim_ents, bool block_name_field_prefix) {
1018  MoFEMFunctionBegin;
1019 
1020  auto mark_dofs = [&](const string block_name, const int &idx_0,
1021  const int &idx_1) {
1022  MoFEMFunctionBeginHot;
1023  if (block_name_field_prefix) {
1024  const string field_block = field_name + "_" + block_name;
1025  CHKERR pushMarkDOFsOnEntities(problem_name, field_block, field_name,
1026  idx_0, idx_1, get_low_dim_ents);
1027  } else {
1028 
1029  CHKERR pushMarkDOFsOnEntities(problem_name, block_name, field_name, idx_0,
1030  idx_1, get_low_dim_ents);
1031  }
1032  MoFEMFunctionReturnHot(0);
1033  };
1034 
1035  // displacement
1036  CHKERR mark_dofs("FIX_X", 0, 0);
1037  CHKERR mark_dofs("FIX_Y", 1, 1);
1038  CHKERR mark_dofs("FIX_Z", 2, 2);
1039  CHKERR mark_dofs("FIX_ALL", 0, MAX_DOFS_ON_ENTITY);
1040 
1041  // rotation
1042  CHKERR mark_dofs("ROTATE_X", 1, 1);
1043  CHKERR mark_dofs("ROTATE_X", 2, 2);
1044  CHKERR mark_dofs("ROTATE_Y", 0, 0);
1045  CHKERR mark_dofs("ROTATE_Y", 2, 2);
1046  CHKERR mark_dofs("ROTATE_Z", 0, 0);
1047  CHKERR mark_dofs("ROTATE_Z", 1, 1);
1048  CHKERR mark_dofs("ROTATE_ALL", 0, MAX_DOFS_ON_ENTITY);
1049 
1050  std::string regex_str;
1051  if (block_name_field_prefix) {
1052  regex_str = (boost::format("%s_%s_%s_(.*)") % problem_name % field_name %
1053  field_name)
1054  .str();
1055  } else {
1056  regex_str = (boost::format("%s_%s_(.*)") % problem_name % field_name).str();
1057  }
1058 
1059  for (auto &m : bcMapByBlockName) {
1060  auto &bc_id = m.first;
1061  if (std::regex_match(bc_id, std::regex(regex_str))) {
1062  auto &bc = m.second;
1063  if (std::regex_match(bc_id, std::regex("(.*)_FIX_X(.*)"))) {
1064  bc->dispBcPtr = boost::make_shared<DisplacementCubitBcData>();
1065  bc->dispBcPtr->data.flag1 = 1;
1066  if (bc->bcAttributes.empty()) {
1067  bc->dispBcPtr->data.value1 = 0;
1068  MOFEM_LOG("BcMngWorld", Sev::warning)
1069  << "Expected one attribute on block but have "
1070  << bc->bcAttributes.size();
1071  } else if (bc->bcAttributes.size() >= 1) {
1072  bc->dispBcPtr->data.value1 = bc->bcAttributes[0];
1073  }
1074  MOFEM_LOG("BcMngWorld", Sev::inform) << "Add X " << bc_id;
1075  MOFEM_LOG("BcMngWorld", Sev::inform) << *bc->dispBcPtr;
1076  } else if (std::regex_match(bc_id, std::regex("(.*)_FIX_Y(.*)"))) {
1077  bc->dispBcPtr = boost::make_shared<DisplacementCubitBcData>();
1078  bc->dispBcPtr->data.flag2 = 1;
1079  if (bc->bcAttributes.empty()) {
1080  bc->dispBcPtr->data.value2 = 0;
1081  MOFEM_LOG("BcMngWorld", Sev::warning)
1082  << "Expected one attribute on block but have "
1083  << bc->bcAttributes.size();
1084  } else if (bc->bcAttributes.size() == 1) {
1085  bc->dispBcPtr->data.value2 = bc->bcAttributes[0];
1086  } else if (bc->bcAttributes.size() >= 2) {
1087  bc->dispBcPtr->data.value2 = bc->bcAttributes[1];
1088  }
1089  MOFEM_LOG("BcMngWorld", Sev::inform) << "Add Y " << bc_id;
1090  MOFEM_LOG("BcMngWorld", Sev::inform) << *(bc->dispBcPtr);
1091  } else if (std::regex_match(bc_id, std::regex("(.*)_FIX_Z(.*)"))) {
1092  bc->dispBcPtr = boost::make_shared<DisplacementCubitBcData>();
1093  bc->dispBcPtr->data.flag3 = 1;
1094  if (bc->bcAttributes.empty()) {
1095  bc->dispBcPtr->data.value3 = 0;
1096  MOFEM_LOG("BcMngWorld", Sev::warning)
1097  << "Expected one attribute on block but have "
1098  << bc->bcAttributes.size();
1099  } else if (bc->bcAttributes.size() == 1) {
1100  bc->dispBcPtr->data.value3 = bc->bcAttributes[0];
1101  } else if (bc->bcAttributes.size() == 3) {
1102  bc->dispBcPtr->data.value3 = bc->bcAttributes[2];
1103  }
1104  MOFEM_LOG("BcMngWorld", Sev::inform) << "Add Z " << bc_id;
1105  MOFEM_LOG("BcMngWorld", Sev::inform) << *(bc->dispBcPtr);
1106  } else if (std::regex_match(bc_id, std::regex("(.*)_FIX_ALL(.*)"))) {
1107  bc->dispBcPtr = boost::make_shared<DisplacementCubitBcData>();
1108  bc->dispBcPtr->data.flag1 = 1;
1109  bc->dispBcPtr->data.flag2 = 1;
1110  bc->dispBcPtr->data.flag3 = 1;
1111  if (bc->bcAttributes.size() >= 1) {
1112  bc->dispBcPtr->data.value1 = bc->bcAttributes[0];
1113  }
1114  if (bc->bcAttributes.size() >= 2) {
1115  bc->dispBcPtr->data.value2 = bc->bcAttributes[1];
1116  }
1117  if (bc->bcAttributes.size() >= 3) {
1118  bc->dispBcPtr->data.value3 = bc->bcAttributes[2];
1119  }
1120  MOFEM_LOG("BcMngWorld", Sev::inform) << "Add ALL " << bc_id;
1121  MOFEM_LOG("BcMngWorld", Sev::inform) << *(bc->dispBcPtr);
1122  } else if (std::regex_match(bc_id, std::regex("(.*)_ROTATE_X(.*)"))) {
1123  bc->dispBcPtr =
1124  boost::make_shared<DisplacementCubitBcDataWithRotation>();
1125  bc->dispBcPtr->data.flag4 = 1;
1126  bc->dispBcPtr->data.flag5 = 0;
1127  bc->dispBcPtr->data.flag6 = 0;
1128  // for the ROTATE_X block the angles can be specified with either one or
1129  // three attributes, e.g. 1, coords or 1,0,0,coords
1130  if (bc->bcAttributes.empty()) {
1131  bc->dispBcPtr->data.value4 = 0;
1132  MOFEM_LOG("BcMngWorld", Sev::warning)
1133  << "Expected one attribute on block on block (angle (1 or 3), "
1134  "center coords(3) but have "
1135  << bc->bcAttributes.size();
1136  } else if (bc->bcAttributes.size() >= 1) {
1137  bc->dispBcPtr->data.value4 = bc->bcAttributes[0];
1138  }
1139  MOFEM_LOG("BcMngWorld", Sev::inform) << "Add X " << bc_id;
1140  MOFEM_LOG("BcMngWorld", Sev::inform) << *bc->dispBcPtr;
1141  if (bc->bcAttributes.size() == 4 || bc->bcAttributes.size() == 6) {
1142  if (auto ext_disp_bc =
1143  dynamic_cast<DisplacementCubitBcDataWithRotation *>(
1144  bc->dispBcPtr.get())) {
1145  auto &o = ext_disp_bc->rotOffset;
1146  for (int a = 0; a != 3; ++a)
1147  o[a] = bc->bcAttributes[bc->bcAttributes.size() - 3 + a];
1148  MOFEM_LOG("BcMngWorld", Sev::inform)
1149  << "Add Rotate X Center: " << o[0] << " " << o[1] << " "
1150  << o[2];
1151  }
1152  }
1153  } else if (std::regex_match(bc_id, std::regex("(.*)_ROTATE_Y(.*)"))) {
1154  bc->dispBcPtr =
1155  boost::make_shared<DisplacementCubitBcDataWithRotation>();
1156  bc->dispBcPtr->data.flag4 = 0;
1157  bc->dispBcPtr->data.flag5 = 1;
1158  bc->dispBcPtr->data.flag6 = 0;
1159  // for the ROTATE_Y block the angles can be specified with either one or
1160  // three attributes, e.g. 1, coords or 0,1,0,coords
1161  if (bc->bcAttributes.empty()) {
1162  bc->dispBcPtr->data.value5 = 0;
1163  MOFEM_LOG("BcMngWorld", Sev::warning)
1164  << "Expected one attribute on block on block (angle (1 or 3), "
1165  "center coords(3) but have "
1166  << bc->bcAttributes.size();
1167  } else if (bc->bcAttributes.size() == 1 ||
1168  bc->bcAttributes.size() == 4) {
1169  bc->dispBcPtr->data.value5 = bc->bcAttributes[0];
1170  } else if (bc->bcAttributes.size() == 6) {
1171  bc->dispBcPtr->data.value5 = bc->bcAttributes[1];
1172  }
1173  MOFEM_LOG("BcMngWorld", Sev::inform) << "Add Y " << bc_id;
1174  MOFEM_LOG("BcMngWorld", Sev::inform) << *(bc->dispBcPtr);
1175  if (bc->bcAttributes.size() == 4 || bc->bcAttributes.size() == 6) {
1176  if (auto ext_disp_bc =
1177  dynamic_cast<DisplacementCubitBcDataWithRotation *>(
1178  bc->dispBcPtr.get())) {
1179  auto &o = ext_disp_bc->rotOffset;
1180  for (int a = 0; a != 3; ++a)
1181  o[a] = bc->bcAttributes[bc->bcAttributes.size() - 3 + a];
1182  MOFEM_LOG("BcMngWorld", Sev::inform)
1183  << "Add Rotate Y Center: " << o[0] << " " << o[1] << " "
1184  << o[2];
1185  }
1186  }
1187  } else if (std::regex_match(bc_id, std::regex("(.*)_ROTATE_Z(.*)"))) {
1188  bc->dispBcPtr =
1189  boost::make_shared<DisplacementCubitBcDataWithRotation>();
1190  bc->dispBcPtr->data.flag4 = 0;
1191  bc->dispBcPtr->data.flag5 = 0;
1192  bc->dispBcPtr->data.flag6 = 1;
1193  // for the ROTATE_Z block the angles can be specified with either one or
1194  // three attributes, e.g. 1, coords or 0,0,1,coords
1195  if (bc->bcAttributes.empty()) {
1196  bc->dispBcPtr->data.value6 = 0;
1197  MOFEM_LOG("BcMngWorld", Sev::warning)
1198  << "Expected one attribute on block (angle (1 or 3), center "
1199  "coords(3) but have "
1200  << bc->bcAttributes.size();
1201  } else if (bc->bcAttributes.size() == 1 ||
1202  bc->bcAttributes.size() == 4) {
1203  bc->dispBcPtr->data.value6 = bc->bcAttributes[0];
1204  } else if (bc->bcAttributes.size() == 3 ||
1205  bc->bcAttributes.size() == 6) {
1206  bc->dispBcPtr->data.value6 = bc->bcAttributes[2];
1207  }
1208  MOFEM_LOG("BcMngWorld", Sev::inform) << "Add Z " << bc_id;
1209  MOFEM_LOG("BcMngWorld", Sev::inform) << *(bc->dispBcPtr);
1210  if (bc->bcAttributes.size() == 4 || bc->bcAttributes.size() == 6) {
1211  if (auto ext_disp_bc =
1212  dynamic_cast<DisplacementCubitBcDataWithRotation *>(
1213  bc->dispBcPtr.get())) {
1214  auto &o = ext_disp_bc->rotOffset;
1215  for (int a = 0; a != 3; ++a)
1216  o[a] = bc->bcAttributes[bc->bcAttributes.size() - 3 + a];
1217  MOFEM_LOG("BcMngWorld", Sev::inform)
1218  << "Add Rotate Z Center: " << o[0] << " " << o[1] << " "
1219  << o[2];
1220  }
1221  }
1222  } else if (std::regex_match(bc_id, std::regex("(.*)_ROTATE_ALL(.*)"))) {
1223  bc->dispBcPtr =
1224  boost::make_shared<DisplacementCubitBcDataWithRotation>();
1225  bc->dispBcPtr->data.flag4 = 1;
1226  bc->dispBcPtr->data.flag5 = 1;
1227  bc->dispBcPtr->data.flag6 = 1;
1228  if (bc->bcAttributes.size() >= 1) {
1229  bc->dispBcPtr->data.value4 = bc->bcAttributes[0];
1230  }
1231  if (bc->bcAttributes.size() >= 2) {
1232  bc->dispBcPtr->data.value5 = bc->bcAttributes[1];
1233  }
1234  if (bc->bcAttributes.size() >= 3) {
1235  bc->dispBcPtr->data.value6 = bc->bcAttributes[2];
1236  }
1237  MOFEM_LOG("BcMngWorld", Sev::inform) << "Add ALL " << bc_id;
1238  MOFEM_LOG("BcMngWorld", Sev::inform) << *(bc->dispBcPtr);
1239  if (bc->bcAttributes.size() > 3) {
1240  if (auto ext_disp_bc =
1241  dynamic_cast<DisplacementCubitBcDataWithRotation *>(
1242  bc->dispBcPtr.get())) {
1243  auto &o = ext_disp_bc->rotOffset;
1244  for (int a = 0; a != 3; ++a)
1245  o[a] = bc->bcAttributes[3 + a];
1246  MOFEM_LOG("BcMngWorld", Sev::inform)
1247  << "Add Rotate ALL Center: " << o[0] << " " << o[1] << " "
1248  << o[2];
1249  }
1250  }
1251  }
1252  }
1253  }
1254 
1255  MoFEMFunctionReturn(0);
1256 }
1257 
1258 template <>
1259 MoFEMErrorCode BcManager::pushMarkDOFsOnEntities<BcScalarMeshsetType<BLOCKSET>>(
1260  const std::string problem_name, const std::string block_name,
1261  const std::string field_name, bool get_low_dim_ents) {
1262  MoFEMFunctionBegin;
1263 
1264  CHKERR pushMarkDOFsOnEntities(problem_name, block_name, field_name, 0,
1265  MAX_DOFS_ON_ENTITY, get_low_dim_ents);
1266 
1267  auto regex_str =
1268  (boost::format("%s_%s_%s(.*)") % problem_name % field_name % block_name)
1269  .str();
1270 
1271  for (auto &m : bcMapByBlockName) {
1272 
1273  auto &bc_id = m.first;
1274 
1275  if (std::regex_match(bc_id, std::regex(regex_str))) {
1276 
1277  auto &bc = m.second;
1278  bc->tempBcPtr = boost::make_shared<TemperatureCubitBcData>();
1279  bc->tempBcPtr->data.flag1 = 1;
1280  if (bc->bcAttributes.empty()) {
1281  bc->tempBcPtr->data.value1 = 0;
1282  MOFEM_LOG("BcMngWorld", Sev::warning)
1283  << "Expected one attribute on block but have "
1284  << bc->bcAttributes.size();
1285  } else if (bc->bcAttributes.size() >= 1) {
1286  bc->tempBcPtr->data.value1 = bc->bcAttributes[0];
1287  }
1288  }
1289  }
1290 
1291  MoFEMFunctionReturn(0);
1292 }
1293 
1294 template <>
1295 MoFEMErrorCode BcManager::pushMarkDOFsOnEntities<DisplacementCubitBcData>(
1296  const std::string problem_name, const std::string field_name,
1297  bool get_low_dim_ents, bool block_name_field_prefix) {
1298  MoFEMFunctionBegin;
1299  // that marks DOFs and create data when are set by cubit nodesets.
1300  CHKERR pushMarkDOFsOnEntities<BcMeshsetType<DISPLACEMENTSET>>(
1301  problem_name, field_name, get_low_dim_ents, block_name_field_prefix);
1302  // that marks DOFs and create data when are set by blocsket.
1303  CHKERR pushMarkDOFsOnEntities<BcDisplacementMeshsetType<BLOCKSET>>(
1304  problem_name, field_name, get_low_dim_ents, block_name_field_prefix);
1305  MoFEMFunctionReturn(0);
1306 }
1307 
1308 template <>
1309 MoFEMErrorCode BcManager::removeBlockDOFsOnEntities<DisplacementCubitBcData>(
1310  const std::string problem_name, const std::string field_name,
1311  bool get_low_dim_ents, bool block_name_field_prefix,
1312  bool is_distributed_mesh) {
1313  MoFEMFunctionBegin;
1314  // that remove DOFs when are set by cubit nodesets.
1315  CHKERR removeBlockDOFsOnEntities<BcMeshsetType<DISPLACEMENTSET>>(
1316  problem_name, field_name, get_low_dim_ents, block_name_field_prefix,
1317  is_distributed_mesh);
1318  // that remove DOFs when are by blocksets
1319  CHKERR removeBlockDOFsOnEntities<BcDisplacementMeshsetType<BLOCKSET>>(
1320  problem_name, field_name, get_low_dim_ents, block_name_field_prefix,
1321  is_distributed_mesh);
1322  // add more ways to remove bcs when appropiate
1323  MoFEMFunctionReturn(0);
1324 }
1325 
1326 template <>
1327 MoFEMErrorCode BcManager::pushMarkDOFsOnEntities<TemperatureCubitBcData>(
1328  const std::string problem_name, const std::string field_name,
1329  bool get_low_dim_ents, bool block_name_field_prefix) {
1330  MoFEMFunctionBegin;
1331  CHKERR pushMarkDOFsOnEntities<BcMeshsetType<TEMPERATURESET>>(
1332  problem_name, field_name, get_low_dim_ents, block_name_field_prefix);
1333 
1334  auto get_block_name = [&]() {
1335  if (block_name_field_prefix)
1336  return (boost::format("%s_FIX_SCALAR") % field_name).str();
1337  else
1338  return field_name;
1339  };
1340 
1341  CHKERR pushMarkDOFsOnEntities<BcScalarMeshsetType<BLOCKSET>>(
1342  problem_name, get_block_name(), field_name, get_low_dim_ents);
1343  MoFEMFunctionReturn(0);
1344 }
1345 
1346 template <>
1347 MoFEMErrorCode BcManager::removeBlockDOFsOnEntities<TemperatureCubitBcData>(
1348  const std::string problem_name, const std::string field_name,
1349  bool get_low_dim_ents, bool block_name_field_prefix,
1350  bool is_distributed_mesh) {
1351  MoFEMFunctionBegin;
1352  CHKERR removeBlockDOFsOnEntities<BcMeshsetType<TEMPERATURESET>>(
1353  problem_name, field_name, get_low_dim_ents, block_name_field_prefix,
1354  is_distributed_mesh);
1355 
1356  auto get_block_name = [&]() {
1357  if (block_name_field_prefix)
1358  return (boost::format("%s_FIX_SCALAR") % field_name).str();
1359  else
1360  return field_name;
1361  };
1362 
1363  CHKERR removeBlockDOFsOnEntities<BcScalarMeshsetType<BLOCKSET>>(
1364  problem_name, get_block_name(), field_name, get_low_dim_ents,
1365  is_distributed_mesh);
1366  MoFEMFunctionReturn(0);
1367 }
1368 
1369 template <>
1370 MoFEMErrorCode BcManager::pushMarkDOFsOnEntities<HeatFluxCubitBcData>(
1371  const std::string problem_name, const std::string field_name,
1372  bool get_low_dim_ents, bool block_name_field_prefix) {
1373  MoFEMFunctionBegin;
1374  CHKERR pushMarkDOFsOnEntities<BcMeshsetType<HEATFLUXSET>>(
1375  problem_name, field_name, get_low_dim_ents, block_name_field_prefix);
1376  MoFEMFunctionReturn(0);
1377 }
1378 
1379 template <>
1380 MoFEMErrorCode BcManager::removeBlockDOFsOnEntities<HeatFluxCubitBcData>(
1381  const std::string problem_name, const std::string field_name,
1382  bool get_low_dim_ents, bool block_name_field_prefix,
1383  bool is_distributed_mesh) {
1384  MoFEMFunctionBegin;
1385  CHKERR removeBlockDOFsOnEntities<BcMeshsetType<HEATFLUXSET>>(
1386  problem_name, field_name, get_low_dim_ents, block_name_field_prefix,
1387  is_distributed_mesh);
1388  MoFEMFunctionReturn(0);
1389 }
1390 
1391 std::pair<std::string, std::string>
1392 BcManager::extractStringFromBlockId(const std::string block_id,
1393  const std::string prb_name) {
1394 
1395  // Assumes that field name is consist with letters and numbers.
1396  // No special characters.
1397  auto field_rgx_str =
1398  (boost::format("%s_([a-zA-Z0-9]*)_(.*)") % prb_name).str();
1399  std::regex field_rgx(field_rgx_str);
1400  std::smatch match_field_name;
1401  std::string field_name;
1402  std::string block_name;
1403 
1404  if (std::regex_search(block_id, match_field_name, field_rgx)) {
1405  field_name = match_field_name[1];
1406  block_name = match_field_name[2];
1407  } else {
1408  CHK_THROW_MESSAGE(MOFEM_DATA_INCONSISTENCY,
1409  "Field name and block name can not be resolved");
1410  }
1411 
1412  return std::make_pair(field_name, block_name);
1413 }
1414 
1415 template <>
1416 MoFEMErrorCode BcManager::pushMarkDOFsOnEntities<BcMeshsetType<FORCESET>>(
1417  const std::string problem_name, const std::string field_name,
1418  bool get_low_dim_ents, bool block_name_field_prefix) {
1419  Interface &m_field = cOre;
1420  auto prb_mng = m_field.getInterface<ProblemsManager>();
1421  MoFEMFunctionBegin;
1422 
1423  if (problem_name.empty())
1424  MOFEM_LOG("BcMngWorld", Sev::warning) << "Argument problem_name is empty";
1425 
1426  if (block_name_field_prefix)
1427  MOFEM_LOG("BcMngWorld", Sev::warning)
1428  << "Argument block_name_field_prefix=true has no effect";
1429 
1430  auto fix_force = [&]() {
1431  MoFEMFunctionBegin;
1432 
1433  auto iterate_meshsets = [&](auto &&meshset_vec_ptr) {
1434  MoFEMFunctionBegin;
1435  for (auto m : meshset_vec_ptr) {
1436  auto bc = boost::make_shared<BCs>();
1437  CHKERR m_field.get_moab().get_entities_by_handle(m->getMeshset(),
1438  bc->bcEnts, true);
1439  bc->forceBcPtr = boost::make_shared<ForceCubitBcData>();
1440  CHKERR m->getBcDataStructure(*(bc->forceBcPtr));
1441 
1442  MOFEM_LOG("BcMngWorld", Sev::verbose)
1443  << "Found block FORCESET id = " << m->getMeshsetId();
1444  MOFEM_LOG("BcMngWorld", Sev::verbose) << *bc->forceBcPtr;
1445 
1446  MOFEM_LOG("BcMngSync", Sev::noisy)
1447  << "Found block FORCESET id = " << m->getMeshsetId()
1448  << " nb. of entities " << bc->bcEnts.size()
1449  << " highest dim of entities "
1450  << BcManagerImplTools::get_dim(bc->bcEnts);
1451  MOFEM_LOG("BcMngSync", Sev::noisy) << *bc->forceBcPtr;
1452  MOFEM_LOG_SEVERITY_SYNC(m_field.get_comm(), Sev::noisy);
1453 
1454  if (problem_name.size()) {
1455 
1456  if (bc->forceBcPtr->data.value2 > 0)
1457  CHKERR prb_mng->modifyMarkDofs(problem_name, ROW, field_name, 0, 0,
1458  ProblemsManager::MarkOP::OR, 1,
1459  bc->bcMarkers);
1460  if (bc->forceBcPtr->data.value3 > 0)
1461  CHKERR prb_mng->modifyMarkDofs(problem_name, ROW, field_name, 1, 1,
1462  ProblemsManager::MarkOP::OR, 1,
1463  bc->bcMarkers);
1464  if (bc->forceBcPtr->data.value4 > 0)
1465  CHKERR prb_mng->modifyMarkDofs(problem_name, ROW, field_name, 2, 2,
1466  ProblemsManager::MarkOP::OR, 1,
1467  bc->bcMarkers);
1468 
1469  if (bc->forceBcPtr->data.value5 > 0) {
1470 
1471  CHKERR prb_mng->modifyMarkDofs(problem_name, ROW, field_name, 1, 1,
1472  ProblemsManager::MarkOP::OR, 1,
1473  bc->bcMarkers);
1474  CHKERR prb_mng->modifyMarkDofs(problem_name, ROW, field_name, 2, 2,
1475  ProblemsManager::MarkOP::OR, 1,
1476  bc->bcMarkers);
1477  }
1478  if (bc->forceBcPtr->data.value5) {
1479 
1480  CHKERR prb_mng->modifyMarkDofs(problem_name, ROW, field_name, 0, 0,
1481  ProblemsManager::MarkOP::OR, 1,
1482  bc->bcMarkers);
1483  CHKERR prb_mng->modifyMarkDofs(problem_name, ROW, field_name, 2, 2,
1484  ProblemsManager::MarkOP::OR, 1,
1485  bc->bcMarkers);
1486  }
1487  if (bc->forceBcPtr->data.value6) {
1488 
1489  CHKERR prb_mng->modifyMarkDofs(problem_name, ROW, field_name, 0, 0,
1490  ProblemsManager::MarkOP::OR, 1,
1491  bc->bcMarkers);
1492  CHKERR prb_mng->modifyMarkDofs(problem_name, ROW, field_name, 1, 1,
1493  ProblemsManager::MarkOP::OR, 1,
1494  bc->bcMarkers);
1495  }
1496  }
1497 
1498  if (get_low_dim_ents) {
1499  auto low_dim_ents =
1500  BcManagerImplTools::get_adj_ents(m_field.get_moab(), bc->bcEnts);
1501  bc->bcEnts.swap(low_dim_ents);
1502  }
1503 
1504  CHKERR m_field.getInterface<CommInterface>()->synchroniseEntities(
1505  bc->bcEnts);
1506  if (problem_name.size())
1507  CHKERR prb_mng->markDofs(problem_name, ROW, ProblemsManager::AND,
1508  bc->bcEnts, bc->bcMarkers);
1509 
1510  const std::string bc_id =
1511  problem_name + "_" + field_name + "_FORCESET" +
1512  boost::lexical_cast<std::string>(m->getMeshsetId());
1513  bcMapByBlockName[bc_id] = bc;
1514  }
1515  MoFEMFunctionReturn(0);
1516  };
1517 
1518  CHKERR iterate_meshsets(
1519 
1520  m_field.getInterface<MeshsetsManager>()->getCubitMeshsetPtr(NODESET |
1521  FORCESET)
1522 
1523  );
1524 
1525  MoFEMFunctionReturn(0);
1526  };
1527 
1528  CHKERR fix_force();
1529 
1530  MoFEMFunctionReturn(0);
1531 }
1532 
1533 template <>
1534 MoFEMErrorCode BcManager::pushMarkDOFsOnEntities<BcForceMeshsetType<BLOCKSET>>(
1535  const std::string problem_name, const std::string field_name,
1536  bool get_low_dim_ents, bool block_name_field_prefix) {
1537  Interface &m_field = cOre;
1538  auto prb_mng = m_field.getInterface<ProblemsManager>();
1539  MoFEMFunctionBegin;
1540 
1541  if (problem_name.size() == 0)
1542  MOFEM_LOG("BcMngWorld", Sev::warning) << "Argument problem_name is empty";
1543 
1544  auto get_force_block = [&](auto block_name) {
1545  MoFEMFunctionBegin;
1546 
1547  for (auto m :
1548  m_field.getInterface<MeshsetsManager>()->getCubitMeshsetPtr(std::regex(
1549 
1550  (boost::format("%s(.*)") % block_name).str()
1551 
1552  ))
1553 
1554  ) {
1555 
1556  const auto block_name = m->getName();
1557 
1558  MOFEM_LOG("BcMngWorld", Sev::inform)
1559  << "Found force block " << block_name;
1560 
1561  auto bc = boost::make_shared<BCs>();
1562  CHKERR m_field.get_moab().get_entities_by_handle(m->getMeshset(),
1563  bc->bcEnts, true);
1564 
1565  CHKERR m->getAttributes(bc->bcAttributes);
1566  if (bc->bcAttributes.size() != 3) {
1567  SETERRQ(m_field.get_comm(), MOFEM_DATA_INCONSISTENCY,
1568  "Expect three block attributes for force block");
1569  }
1570 
1571  bc->forceBcPtr = boost::make_shared<ForceCubitBcData>();
1572  // For details look at ForceCubitBcData in
1573  // mofem/src/multi_indices/BCData.hpp
1574  bc->forceBcPtr->data.value1 = 1;
1575  bc->forceBcPtr->data.value3 = bc->bcAttributes[0];
1576  bc->forceBcPtr->data.value4 = bc->bcAttributes[1];
1577  bc->forceBcPtr->data.value5 = bc->bcAttributes[2];
1578 
1579  MOFEM_LOG("BcMngWorld", Sev::inform) << *bc->forceBcPtr;
1580  MOFEM_LOG("BcMngSync", Sev::noisy)
1581  << "Found block FORCESET id = " << m->getMeshsetId()
1582  << " nb. of entities " << bc->bcEnts.size()
1583  << " highest dim of entities "
1584  << BcManagerImplTools::get_dim(bc->bcEnts);
1585  MOFEM_LOG("BcMngSync", Sev::noisy) << *bc->forceBcPtr;
1586  MOFEM_LOG_SEVERITY_SYNC(m_field.get_comm(), Sev::noisy);
1587 
1588  if (problem_name.size()) {
1589 
1590  if (bc->forceBcPtr->data.value2 > 0)
1591  CHKERR prb_mng->modifyMarkDofs(problem_name, ROW, field_name, 0, 0,
1592  ProblemsManager::MarkOP::OR, 1,
1593  bc->bcMarkers);
1594  if (bc->forceBcPtr->data.value3 > 0)
1595  CHKERR prb_mng->modifyMarkDofs(problem_name, ROW, field_name, 1, 1,
1596  ProblemsManager::MarkOP::OR, 1,
1597  bc->bcMarkers);
1598  if (bc->forceBcPtr->data.value4 > 0)
1599  CHKERR prb_mng->modifyMarkDofs(problem_name, ROW, field_name, 2, 2,
1600  ProblemsManager::MarkOP::OR, 1,
1601  bc->bcMarkers);
1602  }
1603 
1604  if (get_low_dim_ents) {
1605  auto low_dim_ents =
1606  BcManagerImplTools::get_adj_ents(m_field.get_moab(), bc->bcEnts);
1607  bc->bcEnts.swap(low_dim_ents);
1608  }
1609 
1610  CHKERR m_field.getInterface<CommInterface>()->synchroniseEntities(
1611  bc->bcEnts);
1612  if (problem_name.size())
1613  CHKERR prb_mng->markDofs(problem_name, ROW, ProblemsManager::AND,
1614  bc->bcEnts, bc->bcMarkers);
1615 
1616  const std::string bc_id =
1617  problem_name + "_" + field_name + "_" + block_name;
1618  bcMapByBlockName[bc_id] = bc;
1619  }
1620  MoFEMFunctionReturn(0);
1621  };
1622 
1623  CHKERR get_force_block("FORCE");
1624 
1625  MoFEMFunctionReturn(0);
1626 }
1627 
1628 template <>
1629 MoFEMErrorCode BcManager::pushMarkDOFsOnEntities<ForceCubitBcData>(
1630  const std::string problem_name, const std::string field_name,
1631  bool get_low_dim_ents, bool block_name_field_prefix) {
1632  MoFEMFunctionBegin;
1633 
1634  CHKERR pushMarkDOFsOnEntities<BcMeshsetType<FORCESET>>(
1635  problem_name, field_name, get_low_dim_ents, block_name_field_prefix);
1636 
1637  CHKERR pushMarkDOFsOnEntities<BcForceMeshsetType<BLOCKSET>>(
1638  problem_name, field_name, get_low_dim_ents, block_name_field_prefix);
1639 
1640  MoFEMFunctionReturn(0);
1641 }
1642 
1643 template <>
1644 MoFEMErrorCode BcManager::removeBlockDOFsOnEntities<BcMeshsetType<FORCESET>>(
1645  const std::string problem_name, const std::string field_name,
1646  bool get_low_dim_ents, bool block_name_field_prefix,
1647  bool is_distributed_mesh) {
1648  Interface &m_field = cOre;
1649  auto prb_mng = m_field.getInterface<ProblemsManager>();
1650  MoFEMFunctionBegin;
1651 
1652  CHKERR pushMarkDOFsOnEntities<BcMeshsetType<FORCESET>>(
1653  problem_name, field_name, get_low_dim_ents);
1654 
1655  std::array<Range, 3> ents_to_remove;
1656 
1657  for (auto m :
1658 
1659  m_field.getInterface<MeshsetsManager>()->getCubitMeshsetPtr(NODESET |
1660  FORCESET)) {
1661 
1662  const auto block_name = m->getName();
1663  std::string bc_id = problem_name + "_" + field_name + "_" + block_name;
1664 
1665  auto str = boost::format("%s_%s_%s(.*)")
1666 
1667  % problem_name % field_name % block_name;
1668 
1669  if (std::regex_match(bc_id, std::regex(str.str()))) {
1670 
1671  auto bc = bcMapByBlockName.at(bc_id);
1672 
1673  if (auto force_bc = bc->forceBcPtr) {
1674  if (force_bc->data.value3 > 0) {
1675  ents_to_remove[0].merge(bc->bcEnts);
1676  }
1677  if (force_bc->data.value4 > 0) {
1678  ents_to_remove[1].merge(bc->bcEnts);
1679  }
1680  if (force_bc->data.value5 > 0) {
1681  ents_to_remove[2].merge(bc->bcEnts);
1682  }
1683  } else {
1684  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1685  "BC type not implemented");
1686  }
1687  }
1688  }
1689 
1690  auto remove_dofs_on_ents = [&](const Range &ents, const int lo,
1691  const int hi) {
1692  if (is_distributed_mesh)
1693  return prb_mng->removeDofsOnEntities(problem_name, field_name, ents, lo,
1694  hi);
1695  else
1696  return prb_mng->removeDofsOnEntitiesNotDistributed(
1697  problem_name, field_name, ents, lo, hi);
1698  };
1699 
1700  CHKERR remove_dofs_on_ents(ents_to_remove[0], 0, 0);
1701  CHKERR remove_dofs_on_ents(ents_to_remove[1], 1, 1);
1702  CHKERR remove_dofs_on_ents(ents_to_remove[2], 2, 2);
1703 
1704  MoFEMFunctionReturn(0);
1705 }
1706 
1707 template <>
1708 MoFEMErrorCode
1709 BcManager::removeBlockDOFsOnEntities<BcForceMeshsetType<BLOCKSET>>(
1710  const std::string problem_name, const std::string field_name,
1711  bool get_low_dim_ents, bool block_name_field_prefix,
1712  bool is_distributed_mesh) {
1713  Interface &m_field = cOre;
1714  auto prb_mng = m_field.getInterface<ProblemsManager>();
1715  MoFEMFunctionBegin;
1716 
1717  CHKERR pushMarkDOFsOnEntities<BcForceMeshsetType<BLOCKSET>>(
1718  problem_name, field_name, get_low_dim_ents);
1719 
1720  std::array<Range, 3> ents_to_remove;
1721 
1722  for (auto m :
1723 
1724  m_field.getInterface<MeshsetsManager>()->getCubitMeshsetPtr(
1725  BLOCKSET | UNKNOWNNAME)) {
1726 
1727  const auto block_name = m->getName();
1728  std::string bc_id = problem_name + "_" + field_name + "_" + block_name;
1729 
1730  auto str = boost::format("%s_%s_%s(.*)")
1731 
1732  % problem_name % field_name % block_name;
1733 
1734  if (std::regex_match(bc_id, std::regex(str.str()))) {
1735 
1736  auto bc = bcMapByBlockName.at(bc_id);
1737 
1738  if (auto force_bc = bc->forceBcPtr) {
1739  if (force_bc->data.value3 > 0) {
1740  ents_to_remove[0].merge(bc->bcEnts);
1741  }
1742  if (force_bc->data.value4 > 0) {
1743  ents_to_remove[1].merge(bc->bcEnts);
1744  }
1745  if (force_bc->data.value5 > 0) {
1746  ents_to_remove[2].merge(bc->bcEnts);
1747  }
1748  } else {
1749  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1750  "BC type not implemented");
1751  }
1752  }
1753  }
1754 
1755  auto remove_dofs_on_ents = [&](const Range &ents, const int lo,
1756  const int hi) {
1757  if (is_distributed_mesh)
1758  return prb_mng->removeDofsOnEntities(problem_name, field_name, ents, lo,
1759  hi);
1760  else
1761  return prb_mng->removeDofsOnEntitiesNotDistributed(
1762  problem_name, field_name, ents, lo, hi);
1763  };
1764 
1765  CHKERR remove_dofs_on_ents(ents_to_remove[0], 0, 0);
1766  CHKERR remove_dofs_on_ents(ents_to_remove[1], 1, 1);
1767  CHKERR remove_dofs_on_ents(ents_to_remove[2], 2, 2);
1768 
1769  MoFEMFunctionReturn(0);
1770 }
1771 
1772 template <>
1773 MoFEMErrorCode BcManager::removeBlockDOFsOnEntities<ForceCubitBcData>(
1774  const std::string problem_name, const std::string field_name,
1775  bool get_low_dim_ents, bool block_name_field_prefix,
1776  bool is_distributed_mesh) {
1777  MoFEMFunctionBegin;
1778 
1779  CHKERR removeBlockDOFsOnEntities<BcMeshsetType<FORCESET>>(
1780  problem_name, field_name, get_low_dim_ents, block_name_field_prefix,
1781  is_distributed_mesh);
1782 
1783  CHKERR removeBlockDOFsOnEntities<BcForceMeshsetType<BLOCKSET>>(
1784  problem_name, field_name, get_low_dim_ents, block_name_field_prefix,
1785  is_distributed_mesh);
1786 
1787  MoFEMFunctionReturn(0);
1788 };
1789 
1790 MoFEMErrorCode BcManager::pushMarkSideDofs(
1791 
1792  const std::string problem_name, const std::string block_name,
1793  const std::string field_name, int bridge_dim, int lo, int hi
1794 
1795 ) {
1796  Interface &m_field = cOre;
1797  MoFEMFunctionBegin;
1798 
1799  if (problem_name.empty())
1800  MoFEMFunctionReturnHot(0);
1801 
1802  auto iterate_meshsets = [&](auto &&meshset_vec_ptr) {
1803  auto prb_mng = m_field.getInterface<ProblemsManager>();
1804  MoFEMFunctionBegin;
1805  for (auto m : meshset_vec_ptr) {
1806  auto bc = boost::make_shared<BCs>();
1807  CHKERR m_field.get_moab().get_entities_by_handle(m->getMeshset(),
1808  bc->bcEnts, true);
1809  CHKERR m->getAttributes(bc->bcAttributes);
1810 
1811  bc->dofsViewPtr = boost::make_shared<BCs::DofsView>();
1812 
1813  CHKERR m_field.getInterface<CommInterface>()->synchroniseEntities(
1814  bc->bcEnts);
1815  CHKERR prb_mng->getSideDofsOnBrokenSpaceEntities(
1816  *(bc->dofsViewPtr), problem_name, ROW, field_name, bc->bcEnts,
1817  bridge_dim, lo, hi);
1818  CHKERR prb_mng->markDofs(problem_name, ROW, *(bc->dofsViewPtr),
1819  ProblemsManager::OR, bc->bcMarkers);
1820 
1821  MOFEM_LOG("BcMngWorld", Sev::inform)
1822  << "Found block " << m->getName() << " number of attributes "
1823  << bc->bcAttributes.size() << " number of entities "
1824  << bc->bcEnts.size();
1825 
1826  const std::string bc_id =
1827  problem_name + "_" + field_name + "_" + m->getName();
1828  bcMapByBlockName[bc_id] = bc;
1829  }
1830  MoFEMFunctionReturn(0);
1831  };
1832 
1833  CHKERR iterate_meshsets(
1834 
1835  m_field.getInterface<MeshsetsManager>()->getCubitMeshsetPtr(std::regex(
1836 
1837  (boost::format("%s(.*)") % block_name).str()
1838 
1839  ))
1840 
1841  );
1842 
1843  MoFEMFunctionReturn(0);
1844 }
1845 
1846 MoFEMErrorCode BcManager::removeSideDOFs(const std::string problem_name,
1847  const std::string block_name,
1848  const std::string field_name,
1849  int bridge_dim, int lo, int hi,
1850  bool is_distributed_mesh) {
1851  Interface &m_field = cOre;
1852  MoFEMFunctionBegin;
1853 
1854  CHKERR pushMarkSideDofs(problem_name, block_name, field_name, bridge_dim, lo,
1855  hi);
1856 
1857  auto iterate_meshsets = [&](auto &&meshset_vec_ptr) {
1858  MoFEMFunctionBegin;
1859  auto prb_mng = m_field.getInterface<ProblemsManager>();
1860  for (auto m : meshset_vec_ptr) {
1861  const std::string bc_id =
1862  problem_name + "_" + field_name + "_" + m->getName();
1863  auto &bc = bcMapByBlockName.at(bc_id);
1864  CHKERR prb_mng->removeDofs(problem_name, ROW, *(bc->dofsViewPtr), lo, hi);
1865  CHKERR prb_mng->removeDofs(problem_name, COL, *(bc->dofsViewPtr), lo, hi);
1866  }
1867  MoFEMFunctionReturn(0);
1868  };
1869 
1870  if (is_distributed_mesh) {
1871 
1872  CHKERR iterate_meshsets(
1873 
1874  m_field.getInterface<MeshsetsManager>()->getCubitMeshsetPtr(std::regex(
1875 
1876  (boost::format("%s(.*)") % block_name).str()
1877 
1878  ))
1879 
1880  );
1881  } else {
1882  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "Not implemented");
1883  }
1884 
1885  MoFEMFunctionReturn(0);
1886 }
1887 
1888 } // namespace MoFEM
MoFEMFunctionReturnHot
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:460
MoFEM::UnknownInterface::getInterface
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface reference to pointer of interface.
Definition: UnknownInterface.hpp:93
CHK_MOAB_THROW
#define CHK_MOAB_THROW(err, msg)
Check error code of MoAB function and throw MoFEM exception.
Definition: definitions.h:589
MoFEM::BcManager::extractStringFromBlockId
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.
Definition: BcManager.cpp:1392
SIDESET
@ SIDESET
Definition: definitions.h:160
MoFEM::LogManager::checkIfChannelExist
static bool checkIfChannelExist(const std::string channel)
Check if channel exist.
Definition: LogManager.cpp:404
MoFEM::CoreTmp< 0 >
Core (interface) class.
Definition: Core.hpp:82
MoFEM::BcManager::popMarkDOFsOnEntities
boost::shared_ptr< BCs > popMarkDOFsOnEntities(const std::string block_name)
Get bc data and remove element.
Definition: BcManager.cpp:360
MOFEM_LOG_SEVERITY_SYNC
#define MOFEM_LOG_SEVERITY_SYNC(comm, severity)
Synchronise "SYNC" on curtain severity level.
Definition: LogManager.hpp:352
MoFEM::BcManager::query_interface
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
Definition: BcManager.cpp:34
MoFEM::ProblemsManager
Problem manager is used to build and partition problems.
Definition: ProblemsManager.hpp:21
MoFEM::BcManager::pushMarkSideDofs
MoFEMErrorCode pushMarkSideDofs(const std::string problem_name, const std::string block_name, const std::string field_name, int bridge_dim, int lo, int hi)
Mark side DOFs.
Definition: BcManager.cpp:1790
MoFEM::BcManager::BcManager
BcManager(const MoFEM::Core &core)
Definition: BcManager.cpp:41
MoFEM::CoreInterface::get_comm
virtual MPI_Comm & get_comm() const =0
MoFEM::BcManager::cOre
MoFEM::Core & cOre
Definition: BcManager.hpp:368
CHK_THROW_MESSAGE
#define CHK_THROW_MESSAGE(err, msg)
Check and throw MoFEM exception.
Definition: definitions.h:609
MoFEM::Exceptions::MoFEMErrorCode
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
MoFEM::CubitMeshSets
this struct keeps basic methods for moab meshset about material and boundary conditions
Definition: BCMultiIndices.hpp:19
MoFEM::BcManager::getOptions
MoFEMErrorCode getOptions()
get options
Definition: BcManager.cpp:65
MoFEM::BcManager::removeBlockDOFsOnEntities
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.
Definition: BcManager.cpp:73
MoFEM::BcManager::getMergedBlocksMarker
BcMarkerPtr getMergedBlocksMarker(std::vector< std::regex > bc_regex_vec)
Get the Merged Boundary Marker object.
Definition: BcManager.cpp:385
MoFEM::ProblemsManager::OR
@ OR
Definition: ProblemsManager.hpp:416
MoFEM::BcManager::removeSideDOFs
MoFEMErrorCode removeSideDOFs(const std::string problem_name, const std::string block_name, const std::string field_name, int bridge_dim, int lo, int hi, bool is_distributed_mesh=true)
Remove side DOFs.
Definition: BcManager.cpp:1846
MoFEM::LogManager::createSink
static boost::shared_ptr< SinkType > createSink(boost::shared_ptr< std::ostream > stream_ptr, std::string comm_filter)
Create a sink object.
Definition: LogManager.cpp:298
MoFEM::DeprecatedCoreInterface
Deprecated interface functions.
Definition: DeprecatedCoreInterface.hpp:16
UNKNOWNNAME
@ UNKNOWNNAME
Definition: definitions.h:171
ROW
@ ROW
Definition: definitions.h:136
MoFEM::Interface
DeprecatedCoreInterface Interface
Definition: Interface.hpp:2010
MoFEM::BcManager::pushMarkDOFsOnEntities
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:111
NODESET
@ NODESET
Definition: definitions.h:159
CHKERR
#define CHKERR
Inline error check.
Definition: definitions.h:548
MoFEM::CoreInterface::get_moab
virtual moab::Interface & get_moab()=0
MoFEM
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:10
MoFEM::ISManager
Section manager is used to create indexes and sections.
Definition: ISManager.hpp:23
a
constexpr double a
Definition: approx_sphere.cpp:30
MoFEM::BcManager
Simple interface for fast problem set-up.
Definition: BcManager.hpp:25
FORCESET
@ FORCESET
Definition: definitions.h:164
MoFEM::BcManagerImplTools::get_dim
auto get_dim(const Range &ents)
Definition: BcManager.cpp:10
MoFEM::BcManager::getBcMapByBlockName
BcMapByBlockName & getBcMapByBlockName()
Get the bc map.
Definition: BcManager.hpp:243
DISPLACEMENTSET
@ DISPLACEMENTSET
Definition: definitions.h:163
MoFEM::LogManager::getStrmSync
static boost::shared_ptr< std::ostream > getStrmSync()
Get the strm sync object.
Definition: LogManager.cpp:348
COL
@ COL
Definition: definitions.h:136
MoFEM::BcManager::getBlockIS
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:416
MoFEM::LogManager::getStrmWorld
static boost::shared_ptr< std::ostream > getStrmWorld()
Get the strm world object.
Definition: LogManager.cpp:344
MoFEM::LogManager::getStrmSelf
static boost::shared_ptr< std::ostream > getStrmSelf()
Get the strm self object.
Definition: LogManager.cpp:340
MoFEM::BcManager::getMergedBlocksRange
Range getMergedBlocksRange(std::vector< std::regex > bc_regex_vec)
Merge block ranges.
Definition: BcManager.cpp:370
MoFEM::BcManager::BcMarkerPtr
boost::shared_ptr< std::vector< char unsigned > > BcMarkerPtr
Definition: BcManager.hpp:227
MOFEM_LOG_TAG
#define MOFEM_LOG_TAG(channel, tag)
Tag channel.
Definition: LogManager.hpp:339
i
FTensor::Index< 'i', SPACE_DIM > i
Definition: hcurl_divergence_operator_2d.cpp:27
field_name
constexpr auto field_name
Definition: poisson_2d_homogeneous.cpp:13
MoFEM::ProblemsManager::AND
@ AND
Definition: ProblemsManager.hpp:416
MoFEM::BcManager::addBlockDOFsToMPCs
MoFEMErrorCode addBlockDOFsToMPCs(const std::string problem_name, const std::string field_name, bool get_low_dim_ents=false, bool block_name_field_prefix=false, bool is_distributed_mesh=false)
Definition: BcManager.cpp:185
MoFEM::UnknownInterface
base class for all interface classes
Definition: UnknownInterface.hpp:34
MAX_DOFS_ON_ENTITY
#define MAX_DOFS_ON_ENTITY
Maximal number of DOFs on entity.
Definition: definitions.h:249
Range
FTensor::dd
const Tensor2_symmetric_Expr< const ddTensor0< T, Dim, i, j >, typename promote< T, double >::V, Dim, i, j > dd(const Tensor0< T * > &a, const Index< i, Dim > index1, const Index< j, Dim > index2, const Tensor1< int, Dim > &d_ijk, const Tensor1< double, Dim > &d_xyz)
Definition: ddTensor0.hpp:33
MOFEM_LOG
#define MOFEM_LOG(channel, severity)
Log.
Definition: LogManager.hpp:308
MoFEM::CommInterface
Managing BitRefLevels.
Definition: CommInterface.hpp:21
BLOCKSET
@ BLOCKSET
Definition: definitions.h:161
MoFEM::BcManager::bcMapByBlockName
BcMapByBlockName bcMapByBlockName
Definition: BcManager.hpp:370
MoFEM::Exceptions::ierr
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Exceptions.hpp:76
MOFEM_DATA_INCONSISTENCY
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:31
TEMPERATURESET
@ TEMPERATURESET
Definition: definitions.h:168
MoFEM::MeshsetsManager
Interface for managing meshsets containing materials and boundary conditions.
Definition: MeshsetsManager.hpp:104
sdf_hertz_2d_axisymm_plane.d
float d
Definition: sdf_hertz_2d_axisymm_plane.py:4
m
FTensor::Index< 'm', 3 > m
Definition: shallow_wave.cpp:80
MoFEMFunctionBeginHot
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:453
HEATFLUXSET
@ HEATFLUXSET
Definition: definitions.h:169
MoFEM::DisplacementCubitBcDataWithRotation
A specialized version of DisplacementCubitBcData that includes an additional rotation offset.
Definition: EssentialDisplacementCubitBcData.hpp:24
MoFEM::MeshsetsManager::getCubitMeshsetPtr
MoFEMErrorCode getCubitMeshsetPtr(const int ms_id, const CubitBCType cubit_bc_type, const CubitMeshSets **cubit_meshset_ptr) const
get cubit meshset
Definition: MeshsetsManager.cpp:578
MoFEM::LogManager::setLog
static LoggerType & setLog(const std::string channel)
Set ans resset chanel logger.
Definition: LogManager.cpp:389
MoFEM::SmartPetscObj< IS >
MoFEM::MeshsetsManager::checkMeshset
bool checkMeshset(const int ms_id, const CubitBCType cubit_bc_type) const
check for CUBIT Id and CUBIT type
Definition: MeshsetsManager.cpp:360
MoFEMFunctionReturn
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:429
CHKERRG
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:496
MoFEMFunctionBegin
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359
l
FTensor::Index< 'l', 3 > l
Definition: matrix_function.cpp:21
MoFEM::BcManagerImplTools::get_adj_ents
auto get_adj_ents(moab::Interface &moab, const Range &ents)
Definition: BcManager.cpp:17
MoFEM::MPC::COUPLING
@ COUPLING
Generated by Doxygen 1.8.17 and hosted at University of Glasgow