v0.9.1
Classes | Public Types | Public Member Functions | Private Types | Private Member Functions | Private Attributes | List of all members
MoFEM::NodeMergerInterface Struct Reference

Merge node by collapsing edge between them. More...

#include <src/interfaces/NodeMerger.hpp>

Inheritance diagram for MoFEM::NodeMergerInterface:
[legend]
Collaboration diagram for MoFEM::NodeMergerInterface:
[legend]

Classes

struct  FaceMap
 
struct  ParentChild
 

Public Types

typedef multi_index_container< ParentChild, indexed_by< ordered_unique< member< ParentChild, EntityHandle, &ParentChild::pArent > >, ordered_non_unique< member< ParentChild, EntityHandle, &ParentChild::cHild > > > > ParentChildMap
 

Public Member Functions

MoFEMErrorCode query_interface (const MOFEMuuid &uuid, UnknownInterface **iface) const
 
 NodeMergerInterface (const MoFEM::Core &core)
 
MoFEMErrorCode getSubInterfaceOptions ()
 
bool getSuccessMerge ()
 Return true if successful merge. More...
 
void setErrorIfNoCommonEdge (const bool b=true)
 Set error if no common edge. More...
 
MoFEMErrorCode mergeNodes (EntityHandle father, EntityHandle mother, Range &out_tets, Range *tets_ptr=NULL, const bool only_if_improve_quality=false, const double move=0, const int line_search=0, Tag th=NULL, const int verb=0)
 merge nodes which sharing edge More...
 
MoFEMErrorCode mergeNodes (EntityHandle father, EntityHandle mother, BitRefLevel bit, Range *tets_ptr=NULL, const bool only_if_improve_quality=false, const double move=0, Tag th=NULL)
 merge nodes which sharing edge More...
 
MoFEMErrorCode mergeNodes (EntityHandle father, EntityHandle mother, BitRefLevel bit, BitRefLevel tets_from_bit_ref_level, const bool only_if_improve_quality=false, const double move=0, Tag th=NULL)
 merge nodes which sharing edge More...
 
ParentChildMapgetParentChildMap ()
 Get map of parent cand child. More...
 
- Public Member Functions inherited from MoFEM::UnknownInterface
template<class IFACE >
MoFEMErrorCode registerInterface (const MOFEMuuid &uuid, bool error_if_registration_failed=true)
 Register interface. More...
 
template<class IFACE , bool VERIFY = false>
MoFEMErrorCode getInterface (const MOFEMuuid &uuid, IFACE *&iface) const
 Get interface by uuid and return reference to pointer of interface. More...
 
template<class IFACE >
MoFEMErrorCode getInterface (IFACE *&iface) const
 Get interface refernce to pointer of interface. More...
 
template<class IFACE >
MoFEMErrorCode getInterface (IFACE **const iface) const
 Get interface pointer to pointer of interface. More...
 
template<class IFACE , typename boost::enable_if< boost::is_pointer< IFACE >, int >::type = 0>
IFACE getInterface () const
 Get interface pointer to pointer of interface. More...
 
template<class IFACE , typename boost::enable_if< boost::is_reference< IFACE >, int >::type = 0>
IFACE getInterface () const
 Get reference to interface. More...
 
template<class IFACE >
IFACE * getInterface () const
 Function returning pointer to interface. More...
 
virtual ~UnknownInterface ()=default
 
virtual MoFEMErrorCode getLibVersion (Version &version) const
 Get library version. More...
 
virtual const MoFEMErrorCode getFileVersion (moab::Interface &moab, Version &version) const
 Get database major version. More...
 
virtual MoFEMErrorCode getInterfaceVersion (Version &version) const
 Get database major version. More...
 
template<>
MoFEMErrorCode getInterface (const MOFEMuuid &uuid, UnknownInterface *&iface) const
 

Private Types

typedef multi_index_container< FaceMap, indexed_by< hashed_unique< composite_key< FaceMap, member< FaceMap, EntityHandle, &FaceMap::n0 >, member< FaceMap, EntityHandle, &FaceMap::n1 > > > > > FaceMapIdx
 

Private Member Functions

MoFEMErrorCode minQuality (Range &check_tests, EntityHandle father, EntityHandle mother, double *coords_move, double &min_quality, Tag th=NULL, boost::function< double(double, double)> f=[](double a, double b) -> double { return std::min(a, b);})
 Calualte quality if nodes merged. More...
 
MoFEMErrorCode lineSearch (Range &check_tests, EntityHandle father, EntityHandle mother, int line_search, FTensor::Tensor1< double, 3 > &t_move, Tag th=NULL)
 Use bisection method to find point of edge collapse. More...
 

Private Attributes

MoFEM::CorecOre
 
boost::function< double(const double a, const double b)> minQualityFunction
 
bool successMerge
 True if marge is success. More...
 
bool errorIfNoCommonEdge
 Send error if no common edge. More...
 
ParentChildMap parentChildMap
 

Additional Inherited Members

- Protected Member Functions inherited from MoFEM::UnknownInterface
boost::typeindex::type_index getClassIdx (const MOFEMuuid &uid) const
 Get type name for interface Id. More...
 
MOFEMuuid getUId (const boost::typeindex::type_index &class_idx) const
 Get interface Id for class name. More...
 

Detailed Description

Merge node by collapsing edge between them.

Definition at line 30 of file NodeMerger.hpp.

Member Typedef Documentation

◆ FaceMapIdx

typedef multi_index_container< FaceMap, indexed_by<hashed_unique<composite_key< FaceMap, member<FaceMap, EntityHandle, &FaceMap::n0>, member<FaceMap, EntityHandle, &FaceMap::n1> > > > > MoFEM::NodeMergerInterface::FaceMapIdx
private

Definition at line 183 of file NodeMerger.hpp.

◆ ParentChildMap

typedef multi_index_container< ParentChild, indexed_by<ordered_unique< member<ParentChild, EntityHandle, &ParentChild::pArent> >, ordered_non_unique< member<ParentChild, EntityHandle, &ParentChild::cHild> > > > MoFEM::NodeMergerInterface::ParentChildMap

Definition at line 128 of file NodeMerger.hpp.

Constructor & Destructor Documentation

◆ NodeMergerInterface()

MoFEM::NodeMergerInterface::NodeMergerInterface ( const MoFEM::Core core)

Definition at line 37 of file NodeMerger.cpp.

38  : cOre(const_cast<MoFEM::Core &>(core)), minQualityFunction(min_non_abs),
39  successMerge(false), errorIfNoCommonEdge(false) {}
boost::function< double(const double a, const double b)> minQualityFunction
Definition: NodeMerger.hpp:138
static auto min_non_abs(const double a, const double b)
Definition: NodeMerger.cpp:33
bool successMerge
True if marge is success.
Definition: NodeMerger.hpp:140
bool errorIfNoCommonEdge
Send error if no common edge.
Definition: NodeMerger.hpp:141

Member Function Documentation

◆ getParentChildMap()

ParentChildMap& MoFEM::NodeMergerInterface::getParentChildMap ( )

Get map of parent cand child.

Returns

Definition at line 134 of file NodeMerger.hpp.

134 { return parentChildMap; }
ParentChildMap parentChildMap
Definition: NodeMerger.hpp:171

◆ getSubInterfaceOptions()

MoFEMErrorCode MoFEM::NodeMergerInterface::getSubInterfaceOptions ( )

Definition at line 41 of file NodeMerger.cpp.

41  {
42  Interface &m_field = cOre;
44  ierr = PetscOptionsBegin(m_field.get_comm(), "node_merge",
45  "Node merge options", "none");
46  CHKERRQ(ierr);
47 
48  ierr = PetscOptionsEnd();
49  CHKERRQ(ierr);
50 
52 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
CHKERRQ(ierr)
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1791
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412

◆ getSuccessMerge()

bool MoFEM::NodeMergerInterface::getSuccessMerge ( )

Return true if successful merge.

Returns
Error code

Definition at line 44 of file NodeMerger.hpp.

44 { return successMerge; }
bool successMerge
True if marge is success.
Definition: NodeMerger.hpp:140

◆ lineSearch()

MoFEMErrorCode MoFEM::NodeMergerInterface::lineSearch ( Range &  check_tests,
EntityHandle  father,
EntityHandle  mother,
int  line_search,
FTensor::Tensor1< double, 3 > &  t_move,
Tag  th = NULL 
)
private

Use bisection method to find point of edge collapse.

Parameters
check_testsrange of tets to check quality
fatherfirst node of the edge
mothersecond node of the edge
line_searchnumber of iterations
coords_movenode to move
Returns
error code

Definition at line 420 of file NodeMerger.cpp.

422  {
423  Interface &m_field = cOre;
425 
426  EntityHandle conn[] = {father, mother};
427 
428  double coords[6];
429  if (th == NULL) {
430  CHKERR m_field.get_moab().get_coords(conn, 2, coords);
431  } else {
432  CHKERR m_field.get_moab().tag_get_data(th, conn, 2, coords);
433  }
434 
436  FTensor::Tensor1<double, 3> t_coords(coords[0], coords[1], coords[2]);
438  for (int nn = 0; nn != 3; nn++) {
439  t_delta(nn) = coords[3 + nn] - t_coords(nn);
440  }
441 
442  t_move(i) = t_coords(i);
443  double min_quality_i = 1;
444  CHKERR minQuality(check_tests, father, mother, &t_move(0), min_quality_i, th,
446 
447  t_move(i) = t_coords(i) + t_delta(i);
448  double min_quality_k = 1;
449  CHKERR minQuality(check_tests, father, mother, &t_move(0), min_quality_k, th,
451 
452  double alpha_i = 0;
453  double alpha_k = 1;
454 
455  for (int ii = 0; ii != line_search; ii++) {
456  double min_quality = 1;
457  double alpha = (alpha_i + alpha_k) * 0.5;
458  t_move(i) = t_coords(i) + alpha * t_delta(i);
459  CHKERR minQuality(check_tests, father, mother, &t_move(0), min_quality, th,
461  if (min_quality_i >= min_quality_k) {
462  min_quality_k = min_quality;
463  alpha_k = alpha;
464  } else {
465  min_quality_i = min_quality;
466  alpha_i = alpha;
467  }
468  }
469 
470  if (min_quality_i > min_quality_k) {
471  t_move(i) = t_coords(i) + alpha_i * t_delta(i);
472  } else {
473  t_move(i) = t_coords(i) + alpha_k * t_delta(i);
474  }
475 
477 }
boost::function< double(const double a, const double b)> minQualityFunction
Definition: NodeMerger.hpp:138
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
#define CHKERR
Inline error check.
Definition: definitions.h:601
FTensor::Index< 'i', 2 > i
[Common data]
Definition: ContactOps.hpp:24
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1791
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412
MoFEMErrorCode minQuality(Range &check_tests, EntityHandle father, EntityHandle mother, double *coords_move, double &min_quality, Tag th=NULL, boost::function< double(double, double)> f=[](double a, double b) -> double { return std::min(a, b);})
Calualte quality if nodes merged.
Definition: NodeMerger.cpp:356

◆ mergeNodes() [1/3]

MoFEMErrorCode MoFEM::NodeMergerInterface::mergeNodes ( EntityHandle  father,
EntityHandle  mother,
Range &  out_tets,
Range *  tets_ptr = NULL,
const bool  only_if_improve_quality = false,
const double  move = 0,
const int  line_search = 0,
Tag  th = NULL,
const int  verb = 0 
)

merge nodes which sharing edge

Father is sties, mother is merged.

Parameters
fathernode to which mother is merged to.
mothermerged node
out_tetstetrahedra after merge
tets_ptrtest only tets_ptr from range are changed
only_if_improve_qualityDo merge if that improve quality
movefather by fraction of edge length move=[0,1]

Move node on the edge, 0 not move, 1 move to mother side, 0.5 will be in the middle.

Definition at line 54 of file NodeMerger.cpp.

57  {
58  Interface &m_field = cOre;
61 
62  // Get edges adjacent to father and mother, i.e. mother is merged to father.
63  Range father_edges;
64  CHKERR m_field.get_moab().get_adjacencies(&father, 1, 1, false, father_edges);
65  Range mother_edges;
66  CHKERR m_field.get_moab().get_adjacencies(&mother, 1, 1, false, mother_edges);
67 
68  // Get tets adjacent to mother and father
69  Range father_tets;
70  CHKERR m_field.get_moab().get_adjacencies(&father, 1, 3, false, father_tets);
71  Range mother_tets;
72  CHKERR m_field.get_moab().get_adjacencies(&mother, 1, 3, false, mother_tets);
73  if (tets_ptr != NULL) {
74  mother_tets = intersect(mother_tets, *tets_ptr);
75  father_tets = intersect(father_tets, *tets_ptr);
76  }
77 
78  // Find common edge
79  Range common_edge;
80  common_edge = intersect(father_edges, mother_edges);
81  if (tets_ptr != NULL) {
82  Range tets = unite(father_tets, mother_tets);
83  Range tets_edges;
84  CHKERR m_field.get_moab().get_adjacencies(tets, 1, false, tets_edges,
85  moab::Interface::UNION);
86  common_edge = intersect(common_edge, tets_edges);
87  father_edges = intersect(father_edges, tets_edges);
88  mother_edges = intersect(mother_edges, tets_edges);
89  }
90 
91  // No common edge, merge no possible
92  if (errorIfNoCommonEdge && common_edge.empty()) {
93  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
94  "no common edge between nodes");
95  } else if (common_edge.empty()) {
96  Range seed_tets;
97  if (tets_ptr != NULL) {
98  seed_tets.merge(*tets_ptr);
99  }
100  out_tets = seed_tets;
101  successMerge = false;
103  }
104 
105  // Common edge tets, that tests will be squashed
106  Range edge_tets;
107  CHKERR m_field.get_moab().get_adjacencies(common_edge, 3, true, edge_tets);
108  // Intersect with ptr_tets (usually associated with some bit level)
109  if (tets_ptr != NULL) {
110  edge_tets = intersect(edge_tets, *tets_ptr);
111  }
112  // Mother tets, has only one mother vertex and no father vertex.
113  mother_tets = subtract(mother_tets, edge_tets);
114  father_tets = subtract(father_tets, edge_tets);
115 
116  auto get_coords = [&m_field](Tag th, const EntityHandle *conn,
117  const int num_nodes) {
118  VectorDouble coords(3 * num_nodes);
119  if (th == NULL) {
120  CHKERR m_field.get_moab().get_coords(conn, num_nodes, &coords[0]);
121  } else {
122  CHKERR m_field.get_moab().tag_get_data(th, conn, num_nodes, &coords[0]);
123  }
124  return coords;
125  };
126 
127  auto get_tensor = [](VectorDouble &coords, const int shift) {
129  &coords[shift], &coords[shift + 1], &coords[shift + 2]);
130  };
131 
132  // move father coord is move > 0
134  if (move > 0) {
135  EntityHandle conn[] = {father, mother};
136  VectorDouble coords = get_coords(th, conn, 2);
137  auto t_n0 = get_tensor(coords, 0);
138  auto t_n1 = get_tensor(coords, 3);
139  t_move(i) = t_n0(i) + move * (t_n1(i) - t_n0(i));
140  }
141 
142  if (line_search > 0) {
143  Range check_tests = unite(father_tets, mother_tets);
144  CHKERR lineSearch(check_tests, father, mother, line_search, t_move, th);
145  }
146 
147  if (only_if_improve_quality) {
148  Range check_tests = mother_tets;
149  // no need to check father tets since no change of quality for them
150  if (move > 0 || line_search) {
151  check_tests.merge(father_tets);
152  }
153 
154  double min_quality0 = 1;
155  CHKERR minQuality(edge_tets, 0, 0, NULL, min_quality0, th,
157  CHKERR minQuality(check_tests, 0, 0, NULL, min_quality0, th,
159  double min_quality = 1;
160  CHKERR minQuality(check_tests, father, mother,
161  ((move > 0) || line_search) ? &t_move(0) : NULL,
162  min_quality, th, minQualityFunction);
163  if (min_quality < min_quality0) {
164  if (tets_ptr != NULL) {
165  out_tets = *tets_ptr;
166  }
167  successMerge = false;
169  }
170  }
171 
172  // Move node
173  if (move > 0 || line_search) {
174  if (th == NULL) {
175  CHKERR m_field.get_moab().set_coords(&father, 1, &t_move(0));
176  } else {
177  CHKERR m_field.get_moab().tag_set_data(th, &father, 1, &t_move(0));
178  }
179  }
180 
181  auto get_conn = [&m_field](const EntityHandle ent,
182  int *ret_num_nodes = nullptr) {
183  int num_nodes;
184  const EntityHandle *conn;
185  CHKERR m_field.get_moab().get_connectivity(ent, conn, num_nodes, true);
186  if (ret_num_nodes)
187  *ret_num_nodes = num_nodes;
188  return conn;
189  };
190 
191  auto create_tet = [this, &m_field](const EntityHandle *new_conn,
192  const EntityHandle parent) {
193  EntityHandle tet;
194  Range tets;
195  CHKERR m_field.get_moab().get_adjacencies(new_conn, 4, 3, false, tets);
196  bool tet_found = false;
197  for (auto it_tet : tets) {
198  const EntityHandle *tet_conn;
199  int num_nodes;
200  CHKERR m_field.get_moab().get_connectivity(it_tet, tet_conn, num_nodes,
201  true);
202  const EntityHandle *p = std::find(tet_conn, &tet_conn[4], new_conn[0]);
203  if (p != tet_conn + 4) {
204  int s = std::distance(tet_conn, p);
205  int n = 0;
206  for (; n != 4; ++n) {
207  const int idx[] = {0, 1, 2, 3, 0, 1, 2, 3};
208  if (tet_conn[idx[s + n]] != new_conn[n])
209  break;
210  }
211  if (n == 4 && !tet_found) {
212  tet = it_tet;
213  tet_found = true;
214  } else if (n == 4) {
215  THROW_MESSAGE("More that one tet with the same connectivity");
216  }
217  }
218  }
219  if (!tet_found) {
220  // Create tet with new connectivity
221  CHKERR m_field.get_moab().create_element(MBTET, new_conn, 4, tet);
222  CHKERR m_field.get_moab().tag_set_data(cOre.get_th_RefParentHandle(),
223  &tet, 1, &parent);
224  parentChildMap.insert(ParentChild(parent, tet));
225  }
226  return tet;
227  };
228 
229  // clear map
230  parentChildMap.clear();
231 
232  Range created_tets;
233  for (auto m_tet : mother_tets) {
234  const EntityHandle *conn = get_conn(m_tet);
235  EntityHandle new_conn[4];
236  // Replace mother vertices by father vertices
237  int nb_mother_verts = 0;
238  for (int nn = 0; nn != 4; ++nn) {
239  if (conn[nn] == father) {
240  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
241  "Tet has father vertex, impossible but here it is");
242  }
243  if (conn[nn] == mother) {
244  new_conn[nn] = father;
245  nb_mother_verts++;
246  } else {
247  new_conn[nn] = conn[nn];
248  }
249  }
250  if (nb_mother_verts != 1) {
251  SETERRQ1(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
252  "Tet should have only one vertex but have %d", nb_mother_verts);
253  }
254 
255  VectorDouble new_coords = get_coords(th, new_conn, 4);
256 
257  // add tet to range
258  created_tets.insert(create_tet(new_conn, m_tet));
259  }
260 
261  // Loop over mother adjacent entities to use them as parents
262  auto get_adj_ents = [&](const Range &ents, const bool create) {
263  Range adj;
264  for (int dd = 1; dd <= 2; dd++)
265  CHKERR m_field.get_moab().get_adjacencies(ents, dd, create, adj,
266  moab::Interface::UNION);
267  return adj;
268  };
269  auto adj_crated_ents = get_adj_ents(created_tets, true);
270  adj_crated_ents.erase(common_edge[0]);
271 
272  FaceMapIdx face_map;
273  for (auto ent : adj_crated_ents) {
274  int num_nodes;
275  const EntityHandle *conn = get_conn(ent, &num_nodes);
276  EntityHandle small_conn[num_nodes];
277  int ii = 0;
278  bool father_node = false;
279  for (int nn = 0; nn != num_nodes; nn++) {
280  if (conn[nn] == father)
281  father_node = true;
282  else
283  small_conn[ii++] = conn[nn];
284  }
285  if (father_node) {
286  if (ii > 1)
287  std::sort(&small_conn[0], &small_conn[ii]);
288  if (ii == 2)
289  face_map.insert(FaceMap(ent, small_conn[0], small_conn[1]));
290  else
291  face_map.insert(FaceMap(ent, small_conn[0], 0));
292  }
293  }
294 
295  auto adj_mother_ents = get_adj_ents(mother_tets, false);
296  adj_mother_ents.erase(common_edge[0]);
297  for (auto ent : adj_mother_ents) {
298  int num_nodes;
299  const EntityHandle *conn = get_conn(ent, &num_nodes);
300  EntityHandle small_conn[num_nodes];
301  int nb_new_node = 0;
302  int nn = 0;
303  int ii = 0;
304  for (; nn != num_nodes; ++nn) {
305  if (conn[nn] == mother) {
306  nb_new_node++;
307  } else {
308  small_conn[ii++] = conn[nn];
309  }
310  }
311  if (nb_new_node > 0) {
312  if (ii > 1)
313  std::sort(&small_conn[0], &small_conn[ii]);
314 
315  EntityHandle n0 = small_conn[0], n1 = 0;
316  if (ii == 2)
317  n1 = small_conn[1];
318 
319  FaceMapIdx::iterator fit = face_map.find(boost::make_tuple(n0, n1));
320  if (fit != face_map.end()) {
321  const EntityHandle child = fit->e;
322  const EntityHandle parent = ent;
323  if (m_field.get_moab().dimension_from_handle(parent) !=
324  m_field.get_moab().dimension_from_handle(child))
325  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
326  "Huston we have a problem!");
327 
328  // Set parent child relation
329  CHKERR m_field.get_moab().tag_set_data(cOre.get_th_RefParentHandle(),
330  &child, 1, &parent);
331  // Create map
332  parentChildMap.insert(ParentChild(parent, child));
333  }
334  }
335  }
336 
337  // Seed tets to given bit level
338  Range seed_tets;
339  if (tets_ptr != NULL) {
340  seed_tets.merge(*tets_ptr);
341  seed_tets = subtract(seed_tets, edge_tets);
342  seed_tets = subtract(seed_tets, mother_tets);
343  }
344  seed_tets.merge(created_tets);
345  out_tets.swap(seed_tets);
346 
347  successMerge = true;
348 
349  if (verb > VERY_VERBOSE)
350  std::cerr << "Nodes merged" << endl;
351 
353 }
boost::function< double(const double a, const double b)> minQualityFunction
Definition: NodeMerger.hpp:138
Tag get_th_RefParentHandle() const
Definition: Core.hpp:150
FTensor::Index< 'n', 2 > n
Definition: PlasticOps.hpp:33
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
#define THROW_MESSAGE(a)
Throw MoFEM exception.
Definition: definitions.h:625
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:513
bool successMerge
True if marge is success.
Definition: NodeMerger.hpp:140
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
ParentChildMap parentChildMap
Definition: NodeMerger.hpp:171
bool errorIfNoCommonEdge
Send error if no common edge.
Definition: NodeMerger.hpp:141
#define CHKERR
Inline error check.
Definition: definitions.h:601
MoFEMErrorCode lineSearch(Range &check_tests, EntityHandle father, EntityHandle mother, int line_search, FTensor::Tensor1< double, 3 > &t_move, Tag th=NULL)
Use bisection method to find point of edge collapse.
Definition: NodeMerger.cpp:420
FTensor::Index< 'i', 2 > i
[Common data]
Definition: ContactOps.hpp:24
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1791
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412
ublas::vector< double, DoubleAllocator > VectorDouble
Definition: Types.hpp:72
multi_index_container< FaceMap, indexed_by< hashed_unique< composite_key< FaceMap, member< FaceMap, EntityHandle, &FaceMap::n0 >, member< FaceMap, EntityHandle, &FaceMap::n1 > > > > > FaceMapIdx
Definition: NodeMerger.hpp:183
MoFEMErrorCode minQuality(Range &check_tests, EntityHandle father, EntityHandle mother, double *coords_move, double &min_quality, Tag th=NULL, boost::function< double(double, double)> f=[](double a, double b) -> double { return std::min(a, b);})
Calualte quality if nodes merged.
Definition: NodeMerger.cpp:356

◆ mergeNodes() [2/3]

MoFEMErrorCode MoFEM::NodeMergerInterface::mergeNodes ( EntityHandle  father,
EntityHandle  mother,
BitRefLevel  bit,
Range *  tets_ptr = NULL,
const bool  only_if_improve_quality = false,
const double  move = 0,
Tag  th = NULL 
)

merge nodes which sharing edge

Father is sties, mother is merged.

Parameters
fathernode to which mother is merged to.
mothermerged node
bitlevel of mesh merged nodes mesh
testonly tets_ptr from range are changed
only_if_improve_qualityDo merge if that improve quality
movefather by fraction of edge length move=[0,1]

Move node on the edge, 0 not move, 1 move to mother side, 0.5 will be in the middle.

Definition at line 479 of file NodeMerger.cpp.

481  {
482  Interface &m_field = cOre;
484  Range out_tets;
485  CHKERR mergeNodes(father, mother, out_tets, tets_ptr, only_if_improve_quality,
486  move, 0, th);
487  CHKERR m_field.getInterface<BitRefManager>()->setBitRefLevel(out_tets, bit);
489 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
MoFEMErrorCode mergeNodes(EntityHandle father, EntityHandle mother, Range &out_tets, Range *tets_ptr=NULL, const bool only_if_improve_quality=false, const double move=0, const int line_search=0, Tag th=NULL, const int verb=0)
merge nodes which sharing edge
Definition: NodeMerger.cpp:54
#define CHKERR
Inline error check.
Definition: definitions.h:601
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1791
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412

◆ mergeNodes() [3/3]

MoFEMErrorCode MoFEM::NodeMergerInterface::mergeNodes ( EntityHandle  father,
EntityHandle  mother,
BitRefLevel  bit,
BitRefLevel  tets_from_bit_ref_level,
const bool  only_if_improve_quality = false,
const double  move = 0,
Tag  th = NULL 
)

merge nodes which sharing edge

Father is sties, mother is merged.

Parameters
fathernode to which mother is merged to.
mothermerged node
tets_from_bit_ref_levelonly tetrahedrons from bit level are changed
only_if_improve_qualityDo merge if that improve quality
movefather by fraction of edge length move=[0,1]

Move node on the edge, 0 not move, 1 move to mother side, 0.5 will be in the middle.

Definition at line 491 of file NodeMerger.cpp.

494  {
495  Interface &m_field = cOre;
497  Range level_tets;
498  CHKERR m_field.getInterface<BitRefManager>()->getEntitiesByTypeAndRefLevel(
499  tets_from_bit_ref_level, BitRefLevel().set(), MBTET, level_tets);
500  CHKERR mergeNodes(father, mother, bit, &level_tets, only_if_improve_quality,
501  move, th);
503 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
MoFEMErrorCode mergeNodes(EntityHandle father, EntityHandle mother, Range &out_tets, Range *tets_ptr=NULL, const bool only_if_improve_quality=false, const double move=0, const int line_search=0, Tag th=NULL, const int verb=0)
merge nodes which sharing edge
Definition: NodeMerger.cpp:54
#define CHKERR
Inline error check.
Definition: definitions.h:601
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Types.hpp:51
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1791
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412

◆ minQuality()

MoFEMErrorCode MoFEM::NodeMergerInterface::minQuality ( Range &  check_tests,
EntityHandle  father,
EntityHandle  mother,
double coords_move,
double min_quality,
Tag  th = NULL,
boost::function< double(double, double)>  f = [](double a, double b) -> double { return std::min(a, b); } 
)
private

Calualte quality if nodes merged.

Parameters
check_teststets to check
fatherfirst node of the edge
mothersecond node of the edge
coords_movemoved father node
min_qualitycalculated quality
Returns
error code

Definition at line 356 of file NodeMerger.cpp.

359  {
360  Interface &m_field = cOre;
361  double coords[12];
364  for (Range::iterator tit = check_tests.begin(); tit != check_tests.end();
365  tit++) {
366  const EntityHandle *conn;
367  int num_nodes;
368  CHKERR m_field.get_moab().get_connectivity(*tit, conn, num_nodes, true);
369  if (mother > 0) {
370  EntityHandle new_conn[4];
371  // Replace mother vertices by father vertices
372  int nb_mother_verts = 0;
373  int father_nn = 0;
374  for (int nn = 0; nn < 4; nn++) {
375  if (conn[nn] == father) {
376  father_nn = nn;
377  }
378  if (conn[nn] == mother) {
379  new_conn[nn] = father;
380  father_nn = nn;
381  nb_mother_verts++;
382  } else {
383  new_conn[nn] = conn[nn];
384  }
385  }
386  if (nb_mother_verts > 1) {
387  SETERRQ1(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
388  "Tet should have no more than one mother vertex but have %d",
389  nb_mother_verts);
390  }
391  if (th == NULL) {
392  CHKERR m_field.get_moab().get_coords(new_conn, num_nodes, coords);
393  } else {
394  CHKERR m_field.get_moab().tag_get_data(th, new_conn, num_nodes, coords);
395  }
396  if (coords_move) {
397  int shift = 3 * father_nn;
398  for (int nn = 0; nn != 3; nn++) {
399  coords[shift + nn] = coords_move[nn];
400  }
401  }
402  } else {
403  if (th == NULL) {
404  CHKERR m_field.get_moab().get_coords(conn, num_nodes, coords);
405  } else {
406  CHKERR m_field.get_moab().tag_get_data(th, conn, num_nodes, coords);
407  }
408  }
409  double quality = Tools::volumeLengthQuality(coords);
410  if (!std::isnormal(quality)) {
411  min_quality = 0;
413  }
414  min_quality = f(min_quality, quality);
415  }
417 };
static double volumeLengthQuality(const double *coords)
Calculate tetrahedron volume length quality.
Definition: Tools.cpp:34
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:513
#define CHKERR
Inline error check.
Definition: definitions.h:601
FTensor::Index< 'i', 2 > i
[Common data]
Definition: ContactOps.hpp:24
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1791
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412

◆ query_interface()

MoFEMErrorCode MoFEM::NodeMergerInterface::query_interface ( const MOFEMuuid uuid,
UnknownInterface **  iface 
) const
virtual

Implements MoFEM::UnknownInterface.

Definition at line 21 of file NodeMerger.cpp.

22  {
24  *iface = NULL;
25  if (uuid == IDD_MOFEMNodeMerger) {
26  *iface = const_cast<NodeMergerInterface *>(this);
28  }
29  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "unknown interface");
31 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:506
static const MOFEMuuid IDD_MOFEMNodeMerger
Definition: NodeMerger.hpp:23
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:513

◆ setErrorIfNoCommonEdge()

void MoFEM::NodeMergerInterface::setErrorIfNoCommonEdge ( const bool  b = true)

Set error if no common edge.

Parameters
bIf true send error if false no error

Definition at line 50 of file NodeMerger.hpp.

50  {
52  }
bool errorIfNoCommonEdge
Send error if no common edge.
Definition: NodeMerger.hpp:141

Member Data Documentation

◆ cOre

MoFEM::Core& MoFEM::NodeMergerInterface::cOre
private

Definition at line 137 of file NodeMerger.hpp.

◆ errorIfNoCommonEdge

bool MoFEM::NodeMergerInterface::errorIfNoCommonEdge
private

Send error if no common edge.

Definition at line 141 of file NodeMerger.hpp.

◆ minQualityFunction

boost::function<double(const double a, const double b)> MoFEM::NodeMergerInterface::minQualityFunction
private

Definition at line 138 of file NodeMerger.hpp.

◆ parentChildMap

ParentChildMap MoFEM::NodeMergerInterface::parentChildMap
private

Definition at line 171 of file NodeMerger.hpp.

◆ successMerge

bool MoFEM::NodeMergerInterface::successMerge
private

True if marge is success.

Definition at line 140 of file NodeMerger.hpp.


The documentation for this struct was generated from the following files: