Reference documentation for deal.II version GIT relicensing-233-g802318d791 2024-03-28 20:20:02+00:00
\(\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// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 2014 - 2023 by the deal.II authors
5//
6// This file is part of the deal.II library.
7//
8// Part of the source code is dual licensed under Apache-2.0 WITH
9// LLVM-exception OR LGPL-2.1-or-later. Detailed license information
10// governing the source code and code contributions can be found in
11// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II.
12//
13// ------------------------------------------------------------------------
14
15
16#ifndef dealii_occ_utilities_h
17#define dealii_occ_utilities_h
18
19#include <deal.II/base/config.h>
20
21#ifdef DEAL_II_WITH_OPENCASCADE
22
23# include <deal.II/base/point.h>
24
26
27# include <deal.II/grid/tria.h>
28
29# include <string>
30
31// opencascade needs "HAVE_CONFIG_H" to be exported...
32# define HAVE_CONFIG_H
33# include <IFSelect_ReturnStatus.hxx>
34# include <TopoDS_CompSolid.hxx>
35# include <TopoDS_Compound.hxx>
36# include <TopoDS_Edge.hxx>
37# include <TopoDS_Face.hxx>
38# include <TopoDS_Shape.hxx>
39# include <TopoDS_Shell.hxx>
40# include <TopoDS_Solid.hxx>
41# include <TopoDS_Vertex.hxx>
42# include <TopoDS_Wire.hxx>
43# include <gp_Pnt.hxx>
44# undef HAVE_CONFIG_H
45
46
47
49
100namespace OpenCASCADE
101{
109 std::tuple<unsigned int, unsigned int, unsigned int>
110 count_elements(const TopoDS_Shape &shape);
111
119 TopoDS_Shape
120 read_IGES(const std::string &filename, const double scale_factor = 1e-3);
121
125 void
126 write_IGES(const TopoDS_Shape &shape, const std::string &filename);
127
128
134 TopoDS_Shape
135 read_STL(const std::string &filename);
136
152 void
153 write_STL(const TopoDS_Shape &shape,
154 const std::string &filename,
155 const double deflection,
156 const bool sew_different_faces = false,
157 const double sewer_tolerance = 1e-6,
158 const bool is_relative = false,
159 const double angular_deflection = 0.5,
160 const bool in_parallel = false);
161
162
170 TopoDS_Shape
171 read_STEP(const std::string &filename, const double scale_factor = 1e-3);
172
173
177 void
178 write_STEP(const TopoDS_Shape &shape, const std::string &filename);
179
192 double
193 get_shape_tolerance(const TopoDS_Shape &shape);
194
201 TopoDS_Shape
202 intersect_plane(const TopoDS_Shape &in_shape,
203 const double c_x,
204 const double c_y,
205 const double c_z,
206 const double c,
207 const double tolerance = 1e-7);
208
216 TopoDS_Edge
217 join_edges(const TopoDS_Shape &in_shape, const double tolerance = 1e-7);
218
236 template <int dim>
237 TopoDS_Edge
238 interpolation_curve(std::vector<Point<dim>> &curve_points,
239 const Tensor<1, dim> &direction = Tensor<1, dim>(),
240 const bool closed = false,
241 const double tolerance = 1e-7);
242
248 void
249 extract_geometrical_shapes(const TopoDS_Shape &shape,
250 std::vector<TopoDS_Face> &faces,
251 std::vector<TopoDS_Edge> &edges,
252 std::vector<TopoDS_Vertex> &vertices);
253
266 template <int spacedim>
267 void
268 create_triangulation(const TopoDS_Face &face,
270
271
290 template <int spacedim>
291 std::vector<TopoDS_Edge>
294 const Mapping<2, spacedim> &mapping =
296
302 void
303 extract_compound_shapes(const TopoDS_Shape &shape,
304 std::vector<TopoDS_Compound> &compounds,
305 std::vector<TopoDS_CompSolid> &compsolids,
306 std::vector<TopoDS_Solid> &solids,
307 std::vector<TopoDS_Shell> &shells,
308 std::vector<TopoDS_Wire> &wires);
309
324 template <int dim>
325 std::tuple<Point<dim>, TopoDS_Shape, double, double>
326 project_point_and_pull_back(const TopoDS_Shape &in_shape,
327 const Point<dim> &origin,
328 const double tolerance = 1e-7);
329
336 template <int dim>
338 closest_point(const TopoDS_Shape &in_shape,
339 const Point<dim> &origin,
340 const double tolerance = 1e-7);
341
350 template <int dim>
352 push_forward(const TopoDS_Shape &in_shape, const double u, const double v);
353
354
360 std::tuple<Point<3>, Tensor<1, 3>, double, double>
361 push_forward_and_differential_forms(const TopoDS_Face &face,
362 const double u,
363 const double v,
364 const double tolerance = 1e-7);
365
366
374 std::tuple<Point<3>, Tensor<1, 3>, double, double>
375 closest_point_and_differential_forms(const TopoDS_Shape &in_shape,
376 const Point<3> &origin,
377 const double tolerance = 1e-7);
378
379
387 template <int dim>
389 line_intersection(const TopoDS_Shape &in_shape,
390 const Point<dim> &origin,
391 const Tensor<1, dim> &direction,
392 const double tolerance = 1e-7);
393
394
402 template <int spacedim>
404 point(const gp_Pnt &p, const double tolerance = 1e-10);
405
406
410 template <int spacedim>
411 gp_Pnt
412 point(const Point<spacedim> &p);
413
414
421 template <int dim>
422 bool
423 point_compare(const Point<dim> &p1,
424 const Point<dim> &p2,
425 const Tensor<1, dim> &direction = Tensor<1, dim>(),
426 const double tolerance = 1e-10);
427
428
433 template <int dim>
436 << "The point [ " << arg1 << " ] is not on the manifold.");
437
442 template <int dim>
445 << "Projection of point [ " << arg1 << " ] failed.");
446
451 IFSelect_ReturnStatus,
452 << "An OpenCASCADE routine failed with return status "
453 << arg1);
454
459
464} // namespace OpenCASCADE
465
466
468
469#endif // DEAL_II_WITH_OPENCASCADE
470
471#endif
Abstract base class for mapping classes.
Definition mapping.h:316
Definition point.h:111
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:502
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:503
Point< 3 > vertices[4]
#define DeclException0(Exception0)
Definition exceptions.h:471
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:516
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:559
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:855
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:136
Point< dim > push_forward(const TopoDS_Shape &in_shape, const double u, const double v)
Definition utilities.cc:834
std::tuple< unsigned int, unsigned int, unsigned int > count_elements(const TopoDS_Shape &shape)
Definition utilities.cc:91
void write_STEP(const TopoDS_Shape &shape, const std::string &filename)
Definition utilities.cc:387
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:798
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:110
TopoDS_Edge join_edges(const TopoDS_Shape &in_shape, const double tolerance=1e-7)
Definition utilities.cc:445
TopoDS_Shape read_STEP(const std::string &filename, const double scale_factor=1e-3)
Definition utilities.cc:358
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:218
Point< spacedim > point(const gp_Pnt &p, const double tolerance=1e-10)
Definition utilities.cc:191
Point< dim > closest_point(const TopoDS_Shape &in_shape, const Point< dim > &origin, const double tolerance=1e-7)
Definition utilities.cc:788
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:705
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:513
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:293
void create_triangulation(const TopoDS_Face &face, Triangulation< 2, spacedim > &tria)
Definition utilities.cc:892
TopoDS_Shape read_STL(const std::string &filename)
Definition utilities.cc:283
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:427
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:602
double get_shape_tolerance(const TopoDS_Shape &shape)
Definition utilities.cc:403
void write_IGES(const TopoDS_Shape &shape, const std::string &filename)
Definition utilities.cc:270
TopoDS_Shape read_IGES(const std::string &filename, const double scale_factor=1e-3)
Definition utilities.cc:241
const ::parallel::distributed::Triangulation< dim, spacedim > * triangulation