Reference documentation for deal.II version 9.0.0
utilities.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2014 - 2018 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE at
12 // the top level of the deal.II distribution.
13 //
14 // ---------------------------------------------------------------------
15 
16 
17 #ifndef dealii_occ_utilities_h
18 #define dealii_occ_utilities_h
19 
20 #include <deal.II/base/config.h>
21 
22 #ifdef DEAL_II_WITH_OPENCASCADE
23 
24 #include <deal.II/grid/tria.h>
25 #include <deal.II/base/point.h>
26 #include <deal.II/fe/mapping_q1.h>
27 
28 #include <string>
29 
30 // opencascade needs "HAVE_CONFIG_H" to be exported...
31 #define HAVE_CONFIG_H
32 #include <TopoDS_Shape.hxx>
33 #include <TopoDS_Face.hxx>
34 #include <TopoDS_Edge.hxx>
35 #include <TopoDS_Vertex.hxx>
36 #include <TopoDS_Compound.hxx>
37 #include <TopoDS_CompSolid.hxx>
38 #include <TopoDS_Solid.hxx>
39 #include <TopoDS_Shell.hxx>
40 #include <TopoDS_Wire.hxx>
41 #include <IFSelect_ReturnStatus.hxx>
42 #include <gp_Pnt.hxx>
43 #undef HAVE_CONFIG_H
44 
45 
46 
47 DEAL_II_NAMESPACE_OPEN
48 
101 namespace OpenCASCADE
102 {
110  std::tuple<unsigned int, unsigned int, unsigned int>
111  count_elements(const TopoDS_Shape &shape);
112 
120  TopoDS_Shape read_IGES(const std::string &filename,
121  const double scale_factor=1e-3);
122 
126  void write_IGES(const TopoDS_Shape &shape,
127  const std::string &filename);
128 
136  TopoDS_Shape read_STEP(const std::string &filename,
137  const double scale_factor=1e-3);
138 
139 
143  void write_STEP(const TopoDS_Shape &shape,
144  const std::string &filename);
145 
158  double get_shape_tolerance(const TopoDS_Shape &shape);
159 
166  TopoDS_Shape intersect_plane(const TopoDS_Shape &in_shape,
167  const double c_x,
168  const double c_y,
169  const double c_z,
170  const double c,
171  const double tolerance=1e-7);
172 
180  TopoDS_Edge join_edges(const TopoDS_Shape &in_shape,
181  const double tolerance=1e-7);
182 
200  template<int dim>
201  TopoDS_Edge interpolation_curve(std::vector<Point<dim> > &curve_points,
202  const Tensor<1, dim> &direction=Tensor<1,dim>(),
203  const bool closed=false,
204  const double tolerance=1e-7);
205 
211  void extract_geometrical_shapes(const TopoDS_Shape &shape,
212  std::vector<TopoDS_Face> &faces,
213  std::vector<TopoDS_Edge> &edges,
214  std::vector<TopoDS_Vertex> &vertices);
215 
228  template <int spacedim>
229  void create_triangulation(const TopoDS_Face &face,
231 
232 
253  template <int spacedim>
254  std::vector<TopoDS_Edge> create_curves_from_triangulation_boundary
255  (const Triangulation<2,spacedim> &triangulation,
257 
263  void extract_compound_shapes(const TopoDS_Shape &shape,
264  std::vector<TopoDS_Compound> &compounds,
265  std::vector<TopoDS_CompSolid> &compsolids,
266  std::vector<TopoDS_Solid> &solids,
267  std::vector<TopoDS_Shell> &shells,
268  std::vector<TopoDS_Wire> &wires);
269 
284  template <int dim>
285  std::tuple<Point<dim>, TopoDS_Shape, double, double>
286  project_point_and_pull_back(const TopoDS_Shape &in_shape,
287  const Point<dim> &origin,
288  const double tolerance=1e-7);
289 
296  template <int dim>
297  Point<dim> closest_point(const TopoDS_Shape &in_shape,
298  const Point<dim> &origin,
299  const double tolerance=1e-7);
300 
309  template <int dim>
310  Point<dim> push_forward(const TopoDS_Shape &in_shape,
311  const double u,
312  const double v);
313 
314 
320  std::tuple<Point<3>, Tensor<1,3>, double, double>
321  push_forward_and_differential_forms(const TopoDS_Face &face,
322  const double u,
323  const double v,
324  const double tolerance=1e-7);
325 
326 
334  std::tuple<Point<3>, Tensor<1,3>, double, double>
335  closest_point_and_differential_forms(const TopoDS_Shape &in_shape,
336  const Point<3> &origin,
337  const double tolerance=1e-7);
338 
339 
347  template <int dim>
348  Point<dim> line_intersection(const TopoDS_Shape &in_shape,
349  const Point<dim> &origin,
350  const Tensor<1,dim> &direction,
351  const double tolerance=1e-7);
352 
353 
361  template <int spacedim>
362  Point<spacedim> point(const gp_Pnt &p, const double &tolerance=1e-10);
363 
364 
368  template <int spacedim>
369  gp_Pnt point(const Point<spacedim> &p);
370 
371 
378  template <int dim>
379  bool point_compare(const Point<dim> &p1,
380  const Point<dim> &p2,
381  const Tensor<1,dim> &direction = Tensor<1,dim>(),
382  const double tolerance = 1e-10);
383 
384 
389  template <int dim>
391  Point<dim>,
392  <<"The point [ "<<arg1<<" ] is not on the manifold.");
393 
398  template <int dim>
400  Point<dim>,
401  <<"Projection of point [ "<< arg1
402  << " ] failed.");
403 
408  IFSelect_ReturnStatus,
409  <<"An OpenCASCADE routine failed with return status "
410  <<arg1);
411 
416 
421 }
422 
423 
424 DEAL_II_NAMESPACE_CLOSE
425 
426 #endif // DEAL_II_WITH_OPENCASCADE
427 
428 /*------------------------------ occ_utilities.h ------------------------------*/
429 #endif
430 /*------------------------------ occ_utilities.h ------------------------------*/
std::tuple< unsigned int, unsigned int, unsigned int > count_elements(const TopoDS_Shape &shape)
Definition: utilities.cc:89
TopoDS_Shape read_IGES(const std::string &filename, const double scale_factor=1e-3)
Definition: utilities.cc:227
Point< dim > push_forward(const TopoDS_Shape &in_shape, const double u, const double v)
Definition: utilities.cc:719
TopoDS_Shape intersect_plane(const TopoDS_Shape &in_shape, const double c_x, const double c_y, const double c_z, const double c, const double tolerance=1e-7)
Definition: utilities.cc:341
Point< spacedim > point(const gp_Pnt &p, const double &tolerance=1e-10)
Definition: utilities.cc:183
TopoDS_Shape read_STEP(const std::string &filename, const double scale_factor=1e-3)
Definition: utilities.cc:270
TopoDS_Edge interpolation_curve(std::vector< Point< dim > > &curve_points, const Tensor< 1, dim > &direction=Tensor< 1, dim >(), const bool closed=false, const double tolerance=1e-7)
Definition: utilities.cc:459
Point< dim > line_intersection(const TopoDS_Shape &in_shape, const Point< dim > &origin, const Tensor< 1, dim > &direction, const double tolerance=1e-7)
Definition: utilities.cc:418
static ::ExceptionBase & ExcOCCError(IFSelect_ReturnStatus arg1)
bool point_compare(const Point< dim > &p1, const Point< dim > &p2, const Tensor< 1, dim > &direction=Tensor< 1, dim >(), const double tolerance=1e-10)
Definition: utilities.cc:206
std::tuple< Point< 3 >, Tensor< 1, 3 >, double, double > push_forward_and_differential_forms(const TopoDS_Face &face, const double u, const double v, const double tolerance=1e-7)
Definition: utilities.cc:742
static ::ExceptionBase & ExcEdgeIsDegenerate()
void write_IGES(const TopoDS_Shape &shape, const std::string &filename)
Definition: utilities.cc:258
void create_triangulation(const TopoDS_Face &face, Triangulation< 2, spacedim > &tria)
Definition: utilities.cc:774
#define DeclException1(Exception1, type1, outsequence)
Definition: exceptions.h:346
void write_STEP(const TopoDS_Shape &shape, const std::string &filename)
Definition: utilities.cc:301
Abstract base class for mapping classes.
Definition: dof_tools.h:46
#define DeclException0(Exception0)
Definition: exceptions.h:323
Point< dim > closest_point(const TopoDS_Shape &in_shape, const Point< dim > &origin, const double tolerance=1e-7)
Definition: utilities.cc:676
void extract_compound_shapes(const TopoDS_Shape &shape, std::vector< TopoDS_Compound > &compounds, std::vector< TopoDS_CompSolid > &compsolids, std::vector< TopoDS_Solid > &solids, std::vector< TopoDS_Shell > &shells, std::vector< TopoDS_Wire > &wires)
Definition: utilities.cc:130
std::tuple< Point< dim >, TopoDS_Shape, double, double > project_point_and_pull_back(const TopoDS_Shape &in_shape, const Point< dim > &origin, const double tolerance=1e-7)
Definition: utilities.cc:595
double get_shape_tolerance(const TopoDS_Shape &shape)
Definition: utilities.cc:315
TopoDS_Edge join_edges(const TopoDS_Shape &in_shape, const double tolerance=1e-7)
Definition: utilities.cc:354
void extract_geometrical_shapes(const TopoDS_Shape &shape, std::vector< TopoDS_Face > &faces, std::vector< TopoDS_Edge > &edges, std::vector< TopoDS_Vertex > &vertices)
Definition: utilities.cc:105
std::vector< TopoDS_Edge > create_curves_from_triangulation_boundary(const Triangulation< 2, spacedim > &triangulation, const Mapping< 2, spacedim > &mapping=StaticMappingQ1< 2, spacedim >::mapping)
Definition: utilities.cc:497
static ::ExceptionBase & ExcPointNotOnManifold(Point< dim > arg1)
static ::ExceptionBase & ExcProjectionFailed(Point< dim > arg1)
static ::ExceptionBase & ExcUnsupportedShape()
std::tuple< Point< 3 >, Tensor< 1, 3 >, double, double > closest_point_and_differential_forms(const TopoDS_Shape &in_shape, const Point< 3 > &origin, const double tolerance=1e-7)
Definition: utilities.cc:686