Reference documentation for deal.II version 9.2.0
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
utilities.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2014 - 2020 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.md at
12 // the top level directory of deal.II.
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/base/point.h>
25 
26 # include <deal.II/fe/mapping_q1.h>
27 
28 # include <deal.II/grid/tria.h>
29 
30 # include <string>
31 
32 // opencascade needs "HAVE_CONFIG_H" to be exported...
33 # define HAVE_CONFIG_H
34 # include <IFSelect_ReturnStatus.hxx>
35 # include <TopoDS_CompSolid.hxx>
36 # include <TopoDS_Compound.hxx>
37 # include <TopoDS_Edge.hxx>
38 # include <TopoDS_Face.hxx>
39 # include <TopoDS_Shape.hxx>
40 # include <TopoDS_Shell.hxx>
41 # include <TopoDS_Solid.hxx>
42 # include <TopoDS_Vertex.hxx>
43 # include <TopoDS_Wire.hxx>
44 # include <gp_Pnt.hxx>
45 # undef HAVE_CONFIG_H
46 
47 
48 
50 
103 namespace OpenCASCADE
104 {
112  std::tuple<unsigned int, unsigned int, unsigned int>
113  count_elements(const TopoDS_Shape &shape);
114 
122  TopoDS_Shape
123  read_IGES(const std::string &filename, const double scale_factor = 1e-3);
124 
128  void
129  write_IGES(const TopoDS_Shape &shape, const std::string &filename);
130 
131 
137  TopoDS_Shape
138  read_STL(const std::string &filename);
139 
155  void
156  write_STL(const TopoDS_Shape &shape,
157  const std::string & filename,
158  const double deflection,
159  const bool sew_different_faces = false,
160  const double sewer_tolerance = 1e-6,
161  const bool is_relative = false,
162  const double angular_deflection = 0.5,
163  const bool in_parallel = false);
164 
165 
173  TopoDS_Shape
174  read_STEP(const std::string &filename, const double scale_factor = 1e-3);
175 
176 
180  void
181  write_STEP(const TopoDS_Shape &shape, const std::string &filename);
182 
195  double
196  get_shape_tolerance(const TopoDS_Shape &shape);
197 
204  TopoDS_Shape
205  intersect_plane(const TopoDS_Shape &in_shape,
206  const double c_x,
207  const double c_y,
208  const double c_z,
209  const double c,
210  const double tolerance = 1e-7);
211 
219  TopoDS_Edge
220  join_edges(const TopoDS_Shape &in_shape, const double tolerance = 1e-7);
221 
239  template <int dim>
240  TopoDS_Edge
241  interpolation_curve(std::vector<Point<dim>> &curve_points,
242  const Tensor<1, dim> & direction = Tensor<1, dim>(),
243  const bool closed = false,
244  const double tolerance = 1e-7);
245 
251  void
252  extract_geometrical_shapes(const TopoDS_Shape & shape,
253  std::vector<TopoDS_Face> & faces,
254  std::vector<TopoDS_Edge> & edges,
255  std::vector<TopoDS_Vertex> &vertices);
256 
269  template <int spacedim>
270  void
271  create_triangulation(const TopoDS_Face & face,
273 
274 
295  template <int spacedim>
296  std::vector<TopoDS_Edge>
299  const Mapping<2, spacedim> & mapping =
301 
307  void
308  extract_compound_shapes(const TopoDS_Shape & shape,
309  std::vector<TopoDS_Compound> & compounds,
310  std::vector<TopoDS_CompSolid> &compsolids,
311  std::vector<TopoDS_Solid> & solids,
312  std::vector<TopoDS_Shell> & shells,
313  std::vector<TopoDS_Wire> & wires);
314 
329  template <int dim>
330  std::tuple<Point<dim>, TopoDS_Shape, double, double>
331  project_point_and_pull_back(const TopoDS_Shape &in_shape,
332  const Point<dim> & origin,
333  const double tolerance = 1e-7);
334 
341  template <int dim>
342  Point<dim>
343  closest_point(const TopoDS_Shape &in_shape,
344  const Point<dim> & origin,
345  const double tolerance = 1e-7);
346 
355  template <int dim>
356  Point<dim>
357  push_forward(const TopoDS_Shape &in_shape, const double u, const double v);
358 
359 
365  std::tuple<Point<3>, Tensor<1, 3>, double, double>
366  push_forward_and_differential_forms(const TopoDS_Face &face,
367  const double u,
368  const double v,
369  const double tolerance = 1e-7);
370 
371 
379  std::tuple<Point<3>, Tensor<1, 3>, double, double>
380  closest_point_and_differential_forms(const TopoDS_Shape &in_shape,
381  const Point<3> & origin,
382  const double tolerance = 1e-7);
383 
384 
392  template <int dim>
393  Point<dim>
394  line_intersection(const TopoDS_Shape & in_shape,
395  const Point<dim> & origin,
396  const Tensor<1, dim> &direction,
397  const double tolerance = 1e-7);
398 
399 
407  template <int spacedim>
409  point(const gp_Pnt &p, const double tolerance = 1e-10);
410 
411 
415  template <int spacedim>
416  gp_Pnt
417  point(const Point<spacedim> &p);
418 
419 
426  template <int dim>
427  bool
428  point_compare(const Point<dim> & p1,
429  const Point<dim> & p2,
430  const Tensor<1, dim> &direction = Tensor<1, dim>(),
431  const double tolerance = 1e-10);
432 
433 
438  template <int dim>
440  Point<dim>,
441  << "The point [ " << arg1 << " ] is not on the manifold.");
442 
447  template <int dim>
449  Point<dim>,
450  << "Projection of point [ " << arg1 << " ] failed.");
451 
456  IFSelect_ReturnStatus,
457  << "An OpenCASCADE routine failed with return status "
458  << arg1);
459 
464 
469 } // namespace OpenCASCADE
470 
471 
473 
474 # endif // DEAL_II_WITH_OPENCASCADE
475 
476 #endif // dealii_occ_utilities_h
477 /*----------------------------- occ_utilities.h -----------------------------*/
OpenCASCADE::extract_compound_shapes
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:133
OpenCASCADE::push_forward_and_differential_forms
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:849
OpenCASCADE::ExcProjectionFailed
static ::ExceptionBase & ExcProjectionFailed(Point< dim > arg1)
StaticMappingQ1
Definition: mapping_q1.h:88
Triangulation
Definition: tria.h:1109
tria.h
OpenCASCADE::create_curves_from_triangulation_boundary
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:595
OpenCASCADE::ExcOCCError
static ::ExceptionBase & ExcOCCError(IFSelect_ReturnStatus arg1)
mapping_q1.h
OpenCASCADE::join_edges
TopoDS_Edge join_edges(const TopoDS_Shape &in_shape, const double tolerance=1e-7)
Definition: utilities.cc:438
Physics::Elasticity::Kinematics::e
SymmetricTensor< 2, dim, Number > e(const Tensor< 2, dim, Number > &F)
OpenCASCADE::read_STEP
TopoDS_Shape read_STEP(const std::string &filename, const double scale_factor=1e-3)
Definition: utilities.cc:355
OpenCASCADE::write_STEP
void write_STEP(const TopoDS_Shape &shape, const std::string &filename)
Definition: utilities.cc:384
OpenCASCADE::point_compare
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:215
OpenCASCADE::closest_point
Point< dim > closest_point(const TopoDS_Shape &in_shape, const Point< dim > &origin, const double tolerance=1e-7)
Definition: utilities.cc:782
OpenCASCADE::intersect_plane
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:424
OpenCASCADE::get_shape_tolerance
double get_shape_tolerance(const TopoDS_Shape &shape)
Definition: utilities.cc:400
OpenCASCADE::point
Point< spacedim > point(const gp_Pnt &p, const double tolerance=1e-10)
Definition: utilities.cc:188
OpenCASCADE::line_intersection
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:506
OpenCASCADE::read_STL
TopoDS_Shape read_STL(const std::string &filename)
Definition: utilities.cc:280
OpenCASCADE::read_IGES
TopoDS_Shape read_IGES(const std::string &filename, const double scale_factor=1e-3)
Definition: utilities.cc:238
Mapping
Abstract base class for mapping classes.
Definition: mapping.h:302
OpenCASCADE::push_forward
Point< dim > push_forward(const TopoDS_Shape &in_shape, const double u, const double v)
Definition: utilities.cc:828
double
OpenCASCADE::ExcPointNotOnManifold
static ::ExceptionBase & ExcPointNotOnManifold(Point< dim > arg1)
OpenCASCADE::interpolation_curve
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:552
Tensor< 1, dim >
DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:358
OpenCASCADE::create_triangulation
void create_triangulation(const TopoDS_Face &face, Triangulation< 2, spacedim > &tria)
Definition: utilities.cc:886
OpenCASCADE
Definition: boundary_lib.h:33
DeclException1
#define DeclException1(Exception1, type1, outsequence)
Definition: exceptions.h:518
OpenCASCADE::project_point_and_pull_back
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:699
OpenCASCADE::ExcEdgeIsDegenerate
static ::ExceptionBase & ExcEdgeIsDegenerate()
OpenCASCADE::ExcUnsupportedShape
static ::ExceptionBase & ExcUnsupportedShape()
OpenCASCADE::write_STL
void write_STL(const TopoDS_Shape &shape, const std::string &filename, const double deflection, const bool sew_different_faces=false, const double sewer_tolerance=1e-6, const bool is_relative=false, const double angular_deflection=0.5, const bool in_parallel=false)
Definition: utilities.cc:290
vertices
Point< 3 > vertices[4]
Definition: data_out_base.cc:174
OpenCASCADE::extract_geometrical_shapes
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:107
DeclException0
#define DeclException0(Exception0)
Definition: exceptions.h:473
Point< dim >
OpenCASCADE::count_elements
std::tuple< unsigned int, unsigned int, unsigned int > count_elements(const TopoDS_Shape &shape)
Definition: utilities.cc:88
config.h
OpenCASCADE::write_IGES
void write_IGES(const TopoDS_Shape &shape, const std::string &filename)
Definition: utilities.cc:267
triangulation
const typename ::parallel::distributed::Triangulation< dim, spacedim > * triangulation
Definition: p4est_wrappers.cc:69
DEAL_II_NAMESPACE_CLOSE
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:359
OpenCASCADE::closest_point_and_differential_forms
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:792
point.h