Reference documentation for deal.II version 9.5.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\}}\)
Loading...
Searching...
No Matches
utilities.h
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2014 - 2022 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
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
101namespace OpenCASCADE
102{
110 std::tuple<unsigned int, unsigned int, unsigned int>
111 count_elements(const TopoDS_Shape &shape);
112
120 TopoDS_Shape
121 read_IGES(const std::string &filename, const double scale_factor = 1e-3);
122
126 void
127 write_IGES(const TopoDS_Shape &shape, const std::string &filename);
128
129
135 TopoDS_Shape
136 read_STL(const std::string &filename);
137
153 void
154 write_STL(const TopoDS_Shape &shape,
155 const std::string & filename,
156 const double deflection,
157 const bool sew_different_faces = false,
158 const double sewer_tolerance = 1e-6,
159 const bool is_relative = false,
160 const double angular_deflection = 0.5,
161 const bool in_parallel = false);
162
163
171 TopoDS_Shape
172 read_STEP(const std::string &filename, const double scale_factor = 1e-3);
173
174
178 void
179 write_STEP(const TopoDS_Shape &shape, const std::string &filename);
180
193 double
194 get_shape_tolerance(const TopoDS_Shape &shape);
195
202 TopoDS_Shape
203 intersect_plane(const TopoDS_Shape &in_shape,
204 const double c_x,
205 const double c_y,
206 const double c_z,
207 const double c,
208 const double tolerance = 1e-7);
209
217 TopoDS_Edge
218 join_edges(const TopoDS_Shape &in_shape, const double tolerance = 1e-7);
219
237 template <int dim>
238 TopoDS_Edge
239 interpolation_curve(std::vector<Point<dim>> &curve_points,
240 const Tensor<1, dim> & direction = Tensor<1, dim>(),
241 const bool closed = false,
242 const double tolerance = 1e-7);
243
249 void
250 extract_geometrical_shapes(const TopoDS_Shape & shape,
251 std::vector<TopoDS_Face> & faces,
252 std::vector<TopoDS_Edge> & edges,
253 std::vector<TopoDS_Vertex> &vertices);
254
267 template <int spacedim>
268 void
269 create_triangulation(const TopoDS_Face & face,
271
272
291 template <int spacedim>
292 std::vector<TopoDS_Edge>
295 const Mapping<2, spacedim> & mapping =
297
303 void
304 extract_compound_shapes(const TopoDS_Shape & shape,
305 std::vector<TopoDS_Compound> & compounds,
306 std::vector<TopoDS_CompSolid> &compsolids,
307 std::vector<TopoDS_Solid> & solids,
308 std::vector<TopoDS_Shell> & shells,
309 std::vector<TopoDS_Wire> & wires);
310
325 template <int dim>
326 std::tuple<Point<dim>, TopoDS_Shape, double, double>
327 project_point_and_pull_back(const TopoDS_Shape &in_shape,
328 const Point<dim> & origin,
329 const double tolerance = 1e-7);
330
337 template <int dim>
339 closest_point(const TopoDS_Shape &in_shape,
340 const Point<dim> & origin,
341 const double tolerance = 1e-7);
342
351 template <int dim>
353 push_forward(const TopoDS_Shape &in_shape, const double u, const double v);
354
355
361 std::tuple<Point<3>, Tensor<1, 3>, double, double>
362 push_forward_and_differential_forms(const TopoDS_Face &face,
363 const double u,
364 const double v,
365 const double tolerance = 1e-7);
366
367
375 std::tuple<Point<3>, Tensor<1, 3>, double, double>
376 closest_point_and_differential_forms(const TopoDS_Shape &in_shape,
377 const Point<3> & origin,
378 const double tolerance = 1e-7);
379
380
388 template <int dim>
390 line_intersection(const TopoDS_Shape & in_shape,
391 const Point<dim> & origin,
392 const Tensor<1, dim> &direction,
393 const double tolerance = 1e-7);
394
395
403 template <int spacedim>
405 point(const gp_Pnt &p, const double tolerance = 1e-10);
406
407
411 template <int spacedim>
412 gp_Pnt
413 point(const Point<spacedim> &p);
414
415
422 template <int dim>
423 bool
424 point_compare(const Point<dim> & p1,
425 const Point<dim> & p2,
426 const Tensor<1, dim> &direction = Tensor<1, dim>(),
427 const double tolerance = 1e-10);
428
429
434 template <int dim>
437 << "The point [ " << arg1 << " ] is not on the manifold.");
438
443 template <int dim>
446 << "Projection of point [ " << arg1 << " ] failed.");
447
452 IFSelect_ReturnStatus,
453 << "An OpenCASCADE routine failed with return status "
454 << arg1);
455
460
465} // namespace OpenCASCADE
466
467
469
470#endif // DEAL_II_WITH_OPENCASCADE
471
472#endif
Abstract base class for mapping classes.
Definition mapping.h:317
Definition point.h:112
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:472
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
Point< 3 > vertices[4]
#define DeclException0(Exception0)
Definition exceptions.h:465
static ::ExceptionBase & ExcPointNotOnManifold(Point< dim > arg1)
static ::ExceptionBase & ExcEdgeIsDegenerate()
static ::ExceptionBase & ExcProjectionFailed(Point< dim > arg1)
static ::ExceptionBase & ExcOCCError(IFSelect_ReturnStatus arg1)
static ::ExceptionBase & ExcUnsupportedShape()
#define DeclException1(Exception1, type1, outsequence)
Definition exceptions.h:510
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:557
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:858
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:134
Point< dim > push_forward(const TopoDS_Shape &in_shape, const double u, const double v)
Definition utilities.cc:837
std::tuple< unsigned int, unsigned int, unsigned int > count_elements(const TopoDS_Shape &shape)
Definition utilities.cc:92
void write_STEP(const TopoDS_Shape &shape, const std::string &filename)
Definition utilities.cc:385
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:801
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:108
TopoDS_Edge join_edges(const TopoDS_Shape &in_shape, const double tolerance=1e-7)
Definition utilities.cc:443
TopoDS_Shape read_STEP(const std::string &filename, const double scale_factor=1e-3)
Definition utilities.cc:356
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:216
Point< spacedim > point(const gp_Pnt &p, const double tolerance=1e-10)
Definition utilities.cc:189
Point< dim > closest_point(const TopoDS_Shape &in_shape, const Point< dim > &origin, const double tolerance=1e-7)
Definition utilities.cc:791
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:703
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:511
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:291
void create_triangulation(const TopoDS_Face &face, Triangulation< 2, spacedim > &tria)
Definition utilities.cc:895
TopoDS_Shape read_STL(const std::string &filename)
Definition utilities.cc:281
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:425
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:600
double get_shape_tolerance(const TopoDS_Shape &shape)
Definition utilities.cc:401
void write_IGES(const TopoDS_Shape &shape, const std::string &filename)
Definition utilities.cc:268
TopoDS_Shape read_IGES(const std::string &filename, const double scale_factor=1e-3)
Definition utilities.cc:239
const ::parallel::distributed::Triangulation< dim, spacedim > * triangulation
const ::Triangulation< dim, spacedim > & tria