v0.9.0
CutMeshInterface.hpp
Go to the documentation of this file.
1 /** \file CutMeshInterface.hpp
2  * \brief Cut mesh interface
3  *
4  * \ingroup mesh_cut
5  */
6 
7 /*
8  * MoFEM is distributed in the hope that it will be useful, but WITHOUT
9  * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
10  * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
11  * License for more details.
12  *
13  * You should have received a copy of the GNU Lesser General Public
14  * License along with MoFEM. If not, see <http://www.gnu.org/licenses/>
15  */
16 
17 #ifndef __CUTMESHINTERFACE_HPP__
18 #define __CUTMESHINTERFACE_HPP__
19 
20 #include "UnknownInterface.hpp"
21 
22 namespace MoFEM {
23 
26 
27 /**
28  * \brief Interface to cut meshes
29  *
30  * \ingroup mesh_cut
31  */
33 
35  UnknownInterface **iface) const;
36 
38  CutMeshInterface(const MoFEM::Core &core);
40 
44 
46 
47  /**
48  * \brief Get options from command line
49  * @return error code
50  */
53  CHKERR PetscOptionsBegin(PETSC_COMM_WORLD, "cut_", "MOFEM Cut mesh options",
54  "none");
55 
56  CHKERR PetscOptionsInt("-linesearch_steps",
57  "number of bisection steps which line search do to "
58  "find optimal merged nodes position",
59  "", lineSearchSteps, &lineSearchSteps, PETSC_NULL);
60 
61  CHKERR PetscOptionsInt("-max_merging_cycles",
62  "number of maximal merging cycles", "",
64 
65  CHKERR PetscOptionsScalar("-project_entities_quality_trashold",
66  "project entities quality trashold", "",
68  &projectEntitiesQualityTrashold, PETSC_NULL);
69 
70  ierr = PetscOptionsEnd();
71  CHKERRG(ierr);
73  }
74 
75  MoFEMErrorCode setFront(const Range &surface);
76 
77  /**
78  * \brief set surface entities
79  * @param surface entities which going to be added
80  * @return error code
81  */
82  MoFEMErrorCode setSurface(const Range &surface);
83 
84  /**
85  * \brief copy surface entities
86  * @param surface entities which going to be added
87  * @return error code
88  */
89  MoFEMErrorCode copySurface(const Range &surface, Tag th = NULL,
90  double *shift = NULL, double *origin = NULL,
91  double *transform = NULL,
92  const std::string save_mesh = "");
93 
94  /**
95  * \brief set volume entities
96  * @param volume entities which going to be added
97  * @return error code
98  */
99  MoFEMErrorCode setVolume(const Range &volume);
100 
101  /**
102  * \brief merge surface entities
103  * @param surface entities which going to be added
104  * @return error code
105  */
106  MoFEMErrorCode mergeSurface(const Range &surface);
107 
108  /**
109  * \brief merge volume entities
110  * @param volume entities which going to be added
111  * @return error code
112  */
113  MoFEMErrorCode mergeVolumes(const Range &volume);
114 
115  MoFEMErrorCode snapSurfaceSkinToEdges(const Range &fixed_edges,
116  const double rel_tol,
117  const double abs_tol, Tag th = nullptr,
118  const bool debug = false);
119 
120  MoFEMErrorCode snapSurfaceToEdges(const Range &surface_edges,
121  const Range &fixed_edges,
122  const double rel_tol, const double abs_tol,
123  Tag th = nullptr, const bool debug = false);
124 
125  /**
126  * \brief build tree
127  * @return error code
128  */
130 
131  /**
132  * @brief Cut mesh onlu
133  *
134  * @param vol
135  * @param cut_bit
136  * @param th
137  * @param tol_cut
138  * @param tol_cut_close
139  * @param fixed_edges
140  * @param corner_nodes
141  * @param update_meshsets
142  * @param debug
143  * @return MoFEMErrorCode
144  */
145  MoFEMErrorCode cutOnly(Range vol, const BitRefLevel cut_bit, Tag th,
146  const double tol_cut, const double tol_cut_close,
147  Range *fixed_edges = NULL, Range *corner_nodes = NULL,
148  const bool update_meshsets = false,
149  const bool debug = false);
150 
151  /**
152  * @brief Trim mesh only
153  *
154  * @param trim_bit
155  * @param th
156  * @param tol_cut_close
157  * @param fixed_edges
158  * @param corner_nodes
159  * @param update_meshsets
160  * @param debug
161  * @return MoFEMErrorCode
162  */
163  MoFEMErrorCode trimOnly(const BitRefLevel trim_bit, Tag th,
164  const double tol_cut_close, Range *fixed_edges = NULL,
165  Range *corner_nodes = NULL,
166  const bool update_meshsets = false,
167  const bool debug = false);
168 
169  /**
170  * @brief Cut and trim
171  *
172  * @param first_bit
173  * @param th
174  * @param tol_cut
175  * @param tol_cut_close
176  * @param tol_trim_close
177  * @param fixed_edges
178  * @param corner_nodes
179  * @param update_meshsets
180  * @param debug
181  * @return MoFEMErrorCode
182  */
184  cutAndTrim(int &first_bit, Tag th, const double tol_cut,
185  const double tol_cut_close, const double tol_trim_close,
186  Range *fixed_edges = NULL, Range *corner_nodes = NULL,
187  const bool update_meshsets = false, const bool debug = false);
188 
189  /**
190  * @brief Cut, trim and merge
191  *
192  * @param first_bit first bit of bit revel, subsequent set bits are for trim
193  * and merge
194  * @param fraction_level fraction of edges merged at each merge step
195  * @param th tag storring mesh node positions
196  * @param tol_cut tolerance how mesh node should be close to cut surface (mesh
197  * node is moved), should be small
198  * @param tol_cut_close how crack node should be close to mesh (cut surface
199  * node is moved), can be big
200  * @param tol_trim_close how front node should be close to mesh, can be big
201  * @param fixed_edges edges on which nodes can not be moved
202  * @param corner_nodes nodes which can not be moved
203  * @param update_meshsets update meshsets by parents
204  * @param debug swich on debugging
205  * @return MoFEMErrorCode
206  */
207  MoFEMErrorCode cutTrimAndMerge(int &first_bit, const int fraction_level,
208  Tag th, const double tol_cut,
209  const double tol_cut_close,
210  const double tol_trim_close,
211  Range &fixed_edges, Range &corner_nodes,
212  const bool update_meshsets = false,
213  const bool debug = false);
214 
215  /**
216  * @brief Create front from the surface
217  *
218  * @param debug
219  * @return MoFEMErrorCode
220  */
221  MoFEMErrorCode makeFront(const bool debug = false);
222 
223  /**
224  * @brief Calculate distance from mesh nodes to cut surface
225  *
226  * @param intersect_vol
227  * @param verb
228  * @param debug
229  * @return MoFEMErrorCode
230  */
232  const bool debug = false);
233 
234  /**
235  * @brief Calculate distance from mesh nodes to surface front
236  *
237  * @param vol
238  * @param th if not null take coordinates from tag
239  * @param verb
240  * @param debug
241  * @return MoFEMErrorCode
242  */
243  MoFEMErrorCode createFrontLevelSets(Range vol, Tag th = nullptr,
244  int verb = QUIET,
245  const bool debug = false);
246 
247  /**
248  * @brief Create a level sets, i.e. distances from surface and surface front
249  *
250  * @param th
251  * @param vol_edges
252  * @param remove_adj_prims_edges
253  * @param verb
254  * @param debug
255  * @param edges_file_name
256  * @return MoFEMErrorCode
257  */
258  MoFEMErrorCode createLevelSets(Tag th, Range &vol_edges,
259  const bool remove_adj_prims_edges,
260  int verb = QUIET, const bool debug = false,
261  const std::string edges_file_name = string());
262 
263  /**
264  * @brief Create a level sets, i.e. distances from surface and surface front
265  *
266  * @param update_front
267  * @param verb
268  * @param debug
269  * @return MoFEMErrorCode
270  */
271  MoFEMErrorCode createLevelSets(int verb = QUIET, const bool debug = false);
272 
273  /**
274  * @brief Refine and set level sets
275  *
276  * \note Should be run befor cutting
277  *
278  * @param refine_front refine nodes at front
279  * @param update_front update level set at front
280  * @param init_bit_level inital bit ref level to store refined meshes
281  * @param surf_levels number of mesh surface refinement
282  * @param front_levels
283  * @param fixed_edges
284  * @param verb
285  * @param debug
286  * @return MoFEMErrorCode
287  */
288  MoFEMErrorCode refineMesh(const int init_bit_level, const int surf_levels,
289  const int front_levels,
290  Range *fixed_edges = nullptr, int verb = QUIET,
291  const bool debug = false);
292 
293  /**
294  * @brief find edges to cut
295  *
296  * @param vol is tetrahedrons search to cut
297  * @param fixed_edges pointer to fixed edges
298  * @param corner_nodes pointer to corner edges
299  * @param geometry_tol tolerance for close geometry fetchers
300  * @param verb verbosity level
301  * @param debug debugging
302  * @return MoFEMErrorCode
303  */
304  MoFEMErrorCode findEdgesToCut(Range vol, Range *fixed_edges,
305  Range *corner_nodes, const double geometry_tol,
306  int verb = QUIET, const bool debug = false);
307 
308  /**
309  * @brief Find entities on cut surface which can be projected
310  *
311  * @param fixed_edges pointer to fix edges
312  * @param corner_nodes pointer to corner nodes
313  * @param close tolerance is tolerance how close entities has to be
314  * @param verb verbosity level
315  * @param debug true for debuging purposes
316  *
317  */
318  MoFEMErrorCode projectZeroDistanceEnts(Range *fixed_edges,
319  Range *corner_nodes,
320  const double close_tol = 0,
321  const int verb = QUIET,
322  const bool debug = false);
323 
324  /**
325  * \brief cut edges
326  *
327  * For edges to cut (calculated by findEdgesToCut), edges are split in the
328  * middle and then using MoFEM::MeshRefinement interface, tetrahedra mesh
329  * are cut.
330  *
331  * @param bit BitRefLevel of new mesh created by cutting edges
332  * @return error code
333  */
334  MoFEMErrorCode cutEdgesInMiddle(const BitRefLevel bit, Range &cut_vols,
335  Range &cut_surf, Range &cut_verts,
336  const bool debug = false);
337 
338  /**
339  * \brief projecting of mid edge nodes on new mesh on surface
340  * @return error code
341  */
342  MoFEMErrorCode moveMidNodesOnCutEdges(Tag th = NULL);
343 
344  /**
345  * \brief Find edges to trimEdges
346 
347  * To make this work, you need to find edges to cut (findEdgesToCut), then
348  * cut edges in the middle (cutEdgesInMiddle) and finally project edges on
349  * the surface (moveMidNodesOnCutEdges)
350 
351  * @param verb verbosity level
352  * @return error code
353  */
354  MoFEMErrorCode findEdgesToTrim(Range *fixed_edges, Range *corner_nodes,
355  Tag th = NULL, const double tol = 1e-4,
356  const bool debug = false);
357 
358  /**
359  * \brief trim edges
360  * @param bit bit level of the trimmed mesh
361  * @return error code
362  */
364  const bool debug = false);
365 
366  /**
367  * \brief move trimmed edges mid nodes
368  * @return error code
369  */
371 
372  /**
373  * @brief Trim surface from faces beyond front
374  *
375  * @param fixed_edge
376  * @param corner_nodes
377  * @param debug
378  * @return MoFEMErrorCode
379  */
380  MoFEMErrorCode trimSurface(Range *fixed_edge, Range *corner_nodes,
381  const bool debug = false);
382 
383  /**
384  * \brief Remove pathological elements on surface internal front
385  *
386  * Internal surface skin is a set of edges in iterior of the body on boundary
387  * of surface. This set of edges is called surface front. If surface face has
388  * three nodes on surface front, non of the face nodes is split and should be
389  * removed from surface if it is going to be split.
390  *
391  * @param split_bit split bit level
392  * @param bit bit level of split mesh
393  * @param ents ents on the surface which is going to be split
394  * @return error code
395  */
397  Range &ents);
398 
399  /**
400  * \brief split sides
401  * @param split_bit split bit level
402  * @param bit bit level of split mesh
403  * @param ents ents on the surface which is going to be split
404  * @return error code
405  */
406  MoFEMErrorCode splitSides(const BitRefLevel split_bit, const BitRefLevel bit,
407  const Range &ents, Tag th = NULL);
408 
409  /**
410  * @brief Merge edges
411  *
412  * Sort all edges, where sorting is by quality calculated as edge length times
413  * quality of tets adjacent to the edge. Edge is merged if quality if the mesh
414  * is improved.
415  *
416  * @param fraction_level Fraction of edges attemt to be merged at iteration
417  * @param tets Tets of the mesh which edges are merged
418  * @param surface Surface created by edge spliting
419  * @param fixed_edges edges which are geometrical corners of the body
420  * @param corner_nodes vertices on the corners
421  * @param merged_nodes merged nodes
422  * @param out_tets returned test after merge
423  * @param new_surf new surface without merged edges
424  * @param th tag with nodal positons
425  * @param bit_ptr set bit ref level to mesh without merged edges
426  * @param debug
427  * @return MoFEMErrorCode
428  */
429  MoFEMErrorCode mergeBadEdges(const int fraction_level, const Range &tets,
430  const Range &surface, const Range &fixed_edges,
431  const Range &corner_nodes, Range &merged_nodes,
432  Range &out_tets, Range &new_surf, Tag th,
433  const bool update_meshsets = false,
434  const BitRefLevel *bit_ptr = NULL,
435  const bool debug = false);
436 
437  /**
438  * @brief Merge edges
439  *
440  * Sort all edges, where sorting is by quality calculated as edge length times
441  * quality of tets adjacent to the edge. Edge is merged if quality if the mesh
442  * is improved.
443  */
445  mergeBadEdges(const int fraction_level, const BitRefLevel cut_bit,
446  const BitRefLevel trim_bit, const BitRefLevel bit,
447  const Range &surface, const Range &fixed_edges,
448  const Range &corner_nodes, Tag th = NULL,
449  const bool update_meshsets = false, const bool debug = false);
450 
451  /**
452  * \brief set coords to tag
453  * @param th tag handle
454  * @return error code
455  */
456  MoFEMErrorCode setTagData(Tag th, const BitRefLevel bit = BitRefLevel());
457 
458  /**
459  * \brief set coords from tag
460  * @param th tag handle
461  * @return error code
462  */
463  MoFEMErrorCode setCoords(Tag th, const BitRefLevel bit = BitRefLevel(),
464  const BitRefLevel mask = BitRefLevel().set());
465 
466  inline const Range &getVolume() const { return vOlume; }
467  inline const Range &getSurface() const { return sUrface; }
468  inline const Range &getFront() const { return fRont; }
469 
470  inline const Range &getCutEdges() const { return cutEdges; }
471  inline const Range &getNewCutVolumes() const { return cutNewVolumes; }
472  inline const Range &getNewCutSurfaces() const { return cutNewSurfaces; }
473  inline const Range &getNewCutVertices() const { return cutNewVertices; }
474  inline const Range &projectZeroDistanceEnts() const {
475  return zeroDistanceEnts;
476  }
477 
478  inline const Range &getTrimEdges() const { return trimEdges; }
479  inline const Range &getNewTrimVolumes() const { return trimNewVolumes; }
480  inline const Range &getNewTrimSurfaces() const { return trimNewSurfaces; }
481  inline const Range &getNewTrimVertices() const { return trimNewVertices; }
482 
483  inline const Range &getMergedVolumes() const { return mergedVolumes; }
484  inline const Range &getMergedSurfaces() const { return mergedSurfaces; }
485 
486  inline const Range &getCutSurfaceVolumes() const { return cutSurfaceVolumes; }
487  inline const Range &getCutFrontVolumes() const { return cutFrontVolumes; }
488 
489  inline void setTrimFixedEdges(const bool b) { trimFixedEdges = b; };
490 
491  MoFEMErrorCode saveCutEdges(const std::string prefix = "");
492 
494 
495  inline boost::shared_ptr<OrientedBoxTreeTool> &getTreeSurfPtr() {
496  return treeSurfPtr;
497  }
498 
500 
501  struct SaveData {
502  moab::Interface &moab;
503  SaveData(moab::Interface &moab) : moab(moab) {}
504  MoFEMErrorCode operator()(const std::string name, const Range &ents) {
506  EntityHandle meshset;
507  CHKERR moab.create_meshset(MESHSET_SET, meshset);
508  CHKERR moab.add_entities(meshset, ents);
509  CHKERR moab.write_file(name.c_str(), "VTK", "", &meshset, 1);
510  CHKERR moab.delete_entities(&meshset, 1);
512  }
513  };
514 
515  struct LengthMapData {
516  double lEngth;
517  double qUality;
519  bool skip;
520  LengthMapData(const double l, double q, const EntityHandle e)
521  : lEngth(l), qUality(q), eDge(e), skip(false) {}
522  };
523 
524  typedef multi_index_container<
525  LengthMapData,
526  indexed_by<
527 
528  ordered_non_unique<
529  member<LengthMapData, double, &LengthMapData::lEngth>>,
530 
531  hashed_unique<
532  member<LengthMapData, EntityHandle, &LengthMapData::eDge>>,
533 
534  ordered_non_unique<
535  member<LengthMapData, double, &LengthMapData::qUality>>
536 
537  >>
539 
540 private:
541  Range fRont;
542  Range sUrface;
543  Range vOlume;
544 
546 
547  boost::shared_ptr<OrientedBoxTreeTool> treeSurfPtr;
549 
550  Range cutEdges;
556 
560  Range trimEdges;
561 
564 
565  struct TreeData {
566  double dIst;
567  double lEngth;
570  };
571 
572  map<EntityHandle, TreeData> edgesToCut;
573  map<EntityHandle, TreeData> verticesOnCutEdges;
574  map<EntityHandle, TreeData> edgesToTrim;
575  map<EntityHandle, TreeData> verticesOnTrimEdges;
576 
577  double aveLength; ///< Average edge length
578  double maxLength; ///< Maximal edge length
579 
582 };
583 } // namespace MoFEM
584 
585 #endif //__CUTMESHINTERFACE_HPP__
586 
587 /**
588  * \defgroup mesh_cut CutMeshInterface
589  * \brief Interface to mesh cut mesh
590  *
591  * \ingroup mofem
592  */
boost::shared_ptr< OrientedBoxTreeTool > & getTreeSurfPtr()
MoFEMErrorCode buildTree()
build tree
MoFEMErrorCode trimSurface(Range *fixed_edge, Range *corner_nodes, const bool debug=false)
Trim surface from faces beyond front.
map< EntityHandle, TreeData > edgesToTrim
MoFEM interface unique ID.
Interface to cut meshes.
const Range & getMergedVolumes() const
MoFEMErrorCode setTagData(Tag th, const BitRefLevel bit=BitRefLevel())
set coords to tag
MoFEMErrorCode trimEdgesInTheMiddle(const BitRefLevel bit, const bool debug=false)
trim edges
MoFEMErrorCode getSubInterfaceOptions()
double maxLength
Maximal edge length.
MoFEMErrorCode getOptions()
Get options from command line.
MoFEMErrorCode copySurface(const Range &surface, Tag th=NULL, double *shift=NULL, double *origin=NULL, double *transform=NULL, const std::string save_mesh="")
copy surface entities
const Range & getNewCutVertices() const
MoFEMErrorCode snapSurfaceToEdges(const Range &surface_edges, const Range &fixed_edges, const double rel_tol, const double abs_tol, Tag th=nullptr, const bool debug=false)
MoFEMErrorCode setSurface(const Range &surface)
set surface entities
map< EntityHandle, TreeData > verticesOnTrimEdges
const Range & getSurface() const
const Range & getNewTrimVertices() const
MoFEMErrorCode makeFront(const bool debug=false)
Create front from the surface.
base class for all interface classes
boost::shared_ptr< OrientedBoxTreeTool > treeSurfPtr
MoFEMErrorCode snapSurfaceSkinToEdges(const Range &fixed_edges, const double rel_tol, const double abs_tol, Tag th=nullptr, const bool debug=false)
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
MoFEMErrorCode setFront(const Range &surface)
MoFEMErrorCode mergeBadEdges(const int fraction_level, const Range &tets, const Range &surface, const Range &fixed_edges, const Range &corner_nodes, Range &merged_nodes, Range &out_tets, Range &new_surf, Tag th, const bool update_meshsets=false, const BitRefLevel *bit_ptr=NULL, const bool debug=false)
Merge edges.
MoFEMErrorCode createSurfaceLevelSets(int verb=QUIET, const bool debug=false)
Calculate distance from mesh nodes to cut surface.
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:544
MoFEMErrorCode cutTrimAndMerge(int &first_bit, const int fraction_level, Tag th, const double tol_cut, const double tol_cut_close, const double tol_trim_close, Range &fixed_edges, Range &corner_nodes, const bool update_meshsets=false, const bool debug=false)
Cut, trim and merge.
void setTrimFixedEdges(const bool b)
double aveLength
Average edge length.
MoFEMErrorCode operator()(const std::string name, const Range &ents)
const Range & getCutSurfaceVolumes() const
const Range & getNewCutVolumes() const
Core (interface) class.
Definition: Core.hpp:50
MoFEM interface.
const Range & getMergedSurfaces() const
MoFEMErrorCode removePathologicalFrontTris(const BitRefLevel split_bit, Range &ents)
Remove pathological elements on surface internal front.
MoFEMErrorCode cutEdgesInMiddle(const BitRefLevel bit, Range &cut_vols, Range &cut_surf, Range &cut_verts, const bool debug=false)
cut edges
std::bitset< BITINTERFACEUID_SIZE > BitIntefaceId
Definition: Types.hpp:56
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
LengthMapData(const double l, double q, const EntityHandle e)
const Range & getTrimEdges() const
CutMeshInterface(const MoFEM::Core &core)
const Range & getFront() const
const Range & projectZeroDistanceEnts() const
MoFEMErrorCode trimOnly(const BitRefLevel trim_bit, Tag th, const double tol_cut_close, Range *fixed_edges=NULL, Range *corner_nodes=NULL, const bool update_meshsets=false, const bool debug=false)
Trim mesh only.
MoFEMErrorCode saveCutEdges(const std::string prefix="")
Tensor2_Expr< transform_Tensor2< A, B, T, Dim0, Dim1, i, j >, T, Dim0, Dim1, i, j > transform(const Tensor2_Expr< A, T, Dim0, Dim1, i, j > &a, B function)
const Range & getNewTrimSurfaces() const
double tol
const Range & getNewTrimVolumes() const
const Range & getCutEdges() const
MoFEMErrorCode cutAndTrim(int &first_bit, Tag th, const double tol_cut, const double tol_cut_close, const double tol_trim_close, Range *fixed_edges=NULL, Range *corner_nodes=NULL, const bool update_meshsets=false, const bool debug=false)
Cut and trim.
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:66
static const bool debug
const Range & getVolume() const
multi_index_container< LengthMapData, indexed_by< ordered_non_unique< member< LengthMapData, double, &LengthMapData::lEngth > >, hashed_unique< member< LengthMapData, EntityHandle, &LengthMapData::eDge > >, ordered_non_unique< member< LengthMapData, double, &LengthMapData::qUality > > > > LengthMapData_multi_index
MoFEMErrorCode setCoords(Tag th, const BitRefLevel bit=BitRefLevel(), const BitRefLevel mask=BitRefLevel().set())
set coords from tag
MoFEMErrorCode mergeSurface(const Range &surface)
merge surface entities
MoFEMErrorCode moveMidNodesOnCutEdges(Tag th=NULL)
projecting of mid edge nodes on new mesh on surface
MoFEMErrorCode mergeVolumes(const Range &volume)
merge volume entities
MoFEMErrorCode splitSides(const BitRefLevel split_bit, const BitRefLevel bit, const Range &ents, Tag th=NULL)
split sides
MoFEMErrorCode query_interface(const MOFEMuuid &uuid, UnknownInterface **iface) const
const Range & getCutFrontVolumes() const
#define CHKERR
Inline error check.
Definition: definitions.h:596
static const MOFEMuuid IDD_MOFEMCutMesh
MoFEMErrorCode findEdgesToTrim(Range *fixed_edges, Range *corner_nodes, Tag th=NULL, const double tol=1e-4, const bool debug=false)
Find edges to trimEdges.
MoFEMErrorCode refineMesh(const int init_bit_level, const int surf_levels, const int front_levels, Range *fixed_edges=nullptr, int verb=QUIET, const bool debug=false)
Refine and set level sets.
VectorBoundedArray< double, 3 > VectorDouble3
Definition: Types.hpp:85
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Types.hpp:51
MoFEMErrorCode createLevelSets(Tag th, Range &vol_edges, const bool remove_adj_prims_edges, int verb=QUIET, const bool debug=false, const std::string edges_file_name=string())
Create a level sets, i.e. distances from surface and surface front.
MoFEMErrorCode findEdgesToCut(Range vol, Range *fixed_edges, Range *corner_nodes, const double geometry_tol, int verb=QUIET, const bool debug=false)
find edges to cut
map< EntityHandle, TreeData > verticesOnCutEdges
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407
MoFEMErrorCode saveTrimEdges()
const Range & getNewCutSurfaces() const
MoFEMErrorCode createFrontLevelSets(Range vol, Tag th=nullptr, int verb=QUIET, const bool debug=false)
Calculate distance from mesh nodes to surface front.
map< EntityHandle, TreeData > edgesToCut
MoFEMErrorCode setVolume(const Range &volume)
set volume entities
MoFEMErrorCode cutOnly(Range vol, const BitRefLevel cut_bit, Tag th, const double tol_cut, const double tol_cut_close, Range *fixed_edges=NULL, Range *corner_nodes=NULL, const bool update_meshsets=false, const bool debug=false)
Cut mesh onlu.
MoFEMErrorCode moveMidNodesOnTrimmedEdges(Tag th=NULL)
move trimmed edges mid nodes