v0.10.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  /**
103  * @brief Set the constrain surface object
104  *
105  * Add surfaces which are restricted by mesh cut. Example of surface which
106  * needs to be respected is an interface, the boundary between two materials,
107  * or crack surface.
108  *
109  * @param surf
110  * @return MoFEMErrorCode
111  */
112  MoFEMErrorCode setConstrainSurface(const Range surf);
113 
114  /**
115  * \brief merge surface entities
116  * @param surface entities which going to be added
117  * @return error code
118  */
119  MoFEMErrorCode mergeSurface(const Range surface);
120 
121  /**
122  * \brief merge volume entities
123  * @param volume entities which going to be added
124  * @return error code
125  */
126  MoFEMErrorCode mergeVolumes(const Range volume);
127 
128  MoFEMErrorCode snapSurfaceSkinToEdges(const Range fixed_edges,
129  const double rel_tol,
130  const double abs_tol, Tag th = nullptr,
131  const bool debug = false);
132 
133  MoFEMErrorCode snapSurfaceToEdges(const Range surface_edges,
134  const Range fixed_edges,
135  const double rel_tol, const double abs_tol,
136  Tag th = nullptr, const bool debug = false);
137 
138  /**
139  * \brief build tree
140  * @return error code
141  */
143 
144  /**
145  * @brief Cut mesh only
146  *
147  * @param vol
148  * @param cut_bit
149  * @param th
150  * @param tol_geometry
151  * @param tol_cut_close
152  * @param fixed_edges
153  * @param corner_nodes
154  * @param update_meshsets
155  * @param debug
156  * @return MoFEMErrorCode
157  */
158  MoFEMErrorCode cutOnly(Range vol, const BitRefLevel cut_bit, Tag th,
159  const double tol_geometry, const double tol_cut_close,
160  Range *fixed_edges = NULL, Range *corner_nodes = NULL,
161  const bool update_meshsets = false,
162  const bool debug = false);
163 
164  /**
165  * @brief Trim mesh only
166  *
167  * @param trim_bit
168  * @param th
169  * @param tol_cut_close
170  * @param fixed_edges
171  * @param corner_nodes
172  * @param update_meshsets
173  * @param debug
174  * @return MoFEMErrorCode
175  */
176  MoFEMErrorCode trimOnly(const BitRefLevel trim_bit, Tag th,
177  const double tol_trim_close,
178  Range *fixed_edges = NULL, Range *corner_nodes = NULL,
179  const bool update_meshsets = false,
180  const bool debug = false);
181 
182  /**
183  * @brief Cut and trim
184  *
185  * @param first_bit
186  * @param th
187  * @param tol_geometry
188  * @param tol_cut_close
189  * @param tol_trim_close
190  * @param fixed_edges
191  * @param corner_nodes
192  * @param update_meshsets
193  * @param debug
194  * @return MoFEMErrorCode
195  */
197  cutAndTrim(int &first_bit, Tag th, const double tol_geometry,
198  const double tol_cut_close, const double tol_trim_close,
199  Range *fixed_edges = NULL, Range *corner_nodes = NULL,
200  const bool update_meshsets = false, const bool debug = false);
201 
202  /**
203  * @brief Cut, trim and merge
204  *
205  * @param first_bit first bit of bit revel, subsequent set bits are for trim
206  * and merge
207  * @param fraction_level fraction of edges merged at each merge step
208  * @param th tag storring mesh node positions
209  * @param tol_geometry tolerance how mesh node should be close to cut surface (mesh
210  * node is moved), should be small
211  * @param tol_cut_close how crack node should be close to mesh (cut surface
212  * node is moved), can be big
213  * @param tol_trim_close how front node should be close to mesh, can be big
214  * @param fixed_edges edges on which nodes can not be moved
215  * @param corner_nodes nodes which can not be moved
216  * @param update_meshsets update meshsets by parents
217  * @param debug swich on debugging
218  * @return MoFEMErrorCode
219  */
220  MoFEMErrorCode cutTrimAndMerge(int &first_bit, const int fraction_level,
221  Tag th, const double tol_geometry,
222  const double tol_cut_close,
223  const double tol_trim_close,
224  Range &fixed_edges, Range &corner_nodes,
225  const bool update_meshsets = false,
226  const bool debug = false);
227 
228  /**
229  * @brief Create front from the surface
230  *
231  * @param debug
232  * @return MoFEMErrorCode
233  */
234  MoFEMErrorCode makeFront(const bool debug = false);
235 
236  /**
237  * @brief Calculate distance from mesh nodes to cut surface
238  *
239  * @param intersect_vol
240  * @param verb
241  * @param debug
242  * @return MoFEMErrorCode
243  */
245  const bool debug = false);
246 
247  /**
248  * @brief Calculate distance from mesh nodes to surface front
249  *
250  * @param vol
251  * @param th if not null take coordinates from tag
252  * @param verb
253  * @param debug
254  * @return MoFEMErrorCode
255  */
256  MoFEMErrorCode createFrontLevelSets(Range vol, Tag th = nullptr,
257  int verb = QUIET,
258  const bool debug = false);
259 
260  /**
261  * @brief Find level set volumes
262  *
263  * @param th
264  * @param vol_edges
265  * @param remove_adj_prims_edges
266  * @param verb
267  * @param debug
268  * @param edges_file_name
269  * @return MoFEMErrorCode
270  */
271  MoFEMErrorCode findLevelSetVolumes(Tag th, Range &vol_edges,
272  const bool remove_adj_prims_edges,
273  int verb = QUIET, const bool debug = false,
274  const std::string edges_file_name = string());
275 
276  /**
277  * @brief Find level set volumes
278  *
279  * @param update_front
280  * @param verb
281  * @param debug
282  * @return MoFEMErrorCode
283  */
285  const bool debug = false);
286 
287  /**
288  * @brief Refine and set level sets
289  *
290  * \note Should be run befor cutting
291  *
292  * @param refine_front refine nodes at front
293  * @param update_front update level set at front
294  * @param init_bit_level inital bit ref level to store refined meshes
295  * @param surf_levels number of mesh surface refinement
296  * @param front_levels
297  * @param fixed_edges
298  * @param verb
299  * @param debug
300  * @return MoFEMErrorCode
301  */
302  MoFEMErrorCode refineMesh(const int init_bit_level, const int surf_levels,
303  const int front_levels,
304  Range *fixed_edges = nullptr, int verb = QUIET,
305  const bool debug = false);
306 
307  /**
308  * @brief find edges to cut
309  *
310  * @param vol is tetrahedrons search to cut
311  * @param verb verbosity level
312  * @param debug debugging
313  * @return MoFEMErrorCode
314  */
315  MoFEMErrorCode findEdgesToCut(Range vol, int verb = QUIET,
316  const bool debug = false);
317 
318  /**
319  * @brief Find entities on cut surface which can be projected
320  *
321  * @param fixed_edges pointer to fix edges
322  * @param corner_nodes pointer to corner nodes
323  * @param close tolerance is tolerance how close entities has to be
324  * @param verb verbosity level
325  * @param debug true for debuging purposes
326  *
327  */
328  MoFEMErrorCode projectZeroDistanceEnts(Range *fixed_edges,
329  Range *corner_nodes,
330  const double geometry_tol = 0,
331  const double close_tol = 0,
332  const int verb = QUIET,
333  const bool debug = false);
334 
335  /**
336  * \brief cut edges
337  *
338  * For edges to cut (calculated by findEdgesToCut), edges are split in the
339  * middle and then using MoFEM::MeshRefinement interface, tetrahedra mesh
340  * are cut.
341  *
342  * @param bit BitRefLevel of new mesh created by cutting edges
343  * @return error code
344  */
345  MoFEMErrorCode cutEdgesInMiddle(const BitRefLevel bit, Range &cut_vols,
346  Range &cut_surf, Range &cut_verts,
347  const bool debug = false);
348 
349  /**
350  * \brief projecting of mid edge nodes on new mesh on surface
351  * @return error code
352  */
353  MoFEMErrorCode moveMidNodesOnCutEdges(Tag th = NULL);
354 
355  /**
356  * \brief Find edges to trimEdges
357 
358  * To make this work, you need to find edges to cut (findEdgesToCut), then
359  * cut edges in the middle (cutEdgesInMiddle) and finally project edges on
360  * the surface (moveMidNodesOnCutEdges)
361 
362  * @param verb verbosity level
363  * @return error code
364  */
365  MoFEMErrorCode findEdgesToTrim(Range *fixed_edges, Range *corner_nodes,
366  Tag th = NULL,
367  const double tol = 1e-4,
368  const bool debug = false);
369 
370  /**
371  * \brief trim edges
372  * @param bit bit level of the trimmed mesh
373  * @return error code
374  */
376  const bool debug = false);
377 
378  /**
379  * \brief move trimmed edges mid nodes
380  * @return error code
381  */
383 
384  /**
385  * @brief Trim surface from faces beyond front
386  *
387  * @param fixed_edge
388  * @param corner_nodes
389  * @param debug
390  * @return MoFEMErrorCode
391  */
392  MoFEMErrorCode trimSurface(Range *fixed_edge, Range *corner_nodes,
393  const double tol = 1e-4, Tag th = NULL,
394  const bool debug = false);
395 
396  /**
397  * \brief Remove pathological elements on surface internal front
398  *
399  * Internal surface skin is a set of edges in iterior of the body on boundary
400  * of surface. This set of edges is called surface front. If surface face has
401  * three nodes on surface front, non of the face nodes is split and should be
402  * removed from surface if it is going to be split.
403  *
404  * @param split_bit split bit level
405  * @param bit bit level of split mesh
406  * @param ents ents on the surface which is going to be split
407  * @return error code
408  */
410  Range &ents);
411 
412  /**
413  * \brief split sides
414  * @param split_bit split bit level
415  * @param bit bit level of split mesh
416  * @param ents ents on the surface which is going to be split
417  * @return error code
418  */
419  MoFEMErrorCode splitSides(const BitRefLevel split_bit, const BitRefLevel bit,
420  const Range &ents, Tag th = NULL);
421 
422  /**
423  * @brief Merge edges
424  *
425  * Sort all edges, where sorting is by quality calculated as edge length times
426  * quality of tets adjacent to the edge. Edge is merged if quality if the mesh
427  * is improved.
428  *
429  * @param fraction_level Fraction of edges attemt to be merged at iteration
430  * @param tets Tets of the mesh which edges are merged
431  * @param surface Surface created by edge spliting
432  * @param fixed_edges edges which are geometrical corners of the body
433  * @param corner_nodes vertices on the corners
434  * @param merged_nodes merged nodes
435  * @param out_tets returned test after merge
436  * @param new_surf new surface without merged edges
437  * @param th tag with nodal positons
438  * @param bit_ptr set bit ref level to mesh without merged edges
439  * @param debug
440  * @return MoFEMErrorCode
441  */
442  MoFEMErrorCode mergeBadEdges(const int fraction_level, const Range &tets,
443  const Range &surface, const Range &fixed_edges,
444  const Range &corner_nodes, Range &merged_nodes,
445  Range &out_tets, Range &new_surf, Tag th,
446  const bool update_meshsets = false,
447  const BitRefLevel *bit_ptr = NULL,
448  const bool debug = false);
449 
450  /**
451  * @brief Merge edges
452  *
453  * Sort all edges, where sorting is by quality calculated as edge length times
454  * quality of tets adjacent to the edge. Edge is merged if quality if the mesh
455  * is improved.
456  */
458  mergeBadEdges(const int fraction_level, const BitRefLevel cut_bit,
459  const BitRefLevel trim_bit, const BitRefLevel bit,
460  const Range &surface, const Range &fixed_edges,
461  const Range &corner_nodes, Tag th = NULL,
462  const bool update_meshsets = false, const bool debug = false);
463 
464  /**
465  * \brief set coords to tag
466  * @param th tag handle
467  * @return error code
468  */
469  MoFEMErrorCode setTagData(Tag th, const BitRefLevel bit = BitRefLevel());
470 
471  /**
472  * \brief set coords from tag
473  * @param th tag handle
474  * @return error code
475  */
476  MoFEMErrorCode setCoords(Tag th, const BitRefLevel bit = BitRefLevel(),
477  const BitRefLevel mask = BitRefLevel().set());
478 
479  inline const Range &getVolume() const { return vOlume; }
480  inline const Range &getSurface() const { return sUrface; }
481  inline const Range &getFront() const { return fRont; }
482 
483  inline const Range &getCutEdges() const { return cutEdges; }
484  inline const Range &getNewCutVolumes() const { return cutNewVolumes; }
485  inline const Range &getNewCutSurfaces() const { return cutNewSurfaces; }
486  inline const Range &getNewCutVertices() const { return cutNewVertices; }
487  inline const Range &projectZeroDistanceEnts() const {
488  return zeroDistanceEnts;
489  }
490 
491  inline const Range &getTrimEdges() const { return trimEdges; }
492  inline const Range &getNewTrimVolumes() const { return trimNewVolumes; }
493  inline const Range &getNewTrimSurfaces() const { return trimNewSurfaces; }
494  inline const Range &getNewTrimVertices() const { return trimNewVertices; }
495 
496  inline const Range &getMergedVolumes() const { return mergedVolumes; }
497  inline const Range &getMergedSurfaces() const { return mergedSurfaces; }
498 
499  inline const Range &getCutSurfaceVolumes() const { return cutSurfaceVolumes; }
500  inline const Range &getCutFrontVolumes() const { return cutFrontVolumes; }
501 
502  MoFEMErrorCode saveCutEdges(const std::string prefix = "");
503 
505 
506  inline boost::shared_ptr<OrientedBoxTreeTool> &getTreeSurfPtr() {
507  return treeSurfPtr;
508  }
509 
511 
512  struct SaveData {
515  MoFEMErrorCode operator()(const std::string name, const Range &ents) {
517  EntityHandle meshset;
518  CHKERR moab.create_meshset(MESHSET_SET, meshset);
519  CHKERR moab.add_entities(meshset, ents);
520  CHKERR moab.write_file(name.c_str(), "VTK", "", &meshset, 1);
521  CHKERR moab.delete_entities(&meshset, 1);
523  }
524  };
525 
526  struct LengthMapData {
527  double lEngth;
528  double qUality;
530  bool skip;
531  LengthMapData(const double l, double q, const EntityHandle e)
532  : lEngth(l), qUality(q), eDge(e), skip(false) {}
533  };
534 
535  typedef multi_index_container<
536  LengthMapData,
537  indexed_by<
538 
539  ordered_non_unique<
540  member<LengthMapData, double, &LengthMapData::lEngth>>,
541 
542  hashed_unique<
543  member<LengthMapData, EntityHandle, &LengthMapData::eDge>>,
544 
545  ordered_non_unique<
546  member<LengthMapData, double, &LengthMapData::qUality>>
547 
548  >>
550 
551 private:
552  Range fRont;
553  Range sUrface;
554  Range vOlume;
555 
556  boost::shared_ptr<OrientedBoxTreeTool> treeSurfPtr;
558 
559  Range cutEdges;
565 
569  Range trimEdges;
570 
573 
574  struct TreeData {
575  double dIst;
576  double lEngth;
579  };
580 
581  map<EntityHandle, TreeData> edgesToCut;
582  map<EntityHandle, TreeData> verticesOnCutEdges;
583  map<EntityHandle, TreeData> edgesToTrim;
584  map<EntityHandle, TreeData> verticesOnTrimEdges;
585 
586  double aveLength; ///< Average edge length
587  double maxLength; ///< Maximal edge length
588 
592 };
593 } // namespace MoFEM
594 
595 #endif //__CUTMESHINTERFACE_HPP__
596 
597 /**
598  * \defgroup mesh_cut CutMeshInterface
599  * \brief Interface to mesh cut mesh
600  *
601  * \ingroup mofem
602  */
MoFEM::IDD_MOFEMCutMesh
static const MOFEMuuid IDD_MOFEMCutMesh
Definition: CutMeshInterface.hpp:24
MoFEM::CutMeshInterface::nbMaxMergingCycles
int nbMaxMergingCycles
Definition: CutMeshInterface.hpp:42
EntityHandle
MoFEM::CutMeshInterface::moveMidNodesOnTrimmedEdges
MoFEMErrorCode moveMidNodesOnTrimmedEdges(Tag th=NULL)
move trimmed edges mid nodes
Definition: CutMeshInterface.cpp:1484
MoFEM::CutMeshInterface::getTreeSurfPtr
boost::shared_ptr< OrientedBoxTreeTool > & getTreeSurfPtr()
Definition: CutMeshInterface.hpp:506
MoFEM::CutMeshInterface::LengthMapData_multi_index
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
Definition: CutMeshInterface.hpp:549
MoFEM::CutMeshInterface::setCoords
MoFEMErrorCode setCoords(Tag th, const BitRefLevel bit=BitRefLevel(), const BitRefLevel mask=BitRefLevel().set())
set coords from tag
Definition: CutMeshInterface.cpp:3149
MoFEM::CutMeshInterface::cutNewVertices
Range cutNewVertices
Definition: CutMeshInterface.hpp:564
MoFEM::CutMeshInterface::getSurface
const Range & getSurface() const
Definition: CutMeshInterface.hpp:480
MoFEM::CutMeshInterface::copySurface
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
Definition: CutMeshInterface.cpp:64
MoFEM::CutMeshInterface::vOlume
Range vOlume
Definition: CutMeshInterface.hpp:554
MoFEM::Exceptions::MoFEMErrorCode
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:67
MoFEM::CutMeshInterface::cutSurfaceVolumes
Range cutSurfaceVolumes
Definition: CutMeshInterface.hpp:589
MoFEM::CutMeshInterface::getNewCutSurfaces
const Range & getNewCutSurfaces() const
Definition: CutMeshInterface.hpp:485
MoFEM::CutMeshInterface::getMergedVolumes
const Range & getMergedVolumes() const
Definition: CutMeshInterface.hpp:496
MoFEM::CutMeshInterface::snapSurfaceToEdges
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)
Definition: CutMeshInterface.cpp:168
MoFEM::CutMeshInterface::trimSurface
MoFEMErrorCode trimSurface(Range *fixed_edge, Range *corner_nodes, const double tol=1e-4, Tag th=NULL, const bool debug=false)
Trim surface from faces beyond front.
Definition: CutMeshInterface.cpp:2047
MoFEM::CutMeshInterface
Interface to cut meshes.
Definition: CutMeshInterface.hpp:32
MoFEM::CutMeshInterface::cutNewVolumes
Range cutNewVolumes
Definition: CutMeshInterface.hpp:560
MoFEM::CutMeshInterface::zeroDistanceEnts
Range zeroDistanceEnts
Definition: CutMeshInterface.hpp:562
MoFEM::CutMeshInterface::rootSetSurf
EntityHandle rootSetSurf
Definition: CutMeshInterface.hpp:557
MoFEM::CutMeshInterface::saveCutEdges
MoFEMErrorCode saveCutEdges(const std::string prefix="")
Definition: CutMeshInterface.cpp:3166
MoFEM::CutMeshInterface::LengthMapData::qUality
double qUality
Definition: CutMeshInterface.hpp:528
MoFEM::CutMeshInterface::findLevelSetVolumes
MoFEMErrorCode findLevelSetVolumes(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())
Find level set volumes.
Definition: CutMeshInterface.cpp:669
MoFEM::CutMeshInterface::verticesOnTrimEdges
map< EntityHandle, TreeData > verticesOnTrimEdges
Definition: CutMeshInterface.hpp:584
MoFEM::CutMeshInterface::SaveData::SaveData
SaveData(moab::Interface &moab)
Definition: CutMeshInterface.hpp:514
MoFEM::CutMeshInterface::cutOnly
MoFEMErrorCode cutOnly(Range vol, const BitRefLevel cut_bit, Tag th, const double tol_geometry, const double tol_cut_close, Range *fixed_edges=NULL, Range *corner_nodes=NULL, const bool update_meshsets=false, const bool debug=false)
Cut mesh only.
Definition: CutMeshInterface.cpp:242
MoFEM::CutMeshInterface::cutFrontVolumes
Range cutFrontVolumes
Definition: CutMeshInterface.hpp:590
MoFEM::CutMeshInterface::moveMidNodesOnCutEdges
MoFEMErrorCode moveMidNodesOnCutEdges(Tag th=NULL)
projecting of mid edge nodes on new mesh on surface
Definition: CutMeshInterface.cpp:1466
MoFEM::CutMeshInterface::getNewTrimVolumes
const Range & getNewTrimVolumes() const
Definition: CutMeshInterface.hpp:492
MoFEM::CutMeshInterface::clearMap
MoFEMErrorCode clearMap()
Definition: CutMeshInterface.cpp:46
MoFEM::CutMeshInterface::LengthMapData::skip
bool skip
Definition: CutMeshInterface.hpp:530
MoFEM::Interface
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1943
MoFEM::CutMeshInterface::makeFront
MoFEMErrorCode makeFront(const bool debug=false)
Create front from the surface.
Definition: CutMeshInterface.cpp:466
MoFEM::CutMeshInterface::findEdgesToTrim
MoFEMErrorCode findEdgesToTrim(Range *fixed_edges, Range *corner_nodes, Tag th=NULL, const double tol=1e-4, const bool debug=false)
Find edges to trimEdges.
Definition: CutMeshInterface.cpp:1499
MoFEM::CutMeshInterface::getNewTrimVertices
const Range & getNewTrimVertices() const
Definition: CutMeshInterface.hpp:494
MoFEM::CutMeshInterface::lineSearchSteps
int lineSearchSteps
Definition: CutMeshInterface.hpp:41
MoFEM::CutMeshInterface::getTrimEdges
const Range & getTrimEdges() const
Definition: CutMeshInterface.hpp:491
MoFEM::CutMeshInterface::constrainSurface
Range constrainSurface
Definition: CutMeshInterface.hpp:591
CUTMESH_INTERFACE
@ CUTMESH_INTERFACE
Definition: definitions.h:70
MoFEM::CutMeshInterface::splitSides
MoFEMErrorCode splitSides(const BitRefLevel split_bit, const BitRefLevel bit, const Range &ents, Tag th=NULL)
split sides
Definition: CutMeshInterface.cpp:2285
MoFEM::CutMeshInterface::cutTrimAndMerge
MoFEMErrorCode cutTrimAndMerge(int &first_bit, const int fraction_level, Tag th, const double tol_geometry, 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.
Definition: CutMeshInterface.cpp:395
MoFEM::CutMeshInterface::getOptions
MoFEMErrorCode getOptions()
Get options from command line.
Definition: CutMeshInterface.hpp:51
CHKERR
#define CHKERR
Inline error check.
Definition: definitions.h:604
UnknownInterface.hpp
MoFEM interface.
MoFEM::CutMeshInterface::createFrontLevelSets
MoFEMErrorCode createFrontLevelSets(Range vol, Tag th=nullptr, int verb=QUIET, const bool debug=false)
Calculate distance from mesh nodes to surface front.
Definition: CutMeshInterface.cpp:616
MoFEM::CutMeshInterface::saveTrimEdges
MoFEMErrorCode saveTrimEdges()
Definition: CutMeshInterface.cpp:3180
MoFEM::CutMeshInterface::LengthMapData::lEngth
double lEngth
Definition: CutMeshInterface.hpp:527
MoFEM::CutMeshInterface::getNewTrimSurfaces
const Range & getNewTrimSurfaces() const
Definition: CutMeshInterface.hpp:493
MoFEM::CutMeshInterface::removePathologicalFrontTris
MoFEMErrorCode removePathologicalFrontTris(const BitRefLevel split_bit, Range &ents)
Remove pathological elements on surface internal front.
Definition: CutMeshInterface.cpp:2254
MoFEM::CutMeshInterface::trimNewVolumes
Range trimNewVolumes
Definition: CutMeshInterface.hpp:566
MoFEM::CutMeshInterface::~CutMeshInterface
~CutMeshInterface()
Definition: CutMeshInterface.hpp:39
MoFEM::Types::BitIntefaceId
std::bitset< BITINTERFACEUID_SIZE > BitIntefaceId
Definition: Types.hpp:56
MoFEM::CutMeshInterface::zeroDistanceVerts
Range zeroDistanceVerts
Definition: CutMeshInterface.hpp:563
MoFEM::CutMeshInterface::LengthMapData
Definition: CutMeshInterface.hpp:526
MoFEM::CutMeshInterface::getCutFrontVolumes
const Range & getCutFrontVolumes() const
Definition: CutMeshInterface.hpp:500
MoFEM::CutMeshInterface::getSubInterfaceOptions
MoFEMErrorCode getSubInterfaceOptions()
Definition: CutMeshInterface.hpp:45
MoFEM::MOFEMuuid
MoFEM interface unique ID.
Definition: UnknownInterface.hpp:26
MoFEM::CutMeshInterface::LengthMapData::LengthMapData
LengthMapData(const double l, double q, const EntityHandle e)
Definition: CutMeshInterface.hpp:531
MoFEM::CutMeshInterface::trimEdges
Range trimEdges
Definition: CutMeshInterface.hpp:569
MoFEM::CutMeshInterface::cutAndTrim
MoFEMErrorCode cutAndTrim(int &first_bit, Tag th, const double tol_geometry, 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.
Definition: CutMeshInterface.cpp:325
MoFEM::CutMeshInterface::query_interface
MoFEMErrorCode query_interface(const MOFEMuuid &uuid, UnknownInterface **iface) const
Definition: CutMeshInterface.cpp:27
MoFEM
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
MoFEM::Types::BitRefLevel
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Types.hpp:51
MoFEM::CutMeshInterface::TreeData::rayPoint
VectorDouble3 rayPoint
Definition: CutMeshInterface.hpp:578
MoFEM::CutMeshInterface::findEdgesToCut
MoFEMErrorCode findEdgesToCut(Range vol, int verb=QUIET, const bool debug=false)
find edges to cut
Definition: CutMeshInterface.cpp:903
debug
static const bool debug
Definition: dm_create_subdm.cpp:26
MoFEM::CutMeshInterface::getVolume
const Range & getVolume() const
Definition: CutMeshInterface.hpp:479
MoFEM::CutMeshInterface::fRont
Range fRont
Definition: CutMeshInterface.hpp:552
MoFEM::CutMeshInterface::trimEdgesInTheMiddle
MoFEMErrorCode trimEdgesInTheMiddle(const BitRefLevel bit, const bool debug=false)
trim edges
Definition: CutMeshInterface.cpp:1987
MoFEM::CutMeshInterface::TreeData::dIst
double dIst
Definition: CutMeshInterface.hpp:575
MoFEM::CutMeshInterface::LengthMapData::eDge
EntityHandle eDge
Definition: CutMeshInterface.hpp:529
MoFEM::CoreTmp< 0 >
Core (interface) class.
Definition: Core.hpp:77
MoFEM::CutMeshInterface::CutMeshInterface
CutMeshInterface(const MoFEM::Core &core)
Definition: CutMeshInterface.cpp:39
MoFEM::CutMeshInterface::sUrface
Range sUrface
Definition: CutMeshInterface.hpp:553
MoFEM::CutMeshInterface::mergedSurfaces
Range mergedSurfaces
Definition: CutMeshInterface.hpp:572
MoFEM::CutMeshInterface::SaveData::moab
moab::Interface & moab
Definition: CutMeshInterface.hpp:513
MoFEM::CutMeshInterface::getNewCutVertices
const Range & getNewCutVertices() const
Definition: CutMeshInterface.hpp:486
MoFEM::CutMeshInterface::mergedVolumes
Range mergedVolumes
Definition: CutMeshInterface.hpp:571
MoFEM::CutMeshInterface::cOre
MoFEM::Core & cOre
Definition: CutMeshInterface.hpp:37
MoFEM::CutMeshInterface::snapSurfaceSkinToEdges
MoFEMErrorCode snapSurfaceSkinToEdges(const Range fixed_edges, const double rel_tol, const double abs_tol, Tag th=nullptr, const bool debug=false)
Definition: CutMeshInterface.cpp:148
MoFEM::CutMeshInterface::getNewCutVolumes
const Range & getNewCutVolumes() const
Definition: CutMeshInterface.hpp:484
MoFEM::CutMeshInterface::projectEntitiesQualityTrashold
double projectEntitiesQualityTrashold
Definition: CutMeshInterface.hpp:43
MoFEM::CutMeshInterface::TreeData
Definition: CutMeshInterface.hpp:574
MoFEM::CutMeshInterface::mergeVolumes
MoFEMErrorCode mergeVolumes(const Range volume)
merge volume entities
Definition: CutMeshInterface.cpp:142
MoFEM::CutMeshInterface::TreeData::unitRayDir
VectorDouble3 unitRayDir
Definition: CutMeshInterface.hpp:577
MoFEM::CutMeshInterface::setFront
MoFEMErrorCode setFront(const Range surface)
Definition: CutMeshInterface.cpp:52
MoFEM::CutMeshInterface::treeSurfPtr
boost::shared_ptr< OrientedBoxTreeTool > treeSurfPtr
Definition: CutMeshInterface.hpp:556
MoFEM::CutMeshInterface::refineMesh
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.
Definition: CutMeshInterface.cpp:816
MoFEM::CutMeshInterface::edgesToTrim
map< EntityHandle, TreeData > edgesToTrim
Definition: CutMeshInterface.hpp:583
MoFEM::CutMeshInterface::setSurface
MoFEMErrorCode setSurface(const Range surface)
set surface entities
Definition: CutMeshInterface.cpp:58
MoFEM::CutMeshInterface::getMergedSurfaces
const Range & getMergedSurfaces() const
Definition: CutMeshInterface.hpp:497
MoFEM::CutMeshInterface::cutEdgesInMiddle
MoFEMErrorCode cutEdgesInMiddle(const BitRefLevel bit, Range &cut_vols, Range &cut_surf, Range &cut_verts, const bool debug=false)
cut edges
Definition: CutMeshInterface.cpp:1316
MoFEM::CutMeshInterface::setVolume
MoFEMErrorCode setVolume(const Range volume)
set volume entities
Definition: CutMeshInterface.cpp:124
MoFEM::CutMeshInterface::aveLength
double aveLength
Average edge length.
Definition: CutMeshInterface.hpp:586
MoFEM::CutMeshInterface::trimNewVertices
Range trimNewVertices
Definition: CutMeshInterface.hpp:567
MoFEM::CutMeshInterface::trimOnly
MoFEMErrorCode trimOnly(const BitRefLevel trim_bit, Tag th, const double tol_trim_close, Range *fixed_edges=NULL, Range *corner_nodes=NULL, const bool update_meshsets=false, const bool debug=false)
Trim mesh only.
Definition: CutMeshInterface.cpp:279
MoFEM::Types::VectorDouble3
VectorBoundedArray< double, 3 > VectorDouble3
Definition: Types.hpp:92
MoFEM::CutMeshInterface::edgesToCut
map< EntityHandle, TreeData > edgesToCut
Definition: CutMeshInterface.hpp:581
QUIET
@ QUIET
Definition: definitions.h:277
MoFEM::CutMeshInterface::getCutEdges
const Range & getCutEdges() const
Definition: CutMeshInterface.hpp:483
MoFEM::CutMeshInterface::createSurfaceLevelSets
MoFEMErrorCode createSurfaceLevelSets(int verb=QUIET, const bool debug=false)
Calculate distance from mesh nodes to cut surface.
Definition: CutMeshInterface.cpp:484
MoFEM::CutMeshInterface::cutNewSurfaces
Range cutNewSurfaces
Definition: CutMeshInterface.hpp:561
MoFEM::CutMeshInterface::mergeBadEdges
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.
Definition: CutMeshInterface.cpp:2321
MoFEM::CutMeshInterface::getCutSurfaceVolumes
const Range & getCutSurfaceVolumes() const
Definition: CutMeshInterface.hpp:499
FTensor::transform
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)
Definition: Tensor2_transform.hpp:27
MoFEMFunctionReturn
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:485
MoFEM::CutMeshInterface::cutEdges
Range cutEdges
Definition: CutMeshInterface.hpp:559
MoFEM::CutMeshInterface::setTagData
MoFEMErrorCode setTagData(Tag th, const BitRefLevel bit=BitRefLevel())
set coords to tag
Definition: CutMeshInterface.cpp:3133
MoFEM::CutMeshInterface::projectZeroDistanceEnts
const Range & projectZeroDistanceEnts() const
Definition: CutMeshInterface.hpp:487
MoFEM::Exceptions::ierr
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Exceptions.hpp:87
CHKERRG
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:552
MoFEM::CutMeshInterface::SaveData
Definition: CutMeshInterface.hpp:512
MoFEM::CutMeshInterface::mergeSurface
MoFEMErrorCode mergeSurface(const Range surface)
merge surface entities
Definition: CutMeshInterface.cpp:136
MoFEM::CutMeshInterface::trimNewSurfaces
Range trimNewSurfaces
Definition: CutMeshInterface.hpp:568
MoFEMFunctionBegin
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:415
l
FTensor::Index< 'l', 3 > l
Definition: matrix_function.cpp:21
MoFEM::CutMeshInterface::SaveData::operator()
MoFEMErrorCode operator()(const std::string name, const Range &ents)
Definition: CutMeshInterface.hpp:515
MoFEM::CutMeshInterface::buildTree
MoFEMErrorCode buildTree()
build tree
Definition: CutMeshInterface.cpp:231
MoFEM::CutMeshInterface::verticesOnCutEdges
map< EntityHandle, TreeData > verticesOnCutEdges
Definition: CutMeshInterface.hpp:582
MoFEM::UnknownInterface
base class for all interface classes
Definition: UnknownInterface.hpp:84
tol
double tol
Definition: mesh_smoothing.cpp:39
MoFEM::CutMeshInterface::setConstrainSurface
MoFEMErrorCode setConstrainSurface(const Range surf)
Set the constrain surface object.
Definition: CutMeshInterface.cpp:130
MoFEM::CutMeshInterface::TreeData::lEngth
double lEngth
Definition: CutMeshInterface.hpp:576
MoFEM::CutMeshInterface::getFront
const Range & getFront() const
Definition: CutMeshInterface.hpp:481
MoFEM::CutMeshInterface::maxLength
double maxLength
Maximal edge length.
Definition: CutMeshInterface.hpp:587