From 0598086f2127566a0616fadb4186a323df7c4f5d Mon Sep 17 00:00:00 2001 From: ESeNonFossiIo Date: Tue, 8 Mar 2016 21:48:25 +0100 Subject: [PATCH] NURBSPatchManifold --- include/deal.II/opencascade/boundary_lib.h | 69 ++++++++++++++++++++++ source/opencascade/boundary_lib.cc | 67 +++++++++++++++++++++ 2 files changed, 136 insertions(+) diff --git a/include/deal.II/opencascade/boundary_lib.h b/include/deal.II/opencascade/boundary_lib.h index ed69a17c8d..b99cb3231c 100644 --- a/include/deal.II/opencascade/boundary_lib.h +++ b/include/deal.II/opencascade/boundary_lib.h @@ -24,6 +24,7 @@ #include #include +#include // 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 + class NURBSPatchManifold : public ChartManifold + { + 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 &space_point) const; + + /** + * TODO: + */ + virtual Point + 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 + // get_new_point(const Quadrature &quad) const; + + /** + * TODO + */ + std::tuple + 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 cb077f15a7..786c978ff4 100644 --- a/source/opencascade/boundary_lib.cc +++ b/source/opencascade/boundary_lib.cc @@ -10,6 +10,8 @@ #include #include #include +#include +#include #include #include #include @@ -252,7 +254,72 @@ namespace OpenCASCADE return point(P); } + template + NURBSPatchManifold:: + 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 + 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" -- 2.39.5