v0.14.0
Loading...
Searching...
No Matches
TetGenInterface.hpp
Go to the documentation of this file.
1/** \file TetGenInterface.hpp
2 * \brief TetGen interface
3 *
4 * MoFEM TetGen interface
5 *
6 *
7 * \ingroup mesh_tetgen
8 */
9
10#ifndef __TETGENINTERFACE_HPP__
11#define __TETGENINTERFACE_HPP__
12
13#include "UnknownInterface.hpp"
14
15class tetgenio;
16
17namespace MoFEM {
18
19/** \brief TetGen interface
20
21 * \ingroup mesh_tetgen
22 */
24
25 MoFEMErrorCode query_interface(boost::typeindex::type_index type_index,
26 UnknownInterface **iface) const;
27
30 : cOre(const_cast<MoFEM::Core &>(core)) {}
31
32 typedef std::map<EntityHandle, unsigned long> moabTetGen_Map;
33 typedef std::map<unsigned long, EntityHandle> tetGenMoab_Map;
34 typedef std::map<int, Range> idxRange_Map;
35
36 /** \brief create TetGen data structure form range of moab entities
37 *
38 * Move mesh to TetGen data structures
39 *
40 * \param ents range of entities (tetrahedrons or nodes)
41 * \param in tetgen data structure (look to TetGen user manual)
42 * \param moab_tetgen_map mapping moab to TetGen entities
43 * \param tetgen_moab_map mapping tetgen to moab entities
44 */
45 MoFEMErrorCode inData(Range &ents, tetgenio &in,
46 moabTetGen_Map &moab_tetgen_map,
47 tetGenMoab_Map &tetgen_moab_map, Tag th = NULL);
48
54 };
55
56 /** \brief set point tags and type
57
58 Set type of entity, look in TetGen manual for details
59
60 \code
61 std::map<int,Range> types_ents;
62 //RIDGEVERTEX
63 types_ents[TetGenInterface::RIDGEVERTEX].merge(region_tets_skin_without_boundary_nodes);
64 //FREESEGVERTEX
65 types_ents[TetGenInterface::FREESEGVERTEX].merge(crack_surface_tris_skin_nodes);
66 //FREEFACETVERTEX
67 types_ents[TetGenInterface::FREEFACETVERTEX].merge(region_tets_skin_nodes);
68 types_ents[TetGenInterface::FREEFACETVERTEX] =
69 subtract(types_ents[TetGenInterface::FREEFACETVERTEX],types_ents[TetGenInterface::RIDGEVERTEX]);
70 types_ents[TetGenInterface::FREEFACETVERTEX] =
71 subtract(types_ents[TetGenInterface::FREEFACETVERTEX],types_ents[TetGenInterface::FREESEGVERTEX]);
72 //FREEVOLVERTEX
73 types_ents[TetGenInterface::FREEVOLVERTEX].merge(region_nodes);
74 types_ents[TetGenInterface::FREEVOLVERTEX] =
75 subtract(types_ents[TetGenInterface::FREEVOLVERTEX],types_ents[TetGenInterface::RIDGEVERTEX]);
76 types_ents[TetGenInterface::FREEVOLVERTEX] =
77 subtract(types_ents[TetGenInterface::FREEVOLVERTEX],types_ents[TetGenInterface::FREESEGVERTEX]);
78 types_ents[TetGenInterface::FREEVOLVERTEX] =
79 subtract(types_ents[TetGenInterface::FREEVOLVERTEX],types_ents[TetGenInterface::FREEFACETVERTEX]);
80 \endcode
81
82 */
83 MoFEMErrorCode setGeomData(tetgenio &in, moabTetGen_Map &moab_tetgen_map,
84 tetGenMoab_Map &tetgen_moab_map,
85 std::map<int, Range> &type_ents);
86
87 /** \brief get entities for TetGen data structure
88
89 \param ents range of entities (tetrahedrons or nodes)
90 \param in tetgen data structure (look to TetGen user manual)
91 \param moab_tetgen_map mapping MoAB to TetGen entities
92 \param tetgen_moab_map mapping TetGen to moab entities
93 \param ents rerun entities which are in TetGen Data structure
94 \param id_in_tags use tags as entity handles, if that is a case use tag to
95 find MoAB vertex id
96 \param error_if_created throw error if node need to be created
97
98 */
99 MoFEMErrorCode outData(tetgenio &in, tetgenio &out,
100 moabTetGen_Map &moab_tetgen_map,
101 tetGenMoab_Map &tetgen_moab_map, Range *ents = NULL,
102 bool id_in_tags = false, bool error_if_created = false,
103 bool assume_first_nodes_the_same = false,
104 Tag th = nullptr);
105
106 /** \brief get entities for TetGen data structure
107
108 \param ents range of entities (tetrahedrons or nodes)
109 \param in tetgen data structure (look to TetGen user manual)
110 \param moab_tetgen_map mapping MoAB to TetGen entities
111 \param tetgen_moab_map mapping TetGen to MoAB entities
112 \param ents rerun entities which are in TetGen data structure
113 \param bit set level to created entities
114 \param error_if_created throw error if node need to be created
115
116 */
117 MoFEMErrorCode outData(tetgenio &in, tetgenio &out,
118 moabTetGen_Map &moab_tetgen_map,
119 tetGenMoab_Map &tetgen_moab_map, BitRefLevel bit,
120 bool id_in_tags = false, bool error_if_created = false,
121 bool assume_first_nodes_the_same = false,
122 Tag th = nullptr);
123
124 /** \brief set markers to faces
125
126 \param markers data structure with markers
127 \param in tetgen data structure (look to TetGen user manual)
128 \param moab_tetgen_map mapping MoAB to TetGen entities
129 \param tetgen_moab_map mapping TetGen to MoAB entities
130
131 */
132 MoFEMErrorCode setFaceData(std::vector<std::pair<Range, int> > &markers,
133 tetgenio &in, moabTetGen_Map &moab_tetgen_map,
134 tetGenMoab_Map &tetgen_moab_map);
135
136 /** \brief get markers to faces
137
138 \param markers data structure with markers
139 \param in tetgen data structure (look to TetGen user manual)
140 \param moab_tetgen_map mapping MoAB to TetGen entities
141 \param tetgen_moab_map mapping TetGen to MoAB entities
142
143 */
145 tetgenio &out, Range *ents = NULL,
146 idxRange_Map *ents_map = NULL,
147 bool only_non_zero = true);
148
149 /** \brief set region data to tetrahedral
150 */
152 setRegionData(std::vector<std::pair<EntityHandle, int> > &regions, tetgenio &in,
153 Tag th = NULL);
154
155 /** \brief get region data to tetrahedral
156 */
157 MoFEMErrorCode getRegionData(tetGenMoab_Map &tetgen_moab_map, tetgenio &out,
158 Range *ents = NULL,
159 idxRange_Map *ents_map = NULL);
160
161 /** \brief run tetgen
162 */
163 MoFEMErrorCode tetRahedralize(char switches[], tetgenio &in, tetgenio &out);
164
165 /** \brief load poly file
166 */
167 MoFEMErrorCode loadPoly(char file_name[], tetgenio &in);
168
169 // Tools for TetGen, i.e. geometry reconstruction from mesh
170
171 MoFEMErrorCode checkPlanar_Trinagle(double coords[], bool *result,
172 const double eps = 1e-9);
173 MoFEMErrorCode groupPlanar_Triangle(Range &tris, std::vector<Range> &sorted,
174 const double eps = 1e-9, Tag th = NULL);
175
176 /** \brief Group surface triangles in planar regions
177
178 \param tris input triangles
179 \param sorted output sorted planar faces
180 \param eps tolerance
181
182 */
184 std::vector<std::vector<Range> > &sorted,
185 const double eps = 1e-9);
186
187 /** make planar polygon facet
188
189 \param ents surface triangles
190 \param polygons output list of polygons
191 \param reduce_edges reduce edges if on the line
192 \param not_reducable_nodes do not reduce node on edge if in this range
193 \param eps tolerance
194
195 \bug assumes that are no holes
196
197 */
199 bool reduce_edges = false,
200 Range *not_reducable_nodes = NULL,
201 const double eps = 1e-9,Tag th = NULL);
202};
203}
204
205#endif //__TETGENINTERFACE_HPP__
206
207/**
208 * \defgroup mesh_tetgen TetGen interface
209 * \brief Interface to run TetGen
210 *
211 * \ingroup mofem
212 */
MoFEM interface.
static const double eps
auto bit
set bit
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Types.hpp:40
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:10
Core (interface) class.
Definition: Core.hpp:82
std::map< unsigned long, EntityHandle > tetGenMoab_Map
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
MoFEMErrorCode makePolygonFacet(Range &ents, Range &polygons, bool reduce_edges=false, Range *not_reducable_nodes=NULL, const double eps=1e-9, Tag th=NULL)
MoFEMErrorCode groupRegion_Triangle(Range &tris, std::vector< std::vector< Range > > &sorted, const double eps=1e-9)
Group surface triangles in planar regions.
MoFEMErrorCode setGeomData(tetgenio &in, moabTetGen_Map &moab_tetgen_map, tetGenMoab_Map &tetgen_moab_map, std::map< int, Range > &type_ents)
set point tags and type
MoFEMErrorCode inData(Range &ents, tetgenio &in, moabTetGen_Map &moab_tetgen_map, tetGenMoab_Map &tetgen_moab_map, Tag th=NULL)
create TetGen data structure form range of moab entities
MoFEMErrorCode tetRahedralize(char switches[], tetgenio &in, tetgenio &out)
run tetgen
std::map< int, Range > idxRange_Map
TetGenInterface(const MoFEM::Core &core)
MoFEMErrorCode loadPoly(char file_name[], tetgenio &in)
load poly file
MoFEMErrorCode checkPlanar_Trinagle(double coords[], bool *result, const double eps=1e-9)
MoFEMErrorCode outData(tetgenio &in, tetgenio &out, moabTetGen_Map &moab_tetgen_map, tetGenMoab_Map &tetgen_moab_map, Range *ents=NULL, bool id_in_tags=false, bool error_if_created=false, bool assume_first_nodes_the_same=false, Tag th=nullptr)
get entities for TetGen data structure
std::map< EntityHandle, unsigned long > moabTetGen_Map
MoFEMErrorCode groupPlanar_Triangle(Range &tris, std::vector< Range > &sorted, const double eps=1e-9, Tag th=NULL)
MoFEMErrorCode getRegionData(tetGenMoab_Map &tetgen_moab_map, tetgenio &out, Range *ents=NULL, idxRange_Map *ents_map=NULL)
get region data to tetrahedral
MoFEMErrorCode getTriangleMarkers(tetGenMoab_Map &tetgen_moab_map, tetgenio &out, Range *ents=NULL, idxRange_Map *ents_map=NULL, bool only_non_zero=true)
get markers to faces
MoFEMErrorCode setFaceData(std::vector< std::pair< Range, int > > &markers, tetgenio &in, moabTetGen_Map &moab_tetgen_map, tetGenMoab_Map &tetgen_moab_map)
set markers to faces
MoFEMErrorCode setRegionData(std::vector< std::pair< EntityHandle, int > > &regions, tetgenio &in, Tag th=NULL)
set region data to tetrahedral
base class for all interface classes