v0.9.0
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 /*
11  * MoFEM is distributed in the hope that it will be useful, but WITHOUT
12  * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
13  * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
14  * License for more details.
15  *
16  * You should have received a copy of the GNU Lesser General Public
17  * License along with MoFEM. If not, see <http://www.gnu.org/licenses/>
18 */
19 
20 #ifndef __TETGENINTERFACE_HPP__
21 #define __TETGENINTERFACE_HPP__
22 
23 #include "UnknownInterface.hpp"
24 
25 class tetgenio;
26 
27 namespace MoFEM {
28 
31 
32 /** \brief TetGen interface
33 
34  * \ingroup mesh_tetgen
35  */
37 
39  UnknownInterface **iface) const;
40 
43  : cOre(const_cast<MoFEM::Core &>(core)) {}
44 
45  typedef std::map<EntityHandle, unsigned long> moabTetGen_Map;
46  typedef std::map<unsigned long, EntityHandle> tetGenMoab_Map;
47  typedef std::map<int, Range> idxRange_Map;
48 
49  /** \brief create TetGen data structure form range of moab entities
50  *
51  * Move mesh to TetGen data structures
52  *
53  * \param ents range of entities (tetrahedrons or nodes)
54  * \param in tetgen data structure (look to TetGen user manual)
55  * \param moab_tetgen_map mapping moab to TetGen entities
56  * \param tetgen_moab_map mapping tetgen to moab entities
57  */
58  MoFEMErrorCode inData(Range &ents, tetgenio &in,
59  moabTetGen_Map &moab_tetgen_map,
60  tetGenMoab_Map &tetgen_moab_map, Tag th = NULL);
61 
67  };
68 
69  /** \brief set point tags and type
70 
71  Set type of entity, look in TetGen manual for details
72 
73  \code
74  std::map<int,Range> types_ents;
75  //RIDGEVERTEX
76  types_ents[TetGenInterface::RIDGEVERTEX].merge(region_tets_skin_without_boundary_nodes);
77  //FREESEGVERTEX
78  types_ents[TetGenInterface::FREESEGVERTEX].merge(crack_surface_tris_skin_nodes);
79  //FREEFACETVERTEX
80  types_ents[TetGenInterface::FREEFACETVERTEX].merge(region_tets_skin_nodes);
81  types_ents[TetGenInterface::FREEFACETVERTEX] =
82  subtract(types_ents[TetGenInterface::FREEFACETVERTEX],types_ents[TetGenInterface::RIDGEVERTEX]);
83  types_ents[TetGenInterface::FREEFACETVERTEX] =
84  subtract(types_ents[TetGenInterface::FREEFACETVERTEX],types_ents[TetGenInterface::FREESEGVERTEX]);
85  //FREEVOLVERTEX
86  types_ents[TetGenInterface::FREEVOLVERTEX].merge(region_nodes);
87  types_ents[TetGenInterface::FREEVOLVERTEX] =
88  subtract(types_ents[TetGenInterface::FREEVOLVERTEX],types_ents[TetGenInterface::RIDGEVERTEX]);
89  types_ents[TetGenInterface::FREEVOLVERTEX] =
90  subtract(types_ents[TetGenInterface::FREEVOLVERTEX],types_ents[TetGenInterface::FREESEGVERTEX]);
91  types_ents[TetGenInterface::FREEVOLVERTEX] =
92  subtract(types_ents[TetGenInterface::FREEVOLVERTEX],types_ents[TetGenInterface::FREEFACETVERTEX]);
93  \endcode
94 
95  */
96  MoFEMErrorCode setGeomData(tetgenio &in, moabTetGen_Map &moab_tetgen_map,
97  tetGenMoab_Map &tetgen_moab_map,
98  std::map<int, Range> &type_ents);
99 
100  /** \brief get entities for TetGen data structure
101 
102  \param ents range of entities (tetrahedrons or nodes)
103  \param in tetgen data structure (look to TetGen user manual)
104  \param moab_tetgen_map mapping MoAB to TetGen entities
105  \param tetgen_moab_map mapping TetGen to moab entities
106  \param ents rerun entities which are in TetGen Data structure
107  \param id_in_tags use tags as entity handles, if that is a case use tag to
108  find MoAB vertex id
109  \param error_if_created throw error if node need to be created
110 
111  */
112  MoFEMErrorCode outData(tetgenio &in, tetgenio &out,
113  moabTetGen_Map &moab_tetgen_map,
114  tetGenMoab_Map &tetgen_moab_map, Range *ents = NULL,
115  bool id_in_tags = false, bool error_if_created = false,
116  bool assume_first_nodes_the_same = false,
117  Tag th = nullptr);
118 
119  /** \brief get entities for TetGen data structure
120 
121  \param ents range of entities (tetrahedrons or nodes)
122  \param in tetgen data structure (look to TetGen user manual)
123  \param moab_tetgen_map mapping MoAB to TetGen entities
124  \param tetgen_moab_map mapping TetGen to MoAB entities
125  \param ents rerun entities which are in TetGen data structure
126  \param bit set level to created entities
127  \param error_if_created throw error if node need to be created
128 
129  */
130  MoFEMErrorCode outData(tetgenio &in, tetgenio &out,
131  moabTetGen_Map &moab_tetgen_map,
132  tetGenMoab_Map &tetgen_moab_map, BitRefLevel bit,
133  bool id_in_tags = false, bool error_if_created = false,
134  bool assume_first_nodes_the_same = false,
135  Tag th = nullptr);
136 
137  /** \brief set markers to faces
138 
139  \param markers data structure with markers
140  \param in tetgen data structure (look to TetGen user manual)
141  \param moab_tetgen_map mapping MoAB to TetGen entities
142  \param tetgen_moab_map mapping TetGen to MoAB entities
143 
144  */
145  MoFEMErrorCode setFaceData(std::vector<std::pair<Range, int> > &markers,
146  tetgenio &in, moabTetGen_Map &moab_tetgen_map,
147  tetGenMoab_Map &tetgen_moab_map);
148 
149  /** \brief get markers to faces
150 
151  \param markers data structure with markers
152  \param in tetgen data structure (look to TetGen user manual)
153  \param moab_tetgen_map mapping MoAB to TetGen entities
154  \param tetgen_moab_map mapping TetGen to MoAB entities
155 
156  */
158  tetgenio &out, Range *ents = NULL,
159  idxRange_Map *ents_map = NULL,
160  bool only_non_zero = true);
161 
162  /** \brief set region data to tetrahedral
163  */
165  setRegionData(std::vector<std::pair<EntityHandle, int> > &regions, tetgenio &in,
166  Tag th = NULL);
167 
168  /** \brief get region data to tetrahedral
169  */
170  MoFEMErrorCode getRegionData(tetGenMoab_Map &tetgen_moab_map, tetgenio &out,
171  Range *ents = NULL,
172  idxRange_Map *ents_map = NULL);
173 
174  /** \brief run tetgen
175  */
176  MoFEMErrorCode tetRahedralize(char switches[], tetgenio &in, tetgenio &out);
177 
178  /** \brief load poly file
179  */
180  MoFEMErrorCode loadPoly(char file_name[], tetgenio &in);
181 
182  // Tools for TetGen, i.e. geometry reconstruction from mesh
183 
184  MoFEMErrorCode checkPlanar_Trinagle(double coords[], bool *result,
185  const double eps = 1e-9);
186  MoFEMErrorCode groupPlanar_Triangle(Range &tris, std::vector<Range> &sorted,
187  const double eps = 1e-9, Tag th = NULL);
188 
189  /** \brief Group surface triangles in planar regions
190 
191  \param tris input triangles
192  \param sorted output sorted planar faces
193  \param eps tolerance
194 
195  */
197  std::vector<std::vector<Range> > &sorted,
198  const double eps = 1e-9);
199 
200  /** make planar polygon facet
201 
202  \param ents surface triangles
203  \param polygons output list of polygons
204  \param reduce_edges reduce edges if on the line
205  \param not_reducable_nodes do not reduce node on edge if in this range
206  \param eps tolerance
207 
208  \bug assumes that are no holes
209 
210  */
211  MoFEMErrorCode makePolygonFacet(Range &ents, Range &polygons,
212  bool reduce_edges = false,
213  Range *not_reducable_nodes = NULL,
214  const double eps = 1e-9,Tag th = NULL);
215 };
216 }
217 
218 #endif //__TETGENINTERFACE_HPP__
219 
220 /**
221  * \defgroup mesh_tetgen TetGen interface
222  * \brief Interface to run TetGen
223  *
224  * \ingroup mofem
225  */
TetGenInterface(const MoFEM::Core &core)
MoFEMErrorCode groupPlanar_Triangle(Range &tris, std::vector< Range > &sorted, const double eps=1e-9, Tag th=NULL)
MoFEMErrorCode tetRahedralize(char switches[], tetgenio &in, tetgenio &out)
run tetgen
MoFEMErrorCode getRegionData(tetGenMoab_Map &tetgen_moab_map, tetgenio &out, Range *ents=NULL, idxRange_Map *ents_map=NULL)
get region data to tetrahedral
MoFEM interface unique ID.
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
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 makePolygonFacet(Range &ents, Range &polygons, bool reduce_edges=false, Range *not_reducable_nodes=NULL, const double eps=1e-9, Tag th=NULL)
base class for all interface classes
MoFEMErrorCode groupRegion_Triangle(Range &tris, std::vector< std::vector< Range > > &sorted, const double eps=1e-9)
Group surface triangles in planar regions.
Core (interface) class.
Definition: Core.hpp:50
MoFEM interface.
std::bitset< BITINTERFACEUID_SIZE > BitIntefaceId
Definition: Types.hpp:56
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
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
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:66
std::map< EntityHandle, unsigned long > moabTetGen_Map
std::map< int, Range > idxRange_Map
std::map< unsigned long, EntityHandle > tetGenMoab_Map
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
static const MOFEMuuid IDD_MOFEMTetGegInterface
MoFEMErrorCode loadPoly(char file_name[], tetgenio &in)
load poly file
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Types.hpp:51
MoFEMErrorCode checkPlanar_Trinagle(double coords[], bool *result, const double eps=1e-9)
static const double eps
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 query_interface(const MOFEMuuid &uuid, UnknownInterface **iface) const