Reference documentation for deal.II version 9.1.1
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
manifold_lib.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2014 - 2019 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_manifold_lib_h
18 #define dealii_occ_manifold_lib_h
19 
20 #include <deal.II/base/config.h>
21 
22 #ifdef DEAL_II_WITH_OPENCASCADE
23 
24 # include <deal.II/grid/manifold.h>
25 
26 # include <deal.II/opencascade/utilities.h>
27 
28 // opencascade needs "HAVE_CONFIG_H" to be exported...
29 # define HAVE_CONFIG_H
30 # include <Adaptor3d_Curve.hxx>
31 # include <Adaptor3d_HCurve.hxx>
32 # include <BRepAdaptor_Curve.hxx>
33 # undef HAVE_CONFIG_H
34 
35 DEAL_II_NAMESPACE_OPEN
36 
42 namespace OpenCASCADE
43 {
65  template <int dim, int spacedim>
66  class NormalProjectionManifold : public FlatManifold<dim, spacedim>
67  {
68  public:
76  NormalProjectionManifold(const TopoDS_Shape &sh,
77  const double tolerance = 1e-7);
78 
82  virtual std::unique_ptr<Manifold<dim, spacedim>>
83  clone() const override;
84 
94  virtual Point<spacedim>
96  const ArrayView<const Point<spacedim>> &surrounding_points,
97  const Point<spacedim> & candidate) const override;
98 
99 
100  protected:
107  const TopoDS_Shape sh;
108 
112  const double tolerance;
113  };
114 
136  template <int dim, int spacedim>
137  class DirectionalProjectionManifold : public FlatManifold<dim, spacedim>
138  {
139  public:
144  DirectionalProjectionManifold(const TopoDS_Shape & sh,
146  const double tolerance = 1e-7);
147 
151  virtual std::unique_ptr<Manifold<dim, spacedim>>
152  clone() const override;
153 
163  virtual Point<spacedim>
165  const ArrayView<const Point<spacedim>> &surrounding_points,
166  const Point<spacedim> & candidate) const override;
167 
168  protected:
175  const TopoDS_Shape sh;
176 
181 
185  const double tolerance;
186  };
187 
188 
233  template <int dim, int spacedim>
234  class NormalToMeshProjectionManifold : public FlatManifold<dim, spacedim>
235  {
236  public:
242  NormalToMeshProjectionManifold(const TopoDS_Shape &sh,
243  const double tolerance = 1e-7);
244 
248  virtual std::unique_ptr<Manifold<dim, spacedim>>
249  clone() const override;
250 
257  virtual Point<spacedim>
259  const ArrayView<const Point<spacedim>> &surrounding_points,
260  const Point<spacedim> & candidate) const override;
261 
262  protected:
269  const TopoDS_Shape sh;
270 
274  const double tolerance;
275  };
276 
296  template <int dim, int spacedim>
297  class ArclengthProjectionLineManifold : public ChartManifold<dim, spacedim, 1>
298  {
299  public:
303  ArclengthProjectionLineManifold(const TopoDS_Shape &sh,
304  const double tolerance = 1e-7);
305 
309  virtual std::unique_ptr<Manifold<dim, spacedim>>
310  clone() const override;
311 
317  virtual Point<1>
318  pull_back(const Point<spacedim> &space_point) const override;
319 
323  virtual Point<spacedim>
324  push_forward(const Point<1> &chart_point) const override;
325 
326  protected:
330  const TopoDS_Shape sh;
331 
336  Handle_Adaptor3d_HCurve curve;
337 
341  const double tolerance;
342 
347  const double length;
348  };
349 
357  template <int dim, int spacedim>
358  class NURBSPatchManifold : public ChartManifold<dim, spacedim, 2>
359  {
360  public:
366  NURBSPatchManifold(const TopoDS_Face &face, const double tolerance = 1e-7);
367 
371  virtual std::unique_ptr<Manifold<dim, spacedim>>
372  clone() const override;
373 
378  virtual Point<2>
379  pull_back(const Point<spacedim> &space_point) const override;
380 
385  virtual Point<spacedim>
386  push_forward(const Point<2> &chart_point) const override;
387 
401  push_forward_gradient(const Point<2> &chart_point) const override;
402 
403  protected:
408  std::tuple<double, double, double, double>
409  get_uv_bounds() const;
410 
414  TopoDS_Face face;
415 
420  double tolerance;
421  };
422 
423 } // namespace OpenCASCADE
424 
427 DEAL_II_NAMESPACE_CLOSE
428 
429 
430 #endif // DEAL_II_WITH_OPENCASCADE
431 #endif // dealii_occ_manifold_lib_h
NormalToMeshProjectionManifold(const TopoDS_Shape &sh, const double tolerance=1e-7)
virtual Point< 2 > pull_back(const Point< spacedim > &space_point) const override
const Tensor< 1, spacedim > direction
Definition: manifold_lib.h:180
virtual Point< spacedim > push_forward(const Point< 1 > &chart_point) const override
NormalProjectionManifold(const TopoDS_Shape &sh, const double tolerance=1e-7)
Definition: manifold_lib.cc:81
virtual std::unique_ptr< Manifold< dim, spacedim > > clone() const override
virtual Point< spacedim > project_to_manifold(const ArrayView< const Point< spacedim >> &surrounding_points, const Point< spacedim > &candidate) const override
virtual Point< spacedim > push_forward(const Point< 2 > &chart_point) const override
virtual Point< spacedim > project_to_manifold(const ArrayView< const Point< spacedim >> &surrounding_points, const Point< spacedim > &candidate) const override
virtual std::unique_ptr< Manifold< dim, spacedim > > clone() const override
virtual std::unique_ptr< Manifold< dim, spacedim > > clone() const override
virtual std::unique_ptr< Manifold< dim, spacedim > > clone() const override
virtual Point< spacedim > project_to_manifold(const ArrayView< const Point< spacedim >> &surrounding_points, const Point< spacedim > &candidate) const override
virtual std::unique_ptr< Manifold< dim, spacedim > > clone() const override
Definition: manifold_lib.cc:94
ArclengthProjectionLineManifold(const TopoDS_Shape &sh, const double tolerance=1e-7)
virtual Point< 1 > pull_back(const Point< spacedim > &space_point) const override
std::tuple< double, double, double, double > get_uv_bounds() const
virtual DerivativeForm< 1, 2, spacedim > push_forward_gradient(const Point< 2 > &chart_point) const override
DirectionalProjectionManifold(const TopoDS_Shape &sh, const Tensor< 1, spacedim > &direction, const double tolerance=1e-7)