7#ifndef __CUTMESHINTERFACE_HPP__
8 #define __CUTMESHINTERFACE_HPP__
40 PetscOptionsBegin(PETSC_COMM_WORLD,
"cut_",
"MOFEM Cut mesh options",
42 CHKERR PetscOptionsInt(
"-linesearch_steps",
43 "number of bisection steps which line search do to "
44 "find optimal merged nodes position",
49 "-max_merging_cycles",
"number of maximal merging cycles",
"",
52 CHKERR PetscOptionsScalar(
"-project_entities_quality_trashold",
53 "project entities quality trashold",
"",
76 double *shift = NULL,
double *origin = NULL,
77 double *transform = NULL,
78 const std::string save_mesh =
"");
114 const double rel_tol,
115 const double abs_tol,
Tag th =
nullptr,
116 const bool debug =
false);
119 const Range fixed_edges,
120 const double rel_tol,
const double abs_tol,
144 const double tol_geometry,
const double tol_cut_close,
145 Range *fixed_edges = NULL,
Range *corner_nodes = NULL,
146 const bool update_meshsets =
false,
147 const bool debug =
false);
162 const double tol_trim_close,
163 Range *fixed_edges = NULL,
Range *corner_nodes = NULL,
164 const bool update_meshsets =
false,
165 const bool debug =
false);
183 const double tol_cut_close,
const double tol_trim_close,
184 Range *fixed_edges = NULL,
Range *corner_nodes = NULL,
185 const bool update_meshsets =
false,
const bool debug =
false);
206 Tag th,
const double tol_geometry,
207 const double tol_cut_close,
208 const double tol_trim_close,
210 const bool update_meshsets =
false,
211 const bool debug =
false);
230 const bool debug =
false);
243 const bool debug =
false);
258 const bool remove_adj_prims_edges,
int verb =
QUIET,
259 const bool debug =
false,
260 const std::string edges_file_name =
string());
271 const bool debug =
false);
289 const int front_levels,
290 Range *fixed_edges =
nullptr,
int verb =
QUIET,
291 const bool debug =
false);
302 const bool debug =
false);
316 const double geometry_tol = 0,
317 const double close_tol = 0,
318 const int verb =
QUIET,
319 const bool debug =
false);
333 const bool debug =
false);
352 Tag th = NULL,
const double tol = 1e-4,
353 const bool debug =
false);
361 const bool debug =
false);
378 const double tol = 1e-4,
Tag th = NULL,
379 const bool debug =
false);
429 const Range &corner_nodes,
Range &merged_nodes,
431 const bool update_meshsets =
false,
433 const bool debug =
false);
447 const bool update_meshsets =
false,
const bool debug =
false);
503 CHKERR moab.create_meshset(MESHSET_SET, meshset);
505 CHKERR moab.write_file(name.c_str(),
"VTK",
"", &meshset, 1);
520 typedef multi_index_container<
525 member<LengthMapData, double, &LengthMapData::lEngth>>,
528 member<LengthMapData, EntityHandle, &LengthMapData::eDge>>,
531 member<LengthMapData, double, &LengthMapData::qUality>>
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define CHKERR
Inline error check.
FTensor::Index< 'l', 3 > l
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
VectorBoundedArray< double, 3 > VectorDouble3
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
implementation of Data Operators for Forces and Sources
LengthMapData(const double l, double q, const EntityHandle e)
MoFEMErrorCode operator()(const std::string name, const Range &ents)
SaveData(moab::Interface &moab)
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 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 & getNewCutSurfaces() const
MoFEMErrorCode setFront(const Range surface)
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.
MoFEMErrorCode moveMidNodesOnTrimmedEdges(Tag th=NULL)
move trimmed edges mid nodes
const Range & getFront() const
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.
const Range & getCutEdges() const
MoFEMErrorCode getOptions()
Get options from command line.
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.
const Range & getVolume() const
double aveLength
Average edge length.
const Range & getTrimEdges() const
MoFEMErrorCode trimEdgesInTheMiddle(const BitRefLevel bit, const bool debug=false)
trim edges
const Range & getNewCutVolumes() const
MoFEMErrorCode setCoords(Tag th, const BitRefLevel bit=BitRefLevel(), const BitRefLevel mask=BitRefLevel().set())
set coords from tag
const Range & projectZeroDistanceEnts() const
const Range & getCutFrontVolumes() const
MoFEMErrorCode findEdgesToCut(Range vol, int verb=QUIET, const bool debug=false)
find edges to cut
const Range & getMergedSurfaces() const
const Range & getNewTrimVolumes() const
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
MoFEMErrorCode cutEdgesInMiddle(const BitRefLevel bit, Range &cut_vols, Range &cut_surf, Range &cut_verts, const bool debug=false)
cut edges
MoFEMErrorCode setSurface(const Range surface)
set surface entities
map< EntityHandle, TreeData > edgesToTrim
map< EntityHandle, TreeData > verticesOnTrimEdges
MoFEMErrorCode buildTree()
build tree
double projectEntitiesQualityTrashold
double maxLength
Maximal edge length.
MoFEMErrorCode splitSides(const BitRefLevel split_bit, const BitRefLevel bit, const Range &ents, Tag th=NULL)
split sides
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 mergeVolumes(const Range volume)
merge volume entities
MoFEMErrorCode getSubInterfaceOptions()
boost::shared_ptr< OrientedBoxTreeTool > & getTreeSurfPtr()
MoFEMErrorCode moveMidNodesOnCutEdges(Tag th=NULL)
projecting of mid edge nodes on new mesh on surface
map< EntityHandle, TreeData > verticesOnCutEdges
MoFEMErrorCode saveCutEdges(const std::string prefix="")
MoFEMErrorCode saveTrimEdges()
MoFEMErrorCode snapSurfaceSkinToEdges(const Range fixed_edges, const double rel_tol, const double abs_tol, Tag th=nullptr, const bool debug=false)
const Range & getNewTrimVertices() const
const Range & getSurface() const
MoFEMErrorCode clearMap()
MoFEMErrorCode removePathologicalFrontTris(const BitRefLevel split_bit, Range &ents)
Remove pathological elements on surface internal front.
MoFEMErrorCode setConstrainSurface(const Range surf)
Set the constrain surface object.
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.
const Range & getNewTrimSurfaces() const
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.
MoFEMErrorCode setVolume(const Range volume)
set volume entities
const Range & getMergedVolumes() const
MoFEMErrorCode createFrontLevelSets(Range vol, Tag th=nullptr, int verb=QUIET, const bool debug=false)
Calculate distance from mesh nodes to surface front.
MoFEMErrorCode setTagData(Tag th, const BitRefLevel bit=BitRefLevel())
set coords to tag
CutMeshInterface(const MoFEM::Core &core)
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.
map< EntityHandle, TreeData > edgesToCut
MoFEMErrorCode makeFront(const bool debug=false)
Create front from the surface.
boost::shared_ptr< OrientedBoxTreeTool > treeSurfPtr
MoFEMErrorCode mergeSurface(const Range surface)
merge surface entities
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.
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)
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
const Range & getCutSurfaceVolumes() const
MoFEMErrorCode createSurfaceLevelSets(int verb=QUIET, const bool debug=false)
Calculate distance from mesh nodes to cut surface.
base class for all interface classes