NURBSPatchManifold

diff --git a/include/deal.II/opencascade/boundary_lib.h b/include/deal.II/opencascade/boundary_lib.h
index ed69a17c8d3e3d0fb9e6a835ac25adc5b0984824..b99cb3231cb2662d12269474e812021a512a7d05 100644 (file)
--- a/include/deal.II/opencascade/boundary_lib.h
+++ b/include/deal.II/opencascade/boundary_lib.h
@@ -24,6 +24,7 @@
 
 #include <deal.II/opencascade/utilities.h>
 #include <deal.II/grid/tria_boundary.h>
+#include <deal.II/grid/manifold.h>
 
 // opencascade needs "HAVE_CONFIG_H" to be exported...
 #define HAVE_CONFIG_H
@@ -319,6 +320,74 @@ namespace OpenCASCADE
      */
     const double length;
   };
+
+  /**
+   * TODO:
+   *
+   * @ingroup manifold
+   *
+   * @author Luca Heltai, 2014
+   */
+  template <int dim, int spacedim>
+  class NURBSPatchManifold : public ChartManifold<dim, spacedim, 2>
+  {
+  public:
+    /**
+     * TODO:
+     * Asser: count shape!
+     */
+    NURBSPatchManifold(const TopoDS_Face &face, const double tolerance = 1e-7);
+
+    /**
+     * TODO:  3 -> 2
+     */
+    virtual Point<2>
+    pull_back(const Point<spacedim> &space_point) const;
+
+    /**
+     * TODO:
+     */
+    virtual Point<spacedim>
+    push_forward(const Point<2> &chart_point) const;
+
+    /**
+     * Given a point in the spacedim dimensional Euclidean space, this
+     * method returns the derivatives of the function $F$ that maps from
+     * the polar coordinate system to the Euclidean coordinate
+     * system. In other words, it is a matrix of size
+     * $\text{spacedim}\times\text{spacedim}$.
+     *
+     * This function is used in the computations required by the
+     * get_tangent_vector() function.
+     *
+     * Refer to the general documentation of this class for more information.
+     */
+    virtual
+    DerivativeForm<1,dim,spacedim>
+    push_forward_gradient(const Point<2> &chart_point) const;
+
+    // /**
+    //  * Let the new point be the average sum of surrounding vertices.
+    //  *
+    //  * In the two dimensional implementation, we use the pull_back and
+    //  * push_forward mechanism. For three dimensions, this does not work well, so
+    //  * we overload the get_new_point function directly.
+    //  */
+    // virtual Point<spacedim>
+    // get_new_point(const Quadrature<spacedim> &quad) const;
+
+    /**
+     * TODO
+     */
+     std::tuple<double, double, double, double>
+     get_uv_bounds() const;
+     
+  private:
+    TopoDS_Face face;
+
+    double tolerance;
+  };
+
 }
 
 /*@}*/

diff --git a/source/opencascade/boundary_lib.cc b/source/opencascade/boundary_lib.cc
index cb077f15a7406279af6755842226f0c3f146448f..786c978ff4e728c9c4d958420d092a327f5c43b0 100644 (file)
--- a/source/opencascade/boundary_lib.cc
+++ b/source/opencascade/boundary_lib.cc
@@ -10,6 +10,8 @@
 #include <GCPnts_AbscissaPoint.hxx>
 #include <ShapeAnalysis_Curve.hxx>
 #include <BRep_Tool.hxx>
+#include <BRepTools.hxx>
+#include <ShapeAnalysis_Surface.hxx>
 #include <TopoDS.hxx>
 #include <Adaptor3d_HCurve.hxx>
 #include <Handle_Adaptor3d_HCurve.hxx>
@@ -252,7 +254,72 @@ namespace OpenCASCADE
     return point(P);
   }
 
+  template <int dim, int spacedim>
+  NURBSPatchManifold<dim, spacedim>::
+  NURBSPatchManifold( const TopoDS_Face &face,
+                      const double tolerance)
+    :
+    face(face),
+    tolerance(tolerance)
+  {}
+
+  template<>
+  Point<2>
+  NURBSPatchManifold<2, 3>::
+  pull_back(const Point<3> &space_point) const
+  {
+    Handle(Geom_Surface) SurfToProj = BRep_Tool::Surface(face);
+
+    ShapeAnalysis_Surface projector(SurfToProj);
+    gp_Pnt2d proj_params = projector.ValueOfUV(point(space_point), tolerance);
+
+    double u = proj_params.X();
+    double v = proj_params.Y();
+
+    return Point<2>(u,v);
+  }
+
+  template<>
+  Point<3>
+  NURBSPatchManifold<2, 3>::
+  push_forward(const Point<2> &chart_point) const
+  {
+    return ::dealii::OpenCASCADE::push_forward(face, chart_point[0], chart_point[1]);
+  }
+
+  template<>
+  DerivativeForm<1,2,3>
+  NURBSPatchManifold<2, 3>::
+  push_forward_gradient(const Point<2> &chart_point) const
+  {
+    DerivativeForm<1,2,3> DX;
+    Handle(Geom_Surface) surf = BRep_Tool::Surface(face);
+
+    gp_Pnt q;
+    gp_Vec Du, Dv;
+    surf->D1(chart_point[0],chart_point[1], q, Du, Dv);
+
+    DX[0][0] = Du.X();
+    DX[0][1] = Du.Y();
+    DX[0][2] = Du.Z();
+    DX[1][0] = Dv.X();
+    DX[1][1] = Dv.Y();
+    DX[1][2] = Dv.Z();
+
+    return DX;
+  }
+
+  template<>
+  std::tuple<double, double, double, double>
+  NURBSPatchManifold<2, 3>::
+  get_uv_bounds() const
+  {
+    Standard_Real umin, umax, vmin, vmax;
+    BRepTools::UVBounds(face, umin, umax, vmin, vmax);    
+    return std::make_tuple(umin, umax, vmin, vmax);
+  }
 
+template class NURBSPatchManifold<2, 3 >; 
 // Explicit instantiations
 #include "boundary_lib.inst"