v0.14.0
MeshRefinement.cpp
Go to the documentation of this file.
1 /** \file MeshRefinement.cpp
2  * \brief Mesh refinement interface
3  */
4 
5 
6 
7 #include <EntityRefine.hpp>
8 
9 namespace MoFEM {
10 
11 typedef multi_index_container<
12  boost::shared_ptr<RefElement>,
13  // ptrWrapperRefElement,
14  indexed_by<
15  ordered_unique<tag<Ent_mi_tag>,
16  const_mem_fun<RefElement::interface_type_RefEntity,
17  EntityHandle, &RefElement::getEnt>>,
18  ordered_non_unique<
19  tag<Ent_Ent_mi_tag>,
20  const_mem_fun<RefElement::interface_type_RefEntity, EntityHandle,
21  &RefElement::getParentEnt>>,
22  ordered_non_unique<
23  tag<Composite_ParentEnt_And_BitsOfRefinedEdges_mi_tag>,
24  composite_key<
25  RefElement,
26  const_mem_fun<RefElement::interface_type_RefEntity,
27  EntityHandle, &RefElement::getParentEnt>,
28  const_mem_fun<RefElement, int,
29  &RefElement::getBitRefEdgesUlong>>>>>
30  RefElement_multiIndex_parents_view;
31 
32 MoFEMErrorCode
33 MeshRefinement::query_interface(boost::typeindex::type_index type_index,
34  UnknownInterface **iface) const {
35  *iface = const_cast<MeshRefinement *>(this);
36  return 0;
37 }
38 
39 MeshRefinement::MeshRefinement(const Core &core)
40  : cOre(const_cast<Core &>(core)) {}
41 
42 MoFEMErrorCode MeshRefinement::addVerticesInTheMiddleOfEdges(
43  const EntityHandle meshset, const BitRefLevel &bit, const bool recursive,
44  int verb, EntityHandle start_v) {
45  Interface &m_field = cOre;
46  moab::Interface &moab = m_field.get_moab();
47  MoFEMFunctionBegin;
48  Range edges;
49  CHKERR moab.get_entities_by_type(meshset, MBEDGE, edges, recursive);
50  if (edges.empty()) {
51  Range tets;
52  CHKERR moab.get_entities_by_type(meshset, MBTET, tets, recursive);
53  CHKERR moab.get_adjacencies(tets, 1, true, edges, moab::Interface::UNION);
54  if (tets.empty()) {
55  Range prisms;
56  CHKERR moab.get_entities_by_type(meshset, MBPRISM, prisms, recursive);
57  for (Range::iterator pit = prisms.begin(); pit != prisms.end(); pit++) {
58  const EntityHandle *conn;
59  int num_nodes;
60  CHKERR moab.get_connectivity(*pit, conn, num_nodes, true);
61  assert(num_nodes == 6);
62  //
63  Range edge;
64  CHKERR moab.get_adjacencies(&conn[0], 2, 1, true, edge);
65  assert(edge.size() == 1);
66  edges.insert(edge[0]);
67  edge.clear();
68  CHKERR moab.get_adjacencies(&conn[1], 2, 1, true, edge);
69  assert(edge.size() == 1);
70  edges.insert(edge[0]);
71  EntityHandle conn_edge2[] = {conn[2], conn[0]};
72  edge.clear();
73  CHKERR moab.get_adjacencies(conn_edge2, 2, 1, true, edge);
74  assert(edge.size() == 1);
75  edges.insert(edge[0]);
76  //
77  edge.clear();
78  CHKERR moab.get_adjacencies(&conn[3], 2, 1, true, edge);
79  assert(edge.size() == 1);
80  edges.insert(edge[0]);
81  edge.clear();
82  CHKERR moab.get_adjacencies(&conn[4], 2, 1, true, edge);
83  assert(edge.size() == 1);
84  edges.insert(edge[0]);
85  EntityHandle conn_edge8[] = {conn[5], conn[3]};
86  edge.clear();
87  CHKERR moab.get_adjacencies(conn_edge8, 2, 1, true, edge);
88  assert(edge.size() == 1);
89  edges.insert(edge[0]);
90  }
91  }
92  }
93  CHKERR addVerticesInTheMiddleOfEdges(edges, bit, verb, start_v);
94  MoFEMFunctionReturn(0);
95 }
96 MoFEMErrorCode MeshRefinement::addVerticesInTheMiddleOfEdges(
97  const Range &ents, const BitRefLevel &bit, int verb, EntityHandle start_v) {
98  Interface &m_field = cOre;
99  moab::Interface &moab = m_field.get_moab();
100  auto refined_ents_ptr = m_field.get_ref_ents();
101  MoFEMFunctionBegin;
102  auto miit =
103  refined_ents_ptr->get<Composite_EntType_and_ParentEntType_mi_tag>()
104  .lower_bound(boost::make_tuple(MBVERTEX, MBEDGE));
105  auto hi_miit =
106  refined_ents_ptr->get<Composite_EntType_and_ParentEntType_mi_tag>()
107  .upper_bound(boost::make_tuple(MBVERTEX, MBEDGE));
108  RefEntity_multiIndex_view_by_ordered_parent_entity ref_parent_ents_view;
109  ref_parent_ents_view.insert(miit, hi_miit);
110 
111  Range edges = ents.subset_by_type(MBEDGE);
112  if (verb >= VERBOSE) {
113  MOFEM_TAG_AND_LOG("WORLD", Sev::verbose, "MeshRefinemnt")
114  << "ref level " << bit << " nb. edges to refine " << edges.size();
115  }
116 
117  std::array<std::vector<double>, 3> vert_coords;
118  for (auto &vc : vert_coords)
119  vc.reserve(edges.size());
120 
121  std::vector<EntityHandle> parent_edge;
122  parent_edge.reserve(edges.size());
123 
124  std::array<double, 9> coords;
125  FTensor::Tensor1<FTensor::PackPtr<double *, 0>, 3> t_0 = {
126  &coords[0], &coords[1], &coords[2]};
127  FTensor::Tensor1<FTensor::PackPtr<double *, 0>, 3> t_1 = {
128  &coords[3], &coords[4], &coords[5]};
129  FTensor::Tensor1<FTensor::PackPtr<double *, 0>, 3> t_2 = {
130  &coords[6], &coords[7], &coords[8]};
131 
132  FTensor::Index<'i', 3> i;
133 
134  Range add_bit;
135  for (auto p_eit = edges.pair_begin(); p_eit != edges.pair_end(); ++p_eit) {
136  auto miit_view = ref_parent_ents_view.lower_bound(p_eit->first);
137 
138  Range edge_having_parent_vertex;
139  if (miit_view != ref_parent_ents_view.end()) {
140  for (auto hi_miit_view = ref_parent_ents_view.upper_bound(p_eit->second);
141  miit_view != hi_miit_view; ++miit_view) {
142  edge_having_parent_vertex.insert(edge_having_parent_vertex.end(),
143  miit_view->get()->getParentEnt());
144  add_bit.insert(add_bit.end(), miit_view->get()->getEnt());
145  }
146  }
147 
148  Range add_vertex_edges =
149  subtract(Range(p_eit->first, p_eit->second), edge_having_parent_vertex);
150 
151  for (auto e : add_vertex_edges)
152  parent_edge.emplace_back(e);
153 
154  for (auto e : add_vertex_edges) {
155 
156  const EntityHandle *conn;
157  int num_nodes;
158  CHKERR moab.get_connectivity(e, conn, num_nodes, false);
159  if(num_nodes == 2) {
160  CHKERR moab.get_coords(conn, num_nodes, coords.data());
161  t_0(i) += t_1(i);
162  t_0(i) *= 0.5;
163  for (auto j : {0, 1, 2})
164  vert_coords[j].emplace_back(t_0(j));
165 
166  } else if (num_nodes == 3) {
167  for (auto j : {0, 1, 2})
168  vert_coords[j].emplace_back(t_2(j));
169  } else {
170  SETERRQ1(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
171  "edge should have 2 or 3 nodes but have %d", num_nodes);
172  }
173  }
174  }
175 
176  CHKERR m_field.getInterface<BitRefManager>()->addBitRefLevel(add_bit, bit);
177 
178  if (!vert_coords[0].empty()) {
179  ReadUtilIface *read_util;
180  CHKERR moab.query_interface(read_util);
181  int num_nodes = vert_coords[0].size();
182  vector<double *> arrays_coord;
183  CHKERR read_util->get_node_coords(3, num_nodes, 0, start_v, arrays_coord);
184  Range verts(start_v, start_v + num_nodes - 1);
185  for (auto dd : {0, 1, 2}) {
186  std::copy(vert_coords[dd].begin(), vert_coords[dd].end(),
187  arrays_coord[dd]);
188  }
189  CHKERR moab.tag_set_data(cOre.get_th_RefParentHandle(), verts,
190  &*parent_edge.begin());
191  CHKERR m_field.getInterface<BitRefManager>()->setEntitiesBitRefLevel(
192  verts, bit, verb);
193  }
194  MoFEMFunctionReturn(0);
195 }
196 
197 MoFEMErrorCode MeshRefinement::refineTets(const EntityHandle meshset,
198  const BitRefLevel &bit, int verb,
199  const bool debug) {
200  Interface &m_field = cOre;
201  moab::Interface &moab = m_field.get_moab();
202  MoFEMFunctionBegin;
203  Range tets;
204  CHKERR moab.get_entities_by_type(meshset, MBTET, tets, false);
205  CHKERR refineTets(tets, bit, verb, debug);
206  MoFEMFunctionReturn(0);
207 }
208 
209 MoFEMErrorCode MeshRefinement::SetParent::operator()(
210  const EntityHandle ent, const EntityHandle parent,
211  const RefEntity_multiIndex *ref_ents_ptr, MoFEM::Core &cOre) {
212  MoFEM::Interface &m_field = cOre;
213  MoFEMFunctionBegin;
214  auto it = ref_ents_ptr->find(ent);
215  if (it != ref_ents_ptr->end()) {
216  if (it->get()->getParentEnt() != parent && ent != parent) {
217  parentsToChange[ent] = parent;
218  }
219  } else {
220  if (ent != parent) {
221  CHKERR m_field.get_moab().tag_set_data(cOre.get_th_RefParentHandle(),
222  &ent, 1, &parent);
223  }
224  }
225  MoFEMFunctionReturn(0);
226 }
227 
228 MoFEMErrorCode MeshRefinement::SetParent::operator()(
229  const RefEntity_multiIndex *ref_ents_ptr) {
230  MoFEMFunctionBegin;
231  for (auto mit = parentsToChange.begin(); mit != parentsToChange.end();
232  ++mit) {
233  auto it = ref_ents_ptr->find(mit->first);
234  bool success = const_cast<RefEntity_multiIndex *>(ref_ents_ptr)
235  ->modify(it, RefEntity_change_parent(mit->second));
236  if (!success) {
237  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
238  "impossible to set parent");
239  }
240  }
241  MoFEMFunctionReturn(0);
242 }
243 
244 MoFEMErrorCode MeshRefinement::refineTets(const Range &_tets,
245  const BitRefLevel &bit, int verb,
246  const bool debug) {
247 
248  MoFEMFunctionBegin;
249 
250  auto set_edge_bits = [](moab::Interface &moab,
251  RefEntity_multiIndex_view_by_ordered_parent_entity
252  &ref_parent_ents_view,
253  EntityHandle tet, BitRefEdges &parent_edges_bit,
254  EntityHandle *edge_new_nodes, int *split_edges) {
255  MoFEMFunctionBeginHot;
256 
257  for (int ee = 0; ee != 6; ++ee) {
258  EntityHandle edge;
259  CHKERR moab.side_element(tet, 1, ee, edge);
260  auto eit = ref_parent_ents_view.find(edge);
261  if (eit != ref_parent_ents_view.end()) {
262  edge_new_nodes[ee] = (*eit)->getEnt();
263  split_edges[parent_edges_bit.count()] = ee;
264  parent_edges_bit.set(ee, 1);
265  }
266  }
267 
268  MoFEMFunctionReturnHot(0);
269  };
270 
271  CHKERR refineTets(_tets, bit, set_edge_bits, verb, debug);
272 
273  MoFEMFunctionReturn(0);
274 }
275 
276 MoFEMErrorCode MeshRefinement::refineTetsHangingNodes(const Range &_tets,
277  const BitRefLevel &bit,
278  int verb,
279  const bool debug) {
280 
281  MoFEMFunctionBegin;
282 
283  auto set_edge_bits = [](moab::Interface &moab,
284  RefEntity_multiIndex_view_by_ordered_parent_entity
285  &ref_parent_ents_view,
286  EntityHandle tet, BitRefEdges &parent_edges_bit,
287  EntityHandle *edge_new_nodes, int *split_edges) {
288  MoFEMFunctionBeginHot;
289 
290  for (int ee = 0; ee != 6; ++ee) {
291  EntityHandle edge;
292  CHKERR moab.side_element(tet, 1, ee, edge);
293  auto eit = ref_parent_ents_view.find(edge);
294  if (eit != ref_parent_ents_view.end()) {
295  edge_new_nodes[ee] = (*eit)->getEnt();
296  split_edges[parent_edges_bit.count()] = ee;
297  parent_edges_bit.set(ee, 1);
298  }
299  }
300 
301  if(parent_edges_bit.count() < 6)
302  parent_edges_bit.reset();
303 
304  MoFEMFunctionReturnHot(0);
305  };
306 
307  CHKERR refineTets(_tets, bit, set_edge_bits, verb, debug);
308 
309  MoFEMFunctionReturn(0);
310 }
311 
312 MoFEMErrorCode
313 MeshRefinement::refineTetsHangingNodes(const EntityHandle meshset,
314  const BitRefLevel &bit, int verb,
315  const bool debug) {
316  Interface &m_field = cOre;
317  moab::Interface &moab = m_field.get_moab();
318  MoFEMFunctionBegin;
319  Range tets;
320  CHKERR moab.get_entities_by_type(meshset, MBTET, tets, false);
321  CHKERR refineTetsHangingNodes(tets, bit, verb, debug);
322  MoFEMFunctionReturn(0);
323 }
324 
325 MoFEMErrorCode MeshRefinement::refineTets(const Range &_tets,
326  const BitRefLevel &bit,
327  SetEdgeBitsFun set_edge_bits,
328  int verb, const bool debug) {
329 
330  Interface &m_field = cOre;
331  moab::Interface &moab = m_field.get_moab();
332  auto refined_ents_ptr = m_field.get_ref_ents();
333  auto refined_finite_elements_ptr = m_field.get_ref_finite_elements();
334  ReadUtilIface *read_util;
335  MoFEMFunctionBegin;
336 
337  CHKERR m_field.get_moab().query_interface(read_util);
338 
339  Range ents_to_set_bit;
340 
341  // FIXME: refinement is based on entity handlers, should work on global ids of
342  // nodes, this will allow parallelise algorithm in the future
343 
344  auto get_parent_ents_view = [&](const auto type_child,
345  const auto type_parent) {
346  auto &ref_ents =
347  refined_ents_ptr->get<Composite_EntType_and_ParentEntType_mi_tag>();
348  auto range =
349  ref_ents.equal_range(boost::make_tuple(type_child, type_parent));
350 
351  RefEntity_multiIndex_view_by_ordered_parent_entity ref_parent_ents_view;
352 
353  using I = decltype(range.first);
354 
355  boost::function<bool(I it)> tester = [&](I it) {
356  return ((*it)->getBitRefLevel() & bit).any();
357  };
358 
359  boost::function<MoFEMErrorCode(I f, I s)> inserter = [&](I f, I s) {
360  ref_parent_ents_view.insert(f, s);
361  return 0;
362  };
363 
364  rangeInserter(range.first, range.second, tester, inserter);
365 
366  return ref_parent_ents_view;
367  };
368 
369  auto ref_parent_ents_view = get_parent_ents_view(MBVERTEX, MBEDGE);
370 
371  auto tets = _tets.subset_by_type(MBTET);
372 
373  auto get_ele_parent_view = [&]() {
374  auto &ref_ents = refined_finite_elements_ptr->get<Ent_mi_tag>();
375  RefElement_multiIndex_parents_view ref_ele_parent_view;
376  for (auto p = tets.pair_begin(); p != tets.pair_end(); ++p) {
377  auto tmi = ref_ents.lower_bound(p->first);
378  auto tme = ref_ents.upper_bound(p->second);
379  ref_ele_parent_view.insert(tmi, tme);
380  }
381  return ref_ele_parent_view;
382  };
383 
384  auto ref_ele_parent_view = get_ele_parent_view();
385 
386  std::vector<EntityHandle> parent_tets_refinded;
387  std::vector<EntityHandle> vertices_of_created_tets;
388  std::vector<BitRefEdges> parent_edges_bit_vec;
389  std::vector<int> nb_new_tets_vec;
390  std::vector<int> sub_type_vec;
391 
392  parent_tets_refinded.reserve(tets.size());
393  vertices_of_created_tets.reserve(4 * tets.size());
394  parent_edges_bit_vec.reserve(tets.size());
395  nb_new_tets_vec.reserve(tets.size());
396  sub_type_vec.reserve(tets.size());
397 
398  // make loop over all tets which going to be refined
399  for (auto tit : tets) {
400 
401  // get tet connectivity
402  const EntityHandle *conn;
403  int num_nodes;
404  CHKERR moab.get_connectivity(tit, conn, num_nodes, true);
405 
406  // get edges
407  BitRefEdges parent_edges_bit(0);
408  EntityHandle edge_new_nodes[] = {0, 0, 0, 0, 0, 0};
409  int split_edges[] = {-1, -1, -1, -1, -1, -1};
410 
411  CHKERR set_edge_bits(moab, ref_parent_ents_view, tit, parent_edges_bit,
412  edge_new_nodes, split_edges);
413 
414  if (parent_edges_bit.count()) {
415 
416  // swap nodes forward
417  EntityHandle _conn_[4];
418  std::copy(&conn[0], &conn[4], &_conn_[0]);
419 
420  // build connectivity for ref tets
421  EntityHandle new_tets_conns[8 * 4];
422  std::fill(&new_tets_conns[0], &new_tets_conns[8 * 4], no_handle);
423 
424  int sub_type = -1, nb_new_tets = 0;
425  switch (parent_edges_bit.count()) {
426  case 1:
427  sub_type = 0;
428  tet_type_1(_conn_, split_edges[0], edge_new_nodes[split_edges[0]],
429  new_tets_conns);
430  nb_new_tets = 2;
431  break;
432  case 2:
433  sub_type =
434  tet_type_2(_conn_, split_edges, edge_new_nodes, new_tets_conns);
435  if (sub_type & (4 | 8 | 16)) {
436  nb_new_tets = 3;
437  break;
438  } else if (sub_type == 1) {
439  nb_new_tets = 4;
440  break;
441  };
442  SETERRQ(PETSC_COMM_SELF, MOFEM_IMPOSSIBLE_CASE, "Impossible case");
443  break;
444  case 3:
445  sub_type =
446  tet_type_3(_conn_, split_edges, edge_new_nodes, new_tets_conns);
447  if (sub_type <= 4) {
448  nb_new_tets = 4;
449  break;
450  } else if (sub_type <= 7) {
451  nb_new_tets = 5;
452  break;
453  }
454  SETERRQ(PETSC_COMM_SELF, MOFEM_IMPOSSIBLE_CASE, "Impossible case");
455  case 4:
456  sub_type =
457  tet_type_4(_conn_, split_edges, edge_new_nodes, new_tets_conns);
458  if (sub_type == 0) {
459  nb_new_tets = 5;
460  break;
461  } else if (sub_type <= 7) {
462  nb_new_tets = 6;
463  break;
464  }
465  SETERRQ(PETSC_COMM_SELF, MOFEM_IMPOSSIBLE_CASE, "Impossible case");
466  case 5:
467  sub_type = tet_type_5(moab, _conn_, edge_new_nodes, new_tets_conns);
468  nb_new_tets = 7;
469  break;
470  case 6:
471  sub_type = 0;
472  tet_type_6(moab, _conn_, edge_new_nodes, new_tets_conns);
473  nb_new_tets = 8;
474  break;
475  default:
476  SETERRQ(PETSC_COMM_SELF, MOFEM_IMPOSSIBLE_CASE, "Impossible case");
477  }
478 
479  // find that tets
480  auto &ref_ele_by_parent_and_ref_edges =
481  ref_ele_parent_view
482  .get<Composite_ParentEnt_And_BitsOfRefinedEdges_mi_tag>();
483  auto it_by_ref_edges = ref_ele_by_parent_and_ref_edges.equal_range(
484  boost::make_tuple(tit, parent_edges_bit.to_ulong()));
485  // check if tet with this refinement scheme already exits
486  std::vector<EntityHandle> ents_to_set_bit_vec;
487  if (std::distance(it_by_ref_edges.first, it_by_ref_edges.second) ==
488  static_cast<size_t>(nb_new_tets)) {
489  ents_to_set_bit_vec.reserve(nb_new_tets);
490  for (int tt = 0; it_by_ref_edges.first != it_by_ref_edges.second;
491  it_by_ref_edges.first++, tt++) {
492  auto tet = (*it_by_ref_edges.first)->getEnt();
493  ents_to_set_bit_vec.emplace_back(tet);
494  // set ref tets entities
495  if (debug) {
496  // add this tet if exist to this ref level
497  auto ref_tet_it = refined_ents_ptr->find(tet);
498  if (ref_tet_it == refined_ents_ptr->end()) {
499  SETERRQ(m_field.get_comm(), MOFEM_DATA_INCONSISTENCY,
500  "Tetrahedron should be in database");
501  }
502  }
503  }
504  ents_to_set_bit.insert_list(ents_to_set_bit_vec.begin(),
505  ents_to_set_bit_vec.end());
506 
507  } else {
508  // if this element was not refined or was refined with different
509  // patterns of split edges create new elements
510 
511  parent_tets_refinded.emplace_back(tit);
512  for (int tt = 0; tt != nb_new_tets; ++tt) {
513  for (auto nn : {0, 1, 2, 3})
514  vertices_of_created_tets.emplace_back(new_tets_conns[4 * tt + nn]);
515  }
516  parent_edges_bit_vec.emplace_back(parent_edges_bit);
517  nb_new_tets_vec.emplace_back(nb_new_tets);
518  sub_type_vec.emplace_back(sub_type);
519  }
520  } else {
521 
522  if (debug) {
523  RefEntity_multiIndex::iterator tit_miit;
524  tit_miit = refined_ents_ptr->find(tit);
525  if (tit_miit == refined_ents_ptr->end())
526  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
527  "can not find this tet");
528  }
529 
530  ents_to_set_bit.insert(tit);
531  }
532  }
533 
534  // Create tets
535  EntityHandle start_e = 0;
536  EntityHandle *conn_moab;
537  const int nb_new_tets = vertices_of_created_tets.size() / 4;
538  read_util->get_element_connect(nb_new_tets, 4, MBTET, 0, start_e, conn_moab);
539  std::copy(vertices_of_created_tets.begin(), vertices_of_created_tets.end(),
540  conn_moab);
541  CHKERR read_util->update_adjacencies(start_e, nb_new_tets, 4, conn_moab);
542  ents_to_set_bit.insert(start_e, start_e + nb_new_tets - 1);
543 
544  // Create adj entities
545  Range ents;
546  for (auto d : {1, 2}) {
547  CHKERR moab.get_adjacencies(ents_to_set_bit, d, true, ents,
548  moab::Interface::UNION);
549  }
550 
551  SetParent set_parent;
552 
553  // Set parents and adjacencies
554  for (int idx = 0; idx != parent_tets_refinded.size(); ++idx) {
555 
556  const EntityHandle tit = parent_tets_refinded[idx];
557  const BitRefEdges &parent_edges_bit = parent_edges_bit_vec[idx];
558  const int nb_new_tets = nb_new_tets_vec[idx];
559 
560  std::array<EntityHandle, 8> ref_tets;
561  for (int tt = 0; tt != nb_new_tets; ++tt, ++start_e)
562  ref_tets[tt] = start_e;
563 
564  if (nb_new_tets) {
565 
566  CHKERR moab.tag_clear_data(cOre.get_th_RefBitEdge(), &ref_tets[0],
567  nb_new_tets, &parent_edges_bit);
568  CHKERR moab.tag_clear_data(cOre.get_th_RefParentHandle(), &ref_tets[0],
569  nb_new_tets, &tit);
570 
571  // hash map of nodes (RefEntity) by edges (EntityHandle)
572  std::map<EntityHandle /*edge*/, EntityHandle /*node*/>
573  map_ref_nodes_by_edges;
574 
575  Range tit_edges;
576  CHKERR moab.get_adjacencies(&tit, 1, 1, false, tit_edges);
577  for (auto edge : tit_edges) {
578  auto eit = ref_parent_ents_view.find(edge);
579  if (eit != ref_parent_ents_view.end()) {
580  map_ref_nodes_by_edges[(*eit)->getParentEnt()] = eit->get()->getEnt();
581  }
582  }
583 
584  // find parents for new edges and faces
585  // get tet edges and faces
586  Range tit_faces;
587  CHKERR moab.get_adjacencies(&tit, 1, 2, false, tit_faces);
588  Range edges_nodes[6], faces_nodes[4];
589  // for edges - add ref nodes
590  // edges_nodes[ee] - contains all nodes on edge ee including mid nodes if
591  // exist
592  Range::iterator eit = tit_edges.begin();
593  for (int ee = 0; eit != tit_edges.end(); eit++, ee++) {
594  CHKERR moab.get_connectivity(&*eit, 1, edges_nodes[ee], true);
595  auto map_miit = map_ref_nodes_by_edges.find(*eit);
596  if (map_miit != map_ref_nodes_by_edges.end()) {
597  edges_nodes[ee].insert(map_miit->second);
598  }
599  }
600  // for faces - add ref nodes
601  // faces_nodes[ff] - contains all nodes on face ff including mid nodes if
602  // exist
603  Range::iterator fit = tit_faces.begin();
604  for (int ff = 0; fit != tit_faces.end(); fit++, ff++) {
605  CHKERR moab.get_connectivity(&*fit, 1, faces_nodes[ff], true);
606  // Get edges on face and loop over those edges to add mid-nodes to range
607  Range fit_edges;
608  CHKERR moab.get_adjacencies(&*fit, 1, 1, false, fit_edges);
609  for (Range::iterator eit2 = fit_edges.begin(); eit2 != fit_edges.end();
610  eit2++) {
611  auto map_miit = map_ref_nodes_by_edges.find(*eit2);
612  if (map_miit != map_ref_nodes_by_edges.end()) {
613  faces_nodes[ff].insert(map_miit->second);
614  }
615  }
616  }
617  // add ref nodes to tet
618  // tet_nodes contains all nodes on tet including mid edge nodes
619  Range tet_nodes;
620  CHKERR moab.get_connectivity(&tit, 1, tet_nodes, true);
621  for (std::map<EntityHandle, EntityHandle>::iterator map_miit =
622  map_ref_nodes_by_edges.begin();
623  map_miit != map_ref_nodes_by_edges.end(); map_miit++) {
624  tet_nodes.insert(map_miit->second);
625  }
626  Range ref_edges;
627  // Get all all edges of refined tets
628  CHKERR moab.get_adjacencies(ref_tets.data(), nb_new_tets, 1, false,
629  ref_edges, moab::Interface::UNION);
630  // Check for all ref edge and set parents
631  for (Range::iterator reit = ref_edges.begin(); reit != ref_edges.end();
632  reit++) {
633  Range ref_edges_nodes;
634  CHKERR moab.get_connectivity(&*reit, 1, ref_edges_nodes, true);
635 
636  // Check if ref edge is an coarse edge (loop over coarse tet edges)
637  int ee = 0;
638  for (; ee < 6; ee++) {
639  // Two nodes are common (node[0],node[1],ref_node (if exist))
640  // this tests if given edge is contained by edge of refined
641  // tetrahedral
642  if (intersect(edges_nodes[ee], ref_edges_nodes).size() == 2) {
643  EntityHandle edge = tit_edges[ee];
644  CHKERR set_parent(*reit, edge, refined_ents_ptr, cOre);
645  break;
646  }
647  }
648  if (ee < 6)
649  continue; // this refined edge is contained by edge of tetrahedral
650  // check if ref edge is in coarse face
651  int ff = 0;
652  for (; ff < 4; ff++) {
653  // two nodes are common (node[0],node[1],ref_node (if exist))
654  // this tests if given face is contained by face of tetrahedral
655  if (intersect(faces_nodes[ff], ref_edges_nodes).size() == 2) {
656  EntityHandle face = tit_faces[ff];
657  CHKERR set_parent(*reit, face, refined_ents_ptr, cOre);
658  break;
659  }
660  }
661  if (ff < 4)
662  continue; // this refined edge is contained by face of tetrahedral
663 
664  // check if ref edge is in coarse tetrahedral (i.e. that is internal
665  // edge of refined tetrahedral)
666  if (intersect(tet_nodes, ref_edges_nodes).size() == 2) {
667  CHKERR set_parent(*reit, tit, refined_ents_ptr, cOre);
668  continue;
669  }
670 
671  // Refined edge is not child of any edge, face or tetrahedral, this is
672  // Impossible edge
673  SETERRQ(m_field.get_comm(), MOFEM_DATA_INCONSISTENCY,
674  "impossible refined edge");
675  }
676 
677  Range ref_faces;
678  CHKERR moab.get_adjacencies(ref_tets.data(), nb_new_tets, 2, false,
679  ref_faces, moab::Interface::UNION);
680  Tag th_interface_side;
681  const int def_side[] = {0};
682  CHKERR moab.tag_get_handle("INTERFACE_SIDE", 1, MB_TYPE_INTEGER,
683  th_interface_side,
684  MB_TAG_CREAT | MB_TAG_SPARSE, def_side);
685  // Check for all ref faces
686  for (auto rfit : ref_faces) {
687  Range ref_faces_nodes;
688  CHKERR moab.get_connectivity(&rfit, 1, ref_faces_nodes, true);
689  // Check if ref face is in coarse face
690  int ff = 0;
691  for (; ff < 4; ff++) {
692  // Check if refined triangle is contained by face of tetrahedral
693  if (intersect(faces_nodes[ff], ref_faces_nodes).size() == 3) {
694  EntityHandle face = tit_faces[ff];
695  CHKERR set_parent(rfit, face, refined_ents_ptr, cOre);
696  int side = 0;
697  // Set face side if it is on interface
698  CHKERR moab.tag_get_data(th_interface_side, &face, 1, &side);
699  CHKERR moab.tag_set_data(th_interface_side, &rfit, 1, &side);
700  break;
701  }
702  }
703  if (ff < 4)
704  continue; // this face is contained by one of tetrahedrons
705  // check if ref face is in coarse tetrahedral
706  // this is ref face which is contained by tetrahedral volume
707  if (intersect(tet_nodes, ref_faces_nodes).size() == 3) {
708  CHKERR set_parent(rfit, tit, refined_ents_ptr, cOre);
709  continue;
710  }
711  SETERRQ(m_field.get_comm(), MOFEM_DATA_INCONSISTENCY,
712  "impossible refined face");
713  }
714  }
715  }
716 
717  CHKERR set_parent(refined_ents_ptr);
718  CHKERR m_field.getInterface<BitRefManager>()->setBitRefLevel(
719  ents_to_set_bit, bit, true, verb, &ents);
720  MoFEMFunctionReturn(0);
721 }
722 MoFEMErrorCode MeshRefinement::refinePrisms(const EntityHandle meshset,
723  const BitRefLevel &bit, int verb) {
724 
725  Interface &m_field = cOre;
726  moab::Interface &moab = m_field.get_moab();
727  auto refined_ents_ptr = m_field.get_ref_ents();
728  auto refined_finite_elements_ptr = m_field.get_ref_finite_elements();
729 
730  // FIXME: refinement is based on entity handlers, should work on global ids of
731  // nodes, this will allow parallelise algorithm in the future
732 
733  MoFEMFunctionBegin;
734 
735  RefElement_multiIndex_parents_view ref_ele_parent_view;
736  ref_ele_parent_view.insert(
737  refined_finite_elements_ptr->get<Ent_mi_tag>().lower_bound(
738  get_id_for_min_type<MBPRISM>()),
739  refined_finite_elements_ptr->get<Ent_mi_tag>().upper_bound(
740  get_id_for_max_type<MBPRISM>()));
741  auto &ref_ele_by_parent_and_ref_edges =
742  ref_ele_parent_view
743  .get<Composite_ParentEnt_And_BitsOfRefinedEdges_mi_tag>();
744  // find all vertices which parent is edge
745  auto &ref_ents_by_comp =
746  refined_ents_ptr->get<Composite_EntType_and_ParentEntType_mi_tag>();
747  RefEntity_multiIndex_view_by_hashed_parent_entity ref_parent_ents_view;
748  ref_parent_ents_view.insert(
749  ref_ents_by_comp.lower_bound(boost::make_tuple(MBVERTEX, MBEDGE)),
750  ref_ents_by_comp.upper_bound(boost::make_tuple(MBVERTEX, MBEDGE)));
751  Range prisms;
752  CHKERR moab.get_entities_by_type(meshset, MBPRISM, prisms, false);
753  Range::iterator pit = prisms.begin();
754  for (; pit != prisms.end(); pit++) {
755  auto miit_prism = refined_ents_ptr->get<Ent_mi_tag>().find(*pit);
756  if (miit_prism == refined_ents_ptr->end()) {
757  SETERRQ(m_field.get_comm(), MOFEM_DATA_INCONSISTENCY,
758  "this prism is not in ref database");
759  }
760  if (verb >= NOISY) {
761  MOFEM_TAG_AND_LOG("WORLD", Sev::noisy, "MeshRefinment")
762  << "ref prism " << **miit_prism;
763  }
764  // prism connectivity
765  int num_nodes;
766  const EntityHandle *conn;
767  CHKERR moab.get_connectivity(*pit, conn, num_nodes, true);
768  assert(num_nodes == 6);
769  // edges connectivity
770  EntityHandle edges[6];
771  for (int ee = 0; ee < 3; ee++) {
772  CHKERR moab.side_element(*pit, 1, ee, edges[ee]);
773  }
774  for (int ee = 6; ee < 9; ee++) {
775  CHKERR moab.side_element(*pit, 1, ee, edges[ee - 3]);
776  }
777  // detect split edges
778  BitRefEdges split_edges(0);
779  EntityHandle edge_nodes[6];
780  std::fill(&edge_nodes[0], &edge_nodes[6], no_handle);
781  for (int ee = 0; ee < 6; ee++) {
782  RefEntity_multiIndex_view_by_hashed_parent_entity::iterator miit_view =
783  ref_parent_ents_view.find(edges[ee]);
784  if (miit_view != ref_parent_ents_view.end()) {
785  if (((*miit_view)->getBitRefLevel() & bit).any()) {
786  edge_nodes[ee] = (*miit_view)->getEnt();
787  split_edges.set(ee);
788  }
789  }
790  }
791  if (split_edges.count() == 0) {
792  *(const_cast<RefEntity *>(miit_prism->get())->getBitRefLevelPtr()) |= bit;
793  if (verb >= VERY_NOISY)
794  MOFEM_TAG_AND_LOG("WORLD", Sev::noisy, "MeshRefinment")
795  << "no refinement";
796  continue;
797  }
798  // check consistency
799  if (verb >= NOISY) {
800  std::ostringstream ss;
801  ss << "prism split edges " << split_edges << " count "
802  << split_edges.count() << std::endl;
803  PetscPrintf(m_field.get_comm(), ss.str().c_str());
804  }
805  // prism ref
806  EntityHandle new_prism_conn[4 * 6];
807  std::fill(&new_prism_conn[0], &new_prism_conn[4 * 6], no_handle);
808  int nb_new_prisms = 0;
809  switch (split_edges.count()) {
810  case 0:
811  break;
812  case 2:
813  CHKERR prism_type_1(conn, split_edges, edge_nodes, new_prism_conn);
814  nb_new_prisms = 2;
815  break;
816  case 4:
817  CHKERR prism_type_2(conn, split_edges, edge_nodes, new_prism_conn);
818  nb_new_prisms = 3;
819  break;
820  case 6:
821  CHKERR prism_type_3(conn, split_edges, edge_nodes, new_prism_conn);
822  nb_new_prisms = 4;
823  break;
824  default:
825  std::ostringstream ss;
826  ss << split_edges << " : [ " << conn[0] << " " << conn[1] << " "
827  << conn[2] << " " << conn[3] << " " << conn[4] << " " << conn[5]
828  << " ]";
829  SETERRQ(m_field.get_comm(), MOFEM_DATA_INCONSISTENCY, ss.str().c_str());
830  }
831  // find that prism
832  std::bitset<4> ref_prism_bit(0);
833  auto it_by_ref_edges = ref_ele_by_parent_and_ref_edges.lower_bound(
834  boost::make_tuple(*pit, split_edges.to_ulong()));
835  auto hi_it_by_ref_edges = ref_ele_by_parent_and_ref_edges.upper_bound(
836  boost::make_tuple(*pit, split_edges.to_ulong()));
837  auto it_by_ref_edges2 = it_by_ref_edges;
838  for (int pp = 0; it_by_ref_edges2 != hi_it_by_ref_edges;
839  it_by_ref_edges2++, pp++) {
840  // Add this tet to this ref
841  *(const_cast<RefElement *>(it_by_ref_edges2->get())
842  ->getBitRefLevelPtr()) |= bit;
843  ref_prism_bit.set(pp, 1);
844  if (verb > 2) {
845  std::ostringstream ss;
846  ss << "is refined " << *it_by_ref_edges2 << std::endl;
847  PetscPrintf(m_field.get_comm(), ss.str().c_str());
848  }
849  }
850  if (it_by_ref_edges != hi_it_by_ref_edges) {
851  if (ref_prism_bit.count() != (unsigned int)nb_new_prisms)
852  SETERRQ(m_field.get_comm(), MOFEM_DATA_INCONSISTENCY,
853  "data inconsistency");
854  } else {
855  EntityHandle ref_prisms[4];
856  // create prism
857  for (int pp = 0; pp < nb_new_prisms; pp++) {
858  if (verb > 3) {
859  std::ostringstream ss;
860  ss << "ref prism " << ref_prism_bit << std::endl;
861  PetscPrintf(m_field.get_comm(), ss.str().c_str());
862  }
863  if (!ref_prism_bit.test(pp)) {
864  CHKERR moab.create_element(MBPRISM, &new_prism_conn[6 * pp], 6,
865  ref_prisms[pp]);
866  CHKERR moab.tag_set_data(cOre.get_th_RefParentHandle(),
867  &ref_prisms[pp], 1, &*pit);
868  CHKERR moab.tag_set_data(cOre.get_th_RefBitLevel(), &ref_prisms[pp],
869  1, &bit);
870  CHKERR moab.tag_set_data(cOre.get_th_RefBitEdge(), &ref_prisms[pp], 1,
871  &split_edges);
872  auto p_ent =
873  const_cast<RefEntity_multiIndex *>(refined_ents_ptr)
874  ->insert(boost::shared_ptr<RefEntity>(new RefEntity(
875  m_field.get_basic_entity_data_ptr(), ref_prisms[pp])));
876  std::pair<RefElement_multiIndex::iterator, bool> p_fe;
877  try {
878  p_fe =
879  const_cast<RefElement_multiIndex *>(refined_finite_elements_ptr)
880  ->insert(boost::shared_ptr<RefElement>(
881  new RefElement_PRISM(*p_ent.first)));
882  } catch (MoFEMException const &e) {
883  SETERRQ(PETSC_COMM_SELF, e.errorCode, e.errorMessage);
884  }
885  ref_prism_bit.set(pp);
886  CHKERR cOre.addPrismToDatabase(ref_prisms[pp]);
887  if (verb > 2) {
888  std::ostringstream ss;
889  ss << "add prism: " << **p_fe.first << std::endl;
890  if (verb > 7) {
891  for (int nn = 0; nn < 6; nn++) {
892  ss << new_prism_conn[nn] << " ";
893  }
894  ss << std::endl;
895  }
896  PetscPrintf(m_field.get_comm(), ss.str().c_str());
897  }
898  }
899  }
900  }
901  }
902  MoFEMFunctionReturn(0);
903 }
904 MoFEMErrorCode MeshRefinement::refineMeshset(const EntityHandle meshset,
905  const BitRefLevel &bit,
906  const bool recursive, int verb) {
907  Interface &m_field = cOre;
908  auto refined_ents_ptr = m_field.get_ref_ents();
909  MoFEMFunctionBegin;
910  auto miit = refined_ents_ptr->find(meshset);
911  if (miit == refined_ents_ptr->end()) {
912  SETERRQ(m_field.get_comm(), MOFEM_DATA_INCONSISTENCY,
913  "this meshset is not in ref database");
914  }
915  CHKERR m_field.getInterface<BitRefManager>()->updateMeshsetByEntitiesChildren(
916  meshset, bit, meshset, MBMAXTYPE, recursive, verb);
917  *(const_cast<RefEntity *>(miit->get())->getBitRefLevelPtr()) |= bit;
918  MoFEMFunctionReturn(0);
919 }
920 
921 MoFEMErrorCode MeshRefinement::refineTris(const EntityHandle meshset,
922  const BitRefLevel &bit, int verb,
923  const bool debug) {
924  Interface &m_field = cOre;
925  moab::Interface &moab = m_field.get_moab();
926  MoFEMFunctionBegin;
927  Range tris;
928  CHKERR moab.get_entities_by_type(meshset, MBTRI, tris, false);
929  CHKERR refineTris(tris, bit, verb, debug);
930  MoFEMFunctionReturn(0);
931 }
932 
933 MoFEMErrorCode MeshRefinement::refineTris(const Range &_tris,
934  const BitRefLevel &bit, int verb,
935  const bool debug) {
936  MoFEMFunctionBegin;
937 
938  auto set_edge_bits = [](moab::Interface &moab,
939  RefEntity_multiIndex_view_by_ordered_parent_entity
940  &ref_parent_ents_view,
941  EntityHandle tet, BitRefEdges &parent_edges_bit,
942  EntityHandle *edge_new_nodes, int *split_edges) {
943  MoFEMFunctionBeginHot;
944 
945  for (int ee = 0; ee != 3; ++ee) {
946  EntityHandle edge;
947  CHKERR moab.side_element(tet, 1, ee, edge);
948  auto eit = ref_parent_ents_view.find(edge);
949  if (eit != ref_parent_ents_view.end()) {
950  edge_new_nodes[ee] = (*eit)->getEnt();
951  split_edges[parent_edges_bit.count()] = ee;
952  parent_edges_bit.set(ee, 1);
953  }
954  }
955 
956  MoFEMFunctionReturnHot(0);
957  };
958 
959  CHKERR refineTris(_tris, bit, set_edge_bits, verb, debug);
960 
961  MoFEMFunctionReturn(0);
962 }
963 
964 MoFEMErrorCode MeshRefinement::refineTrisHangingNodes(const Range &_tris,
965  const BitRefLevel &bit,
966  int verb,
967  const bool debug) {
968 
969  MoFEMFunctionBegin;
970 
971  auto set_edge_bits = [](moab::Interface &moab,
972  RefEntity_multiIndex_view_by_ordered_parent_entity
973  &ref_parent_ents_view,
974  EntityHandle tet, BitRefEdges &parent_edges_bit,
975  EntityHandle *edge_new_nodes, int *split_edges) {
976  MoFEMFunctionBeginHot;
977 
978  for (int ee = 0; ee != 3; ++ee) {
979  EntityHandle edge;
980  CHKERR moab.side_element(tet, 1, ee, edge);
981  auto eit = ref_parent_ents_view.find(edge);
982  if (eit != ref_parent_ents_view.end()) {
983  edge_new_nodes[ee] = (*eit)->getEnt();
984  split_edges[parent_edges_bit.count()] = ee;
985  parent_edges_bit.set(ee, 1);
986  }
987  }
988 
989  if (parent_edges_bit.count() < 3) {
990  parent_edges_bit.reset();
991  }
992 
993  MoFEMFunctionReturnHot(0);
994  };
995 
996  CHKERR refineTris(_tris, bit, set_edge_bits, verb, debug);
997 
998  MoFEMFunctionReturn(0);
999 };
1000 
1001 MoFEMErrorCode
1002 MeshRefinement::refineTrisHangingNodes(const EntityHandle meshset,
1003  const BitRefLevel &bit, int verb,
1004  const bool debug) {
1005  Interface &m_field = cOre;
1006  moab::Interface &moab = m_field.get_moab();
1007  MoFEMFunctionBegin;
1008  Range tris;
1009  CHKERR moab.get_entities_by_type(meshset, MBTRI, tris, false);
1010  CHKERR refineTrisHangingNodes(tris, bit, verb, debug);
1011  MoFEMFunctionReturn(0);
1012 }
1013 
1014 MoFEMErrorCode MeshRefinement::refineTris(const Range &_tris,
1015  const BitRefLevel &bit,
1016  SetEdgeBitsFun set_edge_bits,
1017  int verb, const bool debug) {
1018  Interface &m_field = cOre;
1019  moab::Interface &moab = m_field.get_moab();
1020  auto refined_ents_ptr = m_field.get_ref_ents();
1021  auto refined_finite_elements_ptr = m_field.get_ref_finite_elements();
1022  ReadUtilIface *read_util;
1023 
1024  MoFEMFunctionBegin;
1025 
1026  auto get_parent_ents_view = [&](const auto type_child,
1027  const auto type_parent) {
1028  auto &ref_ents =
1029  refined_ents_ptr->get<Composite_EntType_and_ParentEntType_mi_tag>();
1030  auto range =
1031  ref_ents.equal_range(boost::make_tuple(type_child, type_parent));
1032 
1033  RefEntity_multiIndex_view_by_ordered_parent_entity ref_parent_ents_view;
1034 
1035  using I = decltype(range.first);
1036 
1037  boost::function<bool(I it)> tester = [&](I it) {
1038  return ((*it)->getBitRefLevel() & bit).any();
1039  };
1040 
1041  boost::function<MoFEMErrorCode(I f, I s)> inserter = [&](I f, I s) {
1042  ref_parent_ents_view.insert(f, s);
1043  return 0;
1044  };
1045 
1046  rangeInserter(range.first, range.second, tester, inserter);
1047 
1048  return ref_parent_ents_view;
1049  };
1050 
1051  auto ref_parent_ents_view = get_parent_ents_view(MBVERTEX, MBEDGE);
1052 
1053  auto tris = _tris.subset_by_type(MBTRI);
1054 
1055  auto get_ele_parent_view = [&]() {
1056  auto &ref_ents = refined_finite_elements_ptr->get<Ent_mi_tag>();
1057  RefElement_multiIndex_parents_view ref_ele_parent_view;
1058  for (auto p = tris.pair_begin(); p != tris.pair_end(); ++p) {
1059  auto tmi = ref_ents.lower_bound(p->first);
1060  auto tme = ref_ents.upper_bound(p->second);
1061  ref_ele_parent_view.insert(tmi, tme);
1062  }
1063  return ref_ele_parent_view;
1064  };
1065 
1066  auto ref_ele_parent_view = get_ele_parent_view();
1067 
1068  Range ents_to_set_bit;
1069 
1070  std::vector<EntityHandle> parent_tris_refined; // entities refined
1071  std::vector<EntityHandle> vertices_of_created_tris;
1072  std::vector<BitRefEdges> parent_edges_bit_vec;
1073  std::vector<int> nb_new_tris_vec;
1074 
1075  parent_tris_refined.reserve(tris.size());
1076  vertices_of_created_tris.reserve(3 * tris.size());
1077  parent_edges_bit_vec.reserve(tris.size());
1078  nb_new_tris_vec.reserve(tris.size());
1079 
1080  // make loop over all tets which going to be refined
1081  for (auto tit : tris) {
1082 
1083  // get tet connectivity
1084  const EntityHandle *conn;
1085  int num_nodes;
1086  CHKERR moab.get_connectivity(tit, conn, num_nodes, true);
1087 
1088  // get edges
1089  BitRefEdges parent_edges_bit(0);
1090  EntityHandle edge_new_nodes[] = {0, 0, 0};
1091  int split_edges[] = {-1, -1, -1};
1092 
1093  CHKERR set_edge_bits(moab, ref_parent_ents_view, tit, parent_edges_bit,
1094  edge_new_nodes, split_edges);
1095 
1096  if (parent_edges_bit.count()) {
1097 
1098  // swap nodes forward
1099  EntityHandle _conn_[3];
1100  std::copy(&conn[0], &conn[3], &_conn_[0]);
1101 
1102  // build connectivity for ref tets
1103  EntityHandle new_tets_conns[4 * 3];
1104  std::fill(&new_tets_conns[0], &new_tets_conns[4 * 3], no_handle);
1105 
1106  int nb_new_tris = 0;
1107  switch (parent_edges_bit.count()) {
1108  case 1:
1109  CHKERR tri_type_1(conn, split_edges[0], edge_new_nodes[split_edges[0]],
1110  new_tets_conns);
1111  nb_new_tris += 2;
1112  break;
1113  case 2:
1114  CHKERR tri_type_2(conn, split_edges, edge_new_nodes, new_tets_conns);
1115  nb_new_tris += 3;
1116  break;
1117  case 3:
1118  CHKERR tri_type_3(conn, edge_new_nodes, new_tets_conns);
1119  nb_new_tris += 4;
1120  break;
1121  default:
1122  SETERRQ(PETSC_COMM_SELF, MOFEM_IMPOSSIBLE_CASE, "Impossible case");
1123  }
1124 
1125  // find that tets
1126  auto &ref_ele_by_parent_and_ref_edges =
1127  ref_ele_parent_view
1128  .get<Composite_ParentEnt_And_BitsOfRefinedEdges_mi_tag>();
1129  auto it_by_ref_edges = ref_ele_by_parent_and_ref_edges.equal_range(
1130  boost::make_tuple(tit, parent_edges_bit.to_ulong()));
1131  // check if tet with this refinement scheme already exits
1132  std::vector<EntityHandle> ents_to_set_bit_vec;
1133  if (std::distance(it_by_ref_edges.first, it_by_ref_edges.second) ==
1134  static_cast<size_t>(nb_new_tris)) {
1135  ents_to_set_bit_vec.reserve(nb_new_tris);
1136  for (int tt = 0; it_by_ref_edges.first != it_by_ref_edges.second;
1137  it_by_ref_edges.first++, tt++) {
1138  auto tet = (*it_by_ref_edges.first)->getEnt();
1139  ents_to_set_bit_vec.emplace_back(tet);
1140  // set ref tets entities
1141  if (debug) {
1142  // add this tet if exist to this ref level
1143  auto ref_tet_it = refined_ents_ptr->find(tet);
1144  if (ref_tet_it == refined_ents_ptr->end()) {
1145  SETERRQ(m_field.get_comm(), MOFEM_DATA_INCONSISTENCY,
1146  "Tetrahedron should be in database");
1147  }
1148  }
1149  }
1150  ents_to_set_bit.insert_list(ents_to_set_bit_vec.begin(),
1151  ents_to_set_bit_vec.end());
1152 
1153  } else {
1154  // if this element was not refined or was refined with different
1155  // patterns of split edges create new elements
1156 
1157  parent_tris_refined.emplace_back(tit);
1158  for (int tt = 0; tt != nb_new_tris; ++tt) {
1159  for (auto nn : {0, 1, 2})
1160  vertices_of_created_tris.emplace_back(new_tets_conns[3 * tt + nn]);
1161  }
1162  parent_edges_bit_vec.emplace_back(parent_edges_bit);
1163  nb_new_tris_vec.emplace_back(nb_new_tris);
1164  }
1165  } else {
1166 
1167  if (debug) {
1168  RefEntity_multiIndex::iterator tit_miit;
1169  tit_miit = refined_ents_ptr->find(tit);
1170  if (tit_miit == refined_ents_ptr->end())
1171  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1172  "can not find this tet");
1173  }
1174 
1175  ents_to_set_bit.insert(tit);
1176  }
1177  }
1178 
1179  const int nb_new_tris = vertices_of_created_tris.size() / 3;
1180  EntityHandle start_e = 0;
1181 
1182  if (nb_new_tris) {
1183 
1184  // Create tets
1185  EntityHandle *conn_moab;
1186  CHKERR m_field.get_moab().query_interface(read_util);
1187  CHKERR read_util->get_element_connect(nb_new_tris, 3, MBTRI, 0, start_e,
1188  conn_moab);
1189  std::copy(vertices_of_created_tris.begin(), vertices_of_created_tris.end(),
1190  conn_moab);
1191  CHKERR read_util->update_adjacencies(start_e, nb_new_tris, 3, conn_moab);
1192  ents_to_set_bit.insert(start_e, start_e + nb_new_tris - 1);
1193 
1194  if (debug) {
1195  auto meshset_ptr = get_temp_meshset_ptr(moab);
1196  CHKERR moab.add_entities(*meshset_ptr, ents_to_set_bit);
1197  CHKERR moab.write_file("out_refinded_debug.vtk", "VTK", "",
1198  meshset_ptr->get_ptr(), 1);
1199  }
1200  }
1201 
1202  // Create adj entities
1203  Range ents;
1204  for (auto d : {1}) {
1205  CHKERR moab.get_adjacencies(ents_to_set_bit, d, true, ents,
1206  moab::Interface::UNION);
1207  }
1208 
1209  SetParent set_parent;
1210 
1211  // Set parrents and adjacencies. Iterate over all parent triangles.
1212  for (int idx = 0; idx != parent_tris_refined.size(); ++idx) {
1213 
1214  const auto tit = parent_tris_refined[idx]; // parent triangle
1215  const auto &parent_edges_bit =
1216  parent_edges_bit_vec[idx]; // marks what edges are refined
1217  const auto nb_new_tris =
1218  nb_new_tris_vec[idx]; // number of triangles created by refinement
1219 
1220  // Iterate over refined triangles
1221  std::array<EntityHandle, 4> ref_tris;
1222  for (int tt = 0; tt != nb_new_tris; ++tt, ++start_e)
1223  ref_tris[tt] = start_e;
1224 
1225  if (nb_new_tris) {
1226 
1227  CHKERR moab.tag_clear_data(cOre.get_th_RefBitEdge(), &ref_tris[0],
1228  nb_new_tris, &parent_edges_bit);
1229  CHKERR moab.tag_clear_data(cOre.get_th_RefParentHandle(), &ref_tris[0],
1230  nb_new_tris, &tit);
1231 
1232  // Hash map of nodes (RefEntity) by edges (EntityHandle)
1233  std::map<EntityHandle /*edge*/, EntityHandle /*node*/>
1234  map_ref_nodes_by_edges;
1235 
1236  // find parents for new edges and faces
1237  // get tet edges and faces
1238  Range tit_edges; // parent edges
1239  CHKERR moab.get_adjacencies(&tit, 1, 1, false, tit_edges);
1240  for (auto edge : tit_edges) {
1241  auto eit = ref_parent_ents_view.find(edge);
1242  if (eit != ref_parent_ents_view.end()) {
1243  map_ref_nodes_by_edges[(*eit)->getParentEnt()] = eit->get()->getEnt();
1244  }
1245  }
1246 
1247  Range edges_nodes[3];
1248  // for edges - add ref nodes
1249  // edges_nodes[ee] - contains all nodes on edge ee including mid nodes if
1250  // exist
1251  Range::iterator eit = tit_edges.begin();
1252  for (int ee = 0; eit != tit_edges.end(); eit++, ee++) {
1253  CHKERR moab.get_connectivity(&*eit, 1, edges_nodes[ee], true);
1254  auto map_miit = map_ref_nodes_by_edges.find(*eit);
1255  if (map_miit != map_ref_nodes_by_edges.end()) {
1256  edges_nodes[ee].insert(map_miit->second);
1257  }
1258  }
1259 
1260  // add ref nodes to tri
1261  // tet_nodes contains all nodes on tet including mid edge nodes
1262  Range tri_nodes;
1263  CHKERR moab.get_connectivity(&tit, 1, tri_nodes, true);
1264  for (auto map_miit = map_ref_nodes_by_edges.begin();
1265  map_miit != map_ref_nodes_by_edges.end(); map_miit++) {
1266  tri_nodes.insert(map_miit->second);
1267  }
1268 
1269  Range ref_edges; // edges on all new tets
1270  // Get all all edges of refined tets
1271  CHKERR moab.get_adjacencies(ref_tris.data(), nb_new_tris, 1, false,
1272  ref_edges, moab::Interface::UNION);
1273 
1274  // Check for all ref edge and set parents
1275  for (auto reit = ref_edges.begin(); reit != ref_edges.end(); ++reit) {
1276  Range ref_edges_nodes;
1277  CHKERR moab.get_connectivity(&*reit, 1, ref_edges_nodes, true);
1278  // Check if ref edge is an coarse edge (loop over coarse tet edges)
1279  int ee = 0;
1280  for (; ee < 3; ee++) {
1281  // Two nodes are common (node[0],node[1],ref_node (if exist))
1282  // this tests if given edge is contained by edge of refined
1283  // triangle
1284  if (intersect(edges_nodes[ee], ref_edges_nodes).size() == 2) {
1285  EntityHandle edge = tit_edges[ee];
1286  CHKERR set_parent(*reit, edge, refined_ents_ptr, cOre);
1287  break;
1288  }
1289  }
1290  if (ee < 3)
1291  continue; // this refined edge is contained by edge of triangle
1292 
1293  // check if ref edge is in coarse tetrahedral (i.e. that is internal
1294  // edge of refined tetrahedral)
1295  if (intersect(tri_nodes, ref_edges_nodes).size() == 2) {
1296  CHKERR set_parent(*reit, tit, refined_ents_ptr, cOre);
1297  continue;
1298  }
1299 
1300  // Refined edge is not child of any edge, face or tetrahedral, this is
1301  // Impossible edge
1302  SETERRQ(m_field.get_comm(), MOFEM_DATA_INCONSISTENCY,
1303  "impossible refined edge");
1304  }
1305  }
1306  }
1307 
1308  CHKERR set_parent(refined_ents_ptr);
1309  CHKERR m_field.getInterface<BitRefManager>()->setBitRefLevel(
1310  ents_to_set_bit, bit, true, verb, &ents);
1311 
1312  MoFEMFunctionReturn(0);
1313 }
1314 
1315 } // 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
MoFEM::Ent_mi_tag
Definition: TagMultiIndices.hpp:21
MoFEM::CoreTmp< 0 >::addPrismToDatabase
MoFEMErrorCode addPrismToDatabase(const EntityHandle prism, int verb=DEFAULT_VERBOSITY)
add prim element
Definition: Core.cpp:522
MoFEM::MeshRefinement::SetParent::parentsToChange
map< EntityHandle, EntityHandle > parentsToChange
Definition: MeshRefinement.hpp:182
MoFEM::MeshRefinement::refineTris
MoFEMErrorCode refineTris(const EntityHandle meshset, const BitRefLevel &bit, int verb=QUIET, const bool debug=false)
refine triangles in the meshset
Definition: MeshRefinement.cpp:921
MoFEM::CoreTmp< 0 >
Core (interface) class.
Definition: Core.hpp:82
MoFEM::MeshRefinement::SetEdgeBitsFun
boost::function< MoFEMErrorCode(moab::Interface &moab, RefEntity_multiIndex_view_by_ordered_parent_entity &ref_parent_ents_view, EntityHandle tet, BitRefEdges &parent_edges_bit, EntityHandle *edge_new_nodes, int *split_edges)> SetEdgeBitsFun
Functions setting edges for refinemnt on enetity level.
Definition: MeshRefinement.hpp:203
FTensor::Tensor1
Definition: Tensor1_value.hpp:8
EntityHandle
MoFEM::MeshRefinement::addVerticesInTheMiddleOfEdges
MoFEMErrorCode addVerticesInTheMiddleOfEdges(const EntityHandle meshset, const BitRefLevel &bit, const bool recursive=false, int verb=QUIET, EntityHandle start_v=0)
make vertices in the middle of edges in meshset and add them to refinement levels defined by bit
Definition: MeshRefinement.cpp:42
NOISY
@ NOISY
Definition: definitions.h:224
MoFEM::CoreInterface::get_comm
virtual MPI_Comm & get_comm() const =0
MoFEM::Exceptions::MoFEMException
Exception to catch.
Definition: Exceptions.hpp:20
MoFEM::CoreTmp< 0 >::get_th_RefBitEdge
Tag get_th_RefBitEdge() const
Definition: Core.hpp:199
MoFEM::RefElement_multiIndex
multi_index_container< boost::shared_ptr< RefElement >, indexed_by< ordered_unique< tag< Ent_mi_tag >, const_mem_fun< RefElement::interface_type_RefEntity, EntityHandle, &RefElement::getEnt > > > > RefElement_multiIndex
Definition: RefElementMultiIndices.hpp:180
MoFEM::MeshRefinement::SetParent
Definition: MeshRefinement.hpp:181
MoFEM::Exceptions::MoFEMErrorCode
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
MoFEM::RefEntityTmp< 0 >
Struct keeps handle to refined handle.
Definition: RefEntsMultiIndices.hpp:141
VERY_NOISY
@ VERY_NOISY
Definition: definitions.h:225
MoFEM::RefElement::getBitRefEdgesUlong
virtual int getBitRefEdgesUlong() const
Definition: RefElementMultiIndices.hpp:28
MoFEM::interface_RefEntity< RefEntity >::getEnt
EntityHandle getEnt() const
Get the entity handle.
Definition: RefEntsMultiIndices.hpp:603
MoFEM::tet_type_2
int tet_type_2(const EntityHandle *conn, const int *split_edges, const EntityHandle *edge_new_nodes, EntityHandle *new_tets_conn)
Definition: EntityRefine.cpp:988
MoFEM::prism_type_2
MoFEMErrorCode prism_type_2(const EntityHandle *conn, const BitRefEdges split_edges, const EntityHandle *edge_new_nodes, EntityHandle *new_prism_conn)
Definition: EntityRefine.cpp:1266
MoFEM::no_handle
const EntityHandle no_handle
No entity handle is indicated by zero handle, i.e. root meshset.
Definition: Common.hpp:12
MoFEM::RefEntity_multiIndex_view_by_ordered_parent_entity
multi_index_container< boost::shared_ptr< RefEntity >, indexed_by< ordered_non_unique< const_mem_fun< RefEntity, EntityHandle, &RefEntity::getParentEnt > >, hashed_unique< tag< Composite_EntType_and_ParentEntType_mi_tag >, composite_key< boost::shared_ptr< RefEntity >, const_mem_fun< RefEntity, EntityHandle, &RefEntity::getEnt >, const_mem_fun< RefEntity, EntityHandle, &RefEntity::getParentEnt > >> > > RefEntity_multiIndex_view_by_ordered_parent_entity
Definition: RefEntsMultiIndices.hpp:792
MoFEM::RefEntity_multiIndex
multi_index_container< boost::shared_ptr< RefEntity >, indexed_by< ordered_unique< tag< Ent_mi_tag >, const_mem_fun< RefEntity, EntityHandle, &RefEntity::getEnt > >, ordered_non_unique< tag< Ent_Ent_mi_tag >, const_mem_fun< RefEntity, EntityHandle, &RefEntity::getParentEnt > >, ordered_non_unique< tag< Composite_EntType_and_ParentEntType_mi_tag >, composite_key< RefEntity, const_mem_fun< RefEntity, EntityType, &RefEntity::getEntType >, const_mem_fun< RefEntity, EntityType, &RefEntity::getParentEntType > > >, ordered_non_unique< tag< Composite_ParentEnt_And_EntType_mi_tag >, composite_key< RefEntity, const_mem_fun< RefEntity, EntityType, &RefEntity::getEntType >, const_mem_fun< RefEntity, EntityHandle, &RefEntity::getParentEnt > > > > > RefEntity_multiIndex
Definition: RefEntsMultiIndices.hpp:760
MoFEM::RefElement::interface_type_RefEntity
interface_RefEntity< RefEntity > interface_type_RefEntity
Definition: RefElementMultiIndices.hpp:18
MoFEM::MeshRefinement::cOre
MoFEM::Core & cOre
Definition: MeshRefinement.hpp:31
MOFEM_IMPOSSIBLE_CASE
@ MOFEM_IMPOSSIBLE_CASE
Definition: definitions.h:35
MoFEM::interface_RefEntity< RefEntity >::getParentEnt
EntityHandle getParentEnt() const
Definition: RefEntsMultiIndices.hpp:615
MoFEM::prism_type_1
MoFEMErrorCode prism_type_1(const EntityHandle *conn, const BitRefEdges split_edges, const EntityHandle *edge_new_nodes, EntityHandle *new_prism_conn)
Definition: EntityRefine.cpp:1230
MoFEM::tri_type_3
MoFEMErrorCode tri_type_3(const EntityHandle *conn, const EntityHandle *edge_new_nodes, EntityHandle *new_tris_conn)
Definition: EntityRefine.cpp:1111
MoFEM::DeprecatedCoreInterface
Deprecated interface functions.
Definition: DeprecatedCoreInterface.hpp:16
I
constexpr IntegrationType I
Definition: operators_tests.cpp:31
MoFEM::Interface
DeprecatedCoreInterface Interface
Definition: Interface.hpp:2010
MoFEM::MeshRefinement::refineTets
MoFEMErrorCode refineTets(const EntityHandle meshset, const BitRefLevel &bit, int verb=QUIET, const bool debug=false)
refine TET in the meshset
Definition: MeshRefinement.cpp:197
MoFEM::tet_type_3
int tet_type_3(const EntityHandle *conn, const int *split_edges, const EntityHandle *edge_new_nodes, EntityHandle *new_tets_conn)
Definition: EntityRefine.cpp:565
MoFEM::Composite_EntType_and_ParentEntType_mi_tag
Definition: TagMultiIndices.hpp:77
VERBOSE
@ VERBOSE
Definition: definitions.h:222
bool
CHKERR
#define CHKERR
Inline error check.
Definition: definitions.h:548
MoFEM::CoreInterface::get_moab
virtual moab::Interface & get_moab()=0
MoFEM::RefElement_multiIndex_parents_view
multi_index_container< boost::shared_ptr< RefElement >, indexed_by< ordered_unique< tag< Ent_mi_tag >, const_mem_fun< RefElement::interface_type_RefEntity, EntityHandle, &RefElement::getEnt > >, ordered_non_unique< tag< Ent_Ent_mi_tag >, const_mem_fun< RefElement::interface_type_RefEntity, EntityHandle, &RefElement::getParentEnt > >, ordered_non_unique< tag< Composite_ParentEnt_And_BitsOfRefinedEdges_mi_tag >, composite_key< RefElement, const_mem_fun< RefElement::interface_type_RefEntity, EntityHandle, &RefElement::getParentEnt >, const_mem_fun< RefElement, int, &RefElement::getBitRefEdgesUlong > > > > > RefElement_multiIndex_parents_view
Definition: MeshRefinement.cpp:30
MoFEM::tet_type_5
int tet_type_5(moab::Interface &moab, const EntityHandle *conn, const EntityHandle *edge_new_nodes, EntityHandle *new_tets_conn)
Definition: EntityRefine.cpp:133
MoFEM
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:10
MoFEM::MeshRefinement::MeshRefinement
MeshRefinement(const MoFEM::Core &core)
Definition: MeshRefinement.cpp:39
MoFEM::RefEntity_multiIndex_view_by_hashed_parent_entity
multi_index_container< boost::shared_ptr< RefEntity >, indexed_by< hashed_non_unique< const_mem_fun< RefEntity, EntityHandle, &RefEntity::getParentEnt > >, hashed_unique< tag< Composite_EntType_and_ParentEntType_mi_tag >, composite_key< boost::shared_ptr< RefEntity >, const_mem_fun< RefEntity, EntityHandle, &RefEntity::getEnt >, const_mem_fun< RefEntity, EntityHandle, &RefEntity::getParentEnt > >> > > RefEntity_multiIndex_view_by_hashed_parent_entity
multi-index view of RefEntity by parent entity
Definition: RefEntsMultiIndices.hpp:777
MoFEM::MeshRefinement::refineTrisHangingNodes
MoFEMErrorCode refineTrisHangingNodes(const EntityHandle meshset, const BitRefLevel &bit, int verb=QUIET, const bool debug=false)
refine TRI in the meshset
Definition: MeshRefinement.cpp:1002
MoFEM::MeshRefinement
Mesh refinement interface.
Definition: MeshRefinement.hpp:26
MoFEM::get_temp_meshset_ptr
auto get_temp_meshset_ptr(moab::Interface &moab)
Create smart pointer to temporary meshset.
Definition: Templates.hpp:1886
bit
auto bit
set bit
Definition: hanging_node_approx.cpp:75
MoFEM::CoreInterface::get_basic_entity_data_ptr
virtual boost::shared_ptr< BasicEntityData > & get_basic_entity_data_ptr()=0
Get pointer to basic entity data.
MoFEM::CoreTmp< 0 >::get_th_RefParentHandle
Tag get_th_RefParentHandle() const
Definition: Core.hpp:197
MoFEM::CoreTmp< 0 >::get_th_RefBitLevel
Tag get_th_RefBitLevel() const
Definition: Core.hpp:198
MoFEM::MeshRefinement::SetParent::operator()
MoFEMErrorCode operator()(const EntityHandle ent, const EntityHandle parent, const RefEntity_multiIndex *ref_ents_ptr, MoFEM::Core &cOre)
Definition: MeshRefinement.cpp:209
MoFEM::RefElement
keeps data about abstract refined finite element
Definition: RefElementMultiIndices.hpp:16
MoFEM::RefEntity
RefEntityTmp< 0 > RefEntity
Definition: RefEntsMultiIndices.hpp:566
debug
static const bool debug
Definition: dm_create_subdm.cpp:12
MoFEM::Exceptions::MoFEMException::errorMessage
char errorMessage[1024]
Definition: Exceptions.hpp:22
i
FTensor::Index< 'i', SPACE_DIM > i
Definition: hcurl_divergence_operator_2d.cpp:27
FTensor::Index< 'i', 3 >
MoFEM::tri_type_2
MoFEMErrorCode tri_type_2(const EntityHandle *conn, const int *split_edges, const EntityHandle *edge_new_nodes, EntityHandle *new_tris_conn)
Definition: EntityRefine.cpp:1177
MoFEM::MeshRefinement::refinePrisms
MoFEMErrorCode refinePrisms(const EntityHandle meshset, const BitRefLevel &bit, int verb=QUIET)
refine PRISM in the meshset
Definition: MeshRefinement.cpp:722
MoFEM::UnknownInterface
base class for all interface classes
Definition: UnknownInterface.hpp:34
Range
MoFEM::Exceptions::MoFEMException::errorCode
const int errorCode
Definition: Exceptions.hpp:21
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_TAG_AND_LOG
#define MOFEM_TAG_AND_LOG(channel, severity, tag)
Tag and log in channel.
Definition: LogManager.hpp:362
MoFEM::Composite_ParentEnt_And_BitsOfRefinedEdges_mi_tag
Definition: TagMultiIndices.hpp:79
MoFEM::RefElement_PRISM
keeps data about abstract PRISM finite element
Definition: RefElementMultiIndices.hpp:66
HenckyOps::f
auto f
Definition: HenckyOps.hpp:15
MoFEM::rangeInserter
auto rangeInserter(const I f, const I s, boost::function< bool(I it)> tester, boost::function< MoFEMErrorCode(I f, I s)> inserter)
Insert ranges.
Definition: Templates.hpp:1945
MoFEM::Types::BitRefEdges
std::bitset< BITREFEDGES_SIZE > BitRefEdges
Definition: Types.hpp:34
j
FTensor::Index< 'j', 3 > j
Definition: matrix_function.cpp:19
MoFEM::RefEntity_change_parent
change parent
Definition: RefEntsMultiIndices.hpp:812
MoFEM::CoreInterface::get_ref_ents
virtual const RefEntity_multiIndex * get_ref_ents() const =0
Get the ref ents object.
MoFEM::tet_type_4
int tet_type_4(const EntityHandle *conn, const int *split_edges, const EntityHandle *edge_new_nodes, EntityHandle *new_tets_conn)
Definition: EntityRefine.cpp:239
MoFEM::MeshRefinement::query_interface
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
Definition: MeshRefinement.cpp:33
MOFEM_DATA_INCONSISTENCY
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:31
MoFEM::CoreInterface::get_ref_finite_elements
virtual const RefElement_multiIndex * get_ref_finite_elements() const =0
Get the ref finite elements object.
sdf_hertz_2d_axisymm_plane.d
float d
Definition: sdf_hertz_2d_axisymm_plane.py:4
MoFEM::prism_type_3
MoFEMErrorCode prism_type_3(const EntityHandle *conn, const BitRefEdges split_edges, const EntityHandle *edge_new_nodes, EntityHandle *new_prism_conn)
Definition: EntityRefine.cpp:1339
MoFEM::tet_type_1
void tet_type_1(const EntityHandle *conn, const int split_edge, const EntityHandle edge_new_node, EntityHandle *new_tets_conn)
Definition: EntityRefine.cpp:1096
EntityRefine.hpp
MoFEM::tet_type_6
void tet_type_6(moab::Interface &moab, const EntityHandle *conn, const EntityHandle *edge_new_nodes, EntityHandle *new_tets_conn)
Definition: EntityRefine.cpp:31
MoFEM::BitRefManager
Managing BitRefLevels.
Definition: BitRefManager.hpp:21
MoFEM::MeshRefinement::refineTetsHangingNodes
MoFEMErrorCode refineTetsHangingNodes(const Range &tets, const BitRefLevel &bit, int verb=QUIET, const bool debug=false)
refine TET in the meshset
Definition: MeshRefinement.cpp:276
MoFEM::Types::BitRefLevel
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Types.hpp:40
MoFEMFunctionBeginHot
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:453
MoFEM::tri_type_1
MoFEMErrorCode tri_type_1(const EntityHandle *conn, const int split_edge, const EntityHandle edge_new_node, EntityHandle *new_tris_conn)
Definition: EntityRefine.cpp:1137
convert.int
int
Definition: convert.py:64
MoFEMFunctionReturn
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:429
MoFEMFunctionBegin
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359
MoFEM::MeshRefinement::refineMeshset
MoFEMErrorCode refineMeshset(const EntityHandle meshset, const BitRefLevel &bit, const bool recursive=false, int verb=QUIET)
refine meshset, i.e. add child of refined entities to meshset
Definition: MeshRefinement.cpp:904
Generated by Doxygen 1.8.17 and hosted at University of Glasgow