]> https://gitweb.dealii.org/ - dealii.git/commitdiff
OpenCASCADE classes and utilities now work also in spacedim=2
authorLuca Heltai <luca.heltai@sissa.it>
Sat, 9 Dec 2017 10:59:19 +0000 (11:59 +0100)
committerLuca Heltai <luca.heltai@sissa.it>
Tue, 12 Dec 2017 12:04:54 +0000 (13:04 +0100)
doc/news/changes/minor/20171208LucaHeltai [new file with mode: 0644]
include/deal.II/opencascade/utilities.h
source/opencascade/CMakeLists.txt
source/opencascade/boundary_lib.cc
source/opencascade/boundary_lib.inst.in
source/opencascade/utilities.cc
source/opencascade/utilities.inst.in [new file with mode: 0644]
tests/opencascade/arclength_boundary_02.cc [new file with mode: 0644]
tests/opencascade/arclength_boundary_02.output [new file with mode: 0644]

diff --git a/doc/news/changes/minor/20171208LucaHeltai b/doc/news/changes/minor/20171208LucaHeltai
new file mode 100644 (file)
index 0000000..be924f8
--- /dev/null
@@ -0,0 +1,3 @@
+New: Now opencascade also works for spacedim == 2. 
+<br> 
+(Luca Heltai, 2017/12/08)
index 8d16be67e2e4c493d024c7e380f347172af7c3fe..2fbe26d6912f221d34252422b2ceddfe04f93a13 100644 (file)
@@ -94,7 +94,11 @@ DEAL_II_NAMESPACE_OPEN
  * cases they can be used in two dimensions as well, and the third dimension
  * will be set to zero.
  *
- * @author Luca Heltai, Andrea Mola, 2011--2014.
+ * If you whish to use these tools when the dimension of the space is two, then
+ * make sure your cad files are actually flat and that all z coordinates are equal
+ * to zero, as otherwise you will get many exceptions.
+ *
+ * @author Luca Heltai, Andrea Mola, 2011--2017.
  */
 namespace OpenCASCADE
 {
@@ -179,9 +183,9 @@ namespace OpenCASCADE
                          const double tolerance=1e-7);
 
   /**
-   * Creates a 3D smooth BSpline curve passing through the points in the
-   * assigned vector, and store it in the returned TopoDS_Shape (which is of
-   * type TopoDS_Edge). The points are reordered internally according to their
+   * Creates a smooth BSpline curve passing through the points in the assigned
+   * vector, and store it in the returned TopoDS_Shape (which is of type
+   * TopoDS_Edge). The points are reordered internally according to their
    * scalar product with the direction, if direction is different from zero,
    * otherwise they are used as passed. Notice that this function changes the
    * input points if required by the algorithm.
@@ -195,8 +199,9 @@ namespace OpenCASCADE
    * points. If the algorithm fails in generating such a curve, an exception
    * is thrown.
    */
-  TopoDS_Edge interpolation_curve(std::vector<Point<3> >  &curve_points,
-                                  const Tensor<1,3> &direction=Tensor<1,3>(),
+  template<int dim>
+  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);
 
@@ -211,16 +216,20 @@ namespace OpenCASCADE
                                   std::vector<TopoDS_Vertex> &vertices);
 
   /**
-   * Create a triangulation from a single face. This class extract the first u
-   * and v parameter of the parametric surface making up this face, and
-   * creates a Triangulation<2,3> containing a single coarse cell reflecting
+   * Create a triangulation from a single face. This class extracts the first u
+   * and v parameter of the parametric surface making up this face, and creates
+   * a Triangulation<2,spacedim> containing a single coarse cell reflecting
    * this face. If the surface is not a trimmed surface, the vertices of this
    * cell will coincide with the TopoDS_Vertex vertices of the original
    * TopoDS_Face. This, however, is often not the case, and the user should be
    * careful on how this mesh is used.
+   *
+   * If you call this function with a Triangulation<2,2>, make sure that the
+   * input face has all z coordinates set to zero, or you'll get an exception.
    */
+  template <int spacedim>
   void create_triangulation(const TopoDS_Face &face,
-                            Triangulation<2,3> &tria);
+                            Triangulation<2,spacedim> &tria);
 
 
   /**
@@ -274,9 +283,10 @@ namespace OpenCASCADE
    * the u coordinate and the v coordinate (which is different from zero only
    * if the resulting shape is a face).
    */
-  std::tuple<Point<3>, TopoDS_Shape, double, double>
+  template <int dim>
+  std::tuple<Point<dim>, TopoDS_Shape, double, double>
   project_point_and_pull_back(const TopoDS_Shape &in_shape,
-                              const Point<3> &origin,
+                              const Point<dim> &origin,
                               const double tolerance=1e-7);
 
   /**
@@ -285,9 +295,10 @@ namespace OpenCASCADE
    * are iterated, faces first, then edges, and the returned point is the
    * closest one to the @p in_shape, regardless of its type.
    */
-  Point<3> closest_point(const TopoDS_Shape &in_shape,
-                         const Point<3> &origin,
-                         const double tolerance=1e-7);
+  template <int dim>
+  Point<dim> closest_point(const TopoDS_Shape &in_shape,
+                           const Point<dim> &origin,
+                           const double tolerance=1e-7);
 
   /**
    * Given an elementary shape @p in_shape and the reference coordinates
@@ -297,9 +308,10 @@ namespace OpenCASCADE
    * used as input to this function. If this is not the case, an Exception is
    * thrown.
    */
-  Point<3> push_forward(const TopoDS_Shape &in_shape,
-                        const double u,
-                        const double v);
+  template <int dim>
+  Point<dim> push_forward(const TopoDS_Shape &in_shape,
+                          const double u,
+                          const double v);
 
 
   /**
@@ -334,22 +346,29 @@ namespace OpenCASCADE
    *
    * The optional @p tolerance parameter is used to compute distances.
    */
-  Point<3> line_intersection(const TopoDS_Shape &in_shape,
-                             const Point<3> &origin,
-                             const Tensor<1,3> &direction,
-                             const double tolerance=1e-7);
+  template <int dim>
+  Point<dim> line_intersection(const TopoDS_Shape &in_shape,
+                               const Point<dim> &origin,
+                               const Tensor<1,dim> &direction,
+                               const double tolerance=1e-7);
 
 
   /**
-   * Convert OpenCASCADE point into a Point<3>.
+   * Convert OpenCASCADE point into a Point<spacedim>.
+   *
+   * The tolerance argument is used to check if the non used components of the
+   * OpenCASCADE point are close to zero. If this is not the case, an assertion
+   * is thrown in debug mode.
    */
-  Point<3> point(const gp_Pnt &p);
+  template <int spacedim>
+  Point<spacedim> point(const gp_Pnt &p, const double &tolerance=1e-10);
 
 
   /**
    * Convert Point<3> into OpenCASCADE point.
    */
-  gp_Pnt point(const Point<3> &p);
+  template <int spacedim>
+  gp_Pnt point(const Point<spacedim> &p);
 
 
   /**
@@ -358,9 +377,10 @@ namespace OpenCASCADE
    * optional parameter is used as a relative tolerance when comparing
    * objects.
    */
-  bool point_compare(const Point<3>    &p1,
-                     const Point<3>    &p2,
-                     const Tensor<1,3> &direction = Tensor<1,3>(),
+  template <int dim>
+  bool point_compare(const Point<dim> &p1,
+                     const Point<dim> &p2,
+                     const Tensor<1,dim> &direction = Tensor<1,dim>(),
                      const double       tolerance = 1e-10);
 
 
@@ -368,16 +388,18 @@ namespace OpenCASCADE
    * Exception thrown when the point specified as argument does not lie
    * between @p tolerance from the given TopoDS_Shape.
    */
+  template <int dim>
   DeclException1 (ExcPointNotOnManifold,
-                  Point<3>,
+                  Point<dim>,
                   <<"The point [ "<<arg1<<" ] is not on the manifold.");
 
   /**
    * Exception thrown when the point specified as argument cannot be projected
    * to the manifold.
    */
+  template <int dim>
   DeclException1 (ExcProjectionFailed,
-                  Point<3>,
+                  Point<dim>,
                   <<"Projection of point [ "<< arg1
                   << " ] failed.");
 
index 9a130423538c81482d737c4c4ca36df34c461cab..0d9905db344c587567f2e7a8b5de2f50238b407c 100644 (file)
@@ -22,6 +22,7 @@ SET(_src
 
 SET(_inst
   boundary_lib.inst.in
+  utilities.inst.in
   )
 
 FILE(GLOB _header
index 8125c5bfd233cdd6a991e72f2676127e9a495097..1545c0f4c32824fc0b657dfb6672017d6f1d7a24 100644 (file)
@@ -99,7 +99,7 @@ namespace OpenCASCADE
       Assert(closest_point(sh, surrounding_points[i], tolerance)
              .distance(surrounding_points[i]) <
              std::max(tolerance*surrounding_points[i].norm(), tolerance),
-             ExcPointNotOnManifold(surrounding_points[i]));
+             ExcPointNotOnManifold<spacedim>(surrounding_points[i]));
 #endif
     return closest_point(sh, candidate,tolerance);
   }
@@ -129,7 +129,7 @@ namespace OpenCASCADE
       Assert(closest_point(sh, surrounding_points[i],tolerance)
              .distance(surrounding_points[i]) <
              std::max(tolerance*surrounding_points[i].norm(), tolerance),
-             ExcPointNotOnManifold(surrounding_points[i]));
+             ExcPointNotOnManifold<spacedim>(surrounding_points[i]));
 #endif
     return line_intersection(sh, candidate, direction, tolerance);
   }
@@ -162,7 +162,7 @@ namespace OpenCASCADE
         Assert(closest_point(sh, surrounding_points[i], tolerance)
                .distance(surrounding_points[i]) <
                std::max(tolerance*surrounding_points[i].norm(), tolerance),
-               ExcPointNotOnManifold(surrounding_points[i]));
+               ExcPointNotOnManifold<spacedim>(surrounding_points[i]));
       }
 #endif
 
@@ -248,9 +248,7 @@ namespace OpenCASCADE
     curve(curve_adaptor(sh)),
     tolerance(tolerance),
     length(shape_length(sh))
-  {
-    Assert(spacedim == 3, ExcNotImplemented());
-  }
+  {}
 
 
   template <int dim, int spacedim>
@@ -261,7 +259,7 @@ namespace OpenCASCADE
     ShapeAnalysis_Curve curve_analysis;
     gp_Pnt proj;
     const double dist = curve_analysis.Project(curve->GetCurve(), point(space_point), tolerance, proj, t, true);
-    Assert(dist < tolerance*length, ExcPointNotOnManifold(space_point));
+    Assert(dist < tolerance*length, ExcPointNotOnManifold<spacedim>(space_point));
     (void)dist; // Silence compiler warning in Release mode.
     return Point<1>(GCPnts_AbscissaPoint::Length(curve->GetCurve(),curve->GetCurve().FirstParameter(),t));
   }
@@ -274,7 +272,7 @@ namespace OpenCASCADE
   {
     GCPnts_AbscissaPoint AP(curve->GetCurve(), chart_point[0], curve->GetCurve().FirstParameter());
     gp_Pnt P = curve->GetCurve().Value(AP.Parameter());
-    return point(P);
+    return point<spacedim>(P);
   }
 
   template <int dim, int spacedim>
@@ -286,10 +284,9 @@ namespace OpenCASCADE
     tolerance(tolerance)
   {}
 
-  template <>
-  Point<2>
-  NURBSPatchManifold<2, 3>::
-  pull_back(const Point<3> &space_point) const
+  template <int dim, int spacedim>  Point<2>
+  NURBSPatchManifold<dim, spacedim>::
+  pull_back(const Point<spacedim> &space_point) const
   {
     Handle(Geom_Surface) SurfToProj = BRep_Tool::Surface(face);
 
@@ -302,20 +299,20 @@ namespace OpenCASCADE
     return Point<2>(u,v);
   }
 
-  template <>
-  Point<3>
-  NURBSPatchManifold<2, 3>::
+  template <int dim, int spacedim>
+  Point<spacedim>
+  NURBSPatchManifold<dim, spacedim>::
   push_forward(const Point<2> &chart_point) const
   {
-    return ::dealii::OpenCASCADE::push_forward(face, chart_point[0], chart_point[1]);
+    return ::dealii::OpenCASCADE::push_forward<spacedim>(face, chart_point[0], chart_point[1]);
   }
 
-  template <>
-  DerivativeForm<1,2,3>
-  NURBSPatchManifold<2, 3>::
+  template <int dim, int spacedim>
+  DerivativeForm<1,2,spacedim>
+  NURBSPatchManifold<dim,spacedim>::
   push_forward_gradient(const Point<2> &chart_point) const
   {
-    DerivativeForm<1,2,3> DX;
+    DerivativeForm<1,2,spacedim> DX;
     Handle(Geom_Surface) surf = BRep_Tool::Surface(face);
 
     gp_Pnt q;
@@ -324,17 +321,24 @@ namespace OpenCASCADE
 
     DX[0][0] = Du.X();
     DX[1][0] = Du.Y();
-    DX[2][0] = Du.Z();
+    if (spacedim>2)
+      DX[2][0] = Du.Z();
+    else
+      Assert(std::abs(Du.Z()) < tolerance,
+             ExcMessage("Expecting derivative along Z to be zero! Bailing out."));
     DX[0][1] = Dv.X();
     DX[1][1] = Dv.Y();
-    DX[2][1] = Dv.Z();
-
+    if (spacedim>2)
+      DX[2][1] = Dv.Z();
+    else
+      Assert(std::abs(Dv.Z()) < tolerance,
+             ExcMessage("Expecting derivative along Z to be zero! Bailing out."));
     return DX;
   }
 
-  template <>
+  template <int dim, int spacedim>
   std::tuple<double, double, double, double>
-  NURBSPatchManifold<2, 3>::
+  NURBSPatchManifold<dim,spacedim>::
   get_uv_bounds() const
   {
     Standard_Real umin, umax, vmin, vmax;
@@ -343,7 +347,6 @@ namespace OpenCASCADE
   }
 
 // Explicit instantiations
-  template class NURBSPatchManifold<2, 3 >;
 #include "boundary_lib.inst"
 
 } // end namespace OpenCASCADE
index a43f51b6888f3732cc070c0974f52f7aee55829b..6514422184a3af9d1f583febf0f91dd26fb9e68e 100644 (file)
@@ -21,6 +21,12 @@ for (deal_II_dimension : DIMENSIONS)
     template class DirectionalProjectionBoundary<deal_II_dimension, 3>;
     template class NormalToMeshProjectionBoundary<deal_II_dimension, 3>;
     template class ArclengthProjectionLineManifold<deal_II_dimension, 3>;
+    template class NURBSPatchManifold<deal_II_dimension, 3>;
+#if deal_II_dimension <= 2
+    template class DirectionalProjectionBoundary<deal_II_dimension, 2>;
+    template class ArclengthProjectionLineManifold<deal_II_dimension, 2>;
+    template class NURBSPatchManifold<deal_II_dimension, 2>;
+#endif
 }
 
 
index eefb3b44f35b55031b16f1271d4866e84b6b7744..bd83e5f57104cba9db01374581463f1900364e4b 100644 (file)
@@ -165,27 +165,52 @@ namespace OpenCASCADE
       }
   }
 
-  gp_Pnt point(const Point<3> &p)
+  template <int spacedim>
+  gp_Pnt point(const Point<spacedim> &p)
   {
-    return gp_Pnt(p(0), p(1), p(2));
+    switch (spacedim)
+      {
+      case 1:
+        return gp_Pnt(p[0], 0, 0);
+      case 2:
+        return gp_Pnt(p[0], p[1], 0);
+      case 3:
+        return gp_Pnt(p[0], p[1], p[2]);
+      }
   }
 
-
-  Point<3> point(const gp_Pnt &p)
+  template <int spacedim>
+  Point<spacedim> point(const gp_Pnt &p, const double &tolerance)
   {
-    return Point<3>(p.X(), p.Y(), p.Z());
+    (void) tolerance;
+    switch (spacedim)
+      {
+      case 1:
+        Assert(std::abs(p.Y()) <= tolerance,
+               ExcMessage("Cannot convert OpenCASCADE point to 1d if p.Y() != 0."));
+        Assert(std::abs(p.Z()) <= tolerance,
+               ExcMessage("Cannot convert OpenCASCADE point to 1d if p.Z() != 0."));
+        return Point<spacedim>(p.X());
+      case 2:
+        Assert(std::abs(p.Z()) <= tolerance,
+               ExcMessage("Cannot convert OpenCASCADE point to 2d if p.Z() != 0."));
+        return Point<spacedim>(p.X(), p.Y());
+      case 3:
+        return Point<spacedim>(p.X(), p.Y(), p.Z());
+      }
   }
 
-  bool point_compare(const Point<3>    &p1,
-                     const Point<3>    &p2,
-                     const Tensor<1,3> &direction,
+  template<int dim>
+  bool point_compare(const Point<dim>    &p1,
+                     const Point<dim>    &p2,
+                     const Tensor<1,dim> &direction,
                      const double       tolerance)
   {
     const double rel_tol=std::max(tolerance, std::max(p1.norm(), p2.norm())*tolerance);
     if (direction.norm() > 0.0)
       return (p1*direction < p2*direction-rel_tol);
     else
-      for (int d=2; d>=0; --d)
+      for (int d=dim; d>=0; --d)
         if (p1[d] < p2[d]-rel_tol)
           return true;
         else if (p2[d] < p1[d]-rel_tol)
@@ -387,16 +412,16 @@ namespace OpenCASCADE
     return out_shape;
   }
 
-
-  Point<3> line_intersection(const TopoDS_Shape &in_shape,
-                             const Point<3> &origin,
-                             const Tensor<1,3> &direction,
-                             const double tolerance)
+  template <int dim>
+  Point<dim> line_intersection(const TopoDS_Shape &in_shape,
+                               const Point<dim> &origin,
+                               const Tensor<1,dim> &direction,
+                               const double tolerance)
   {
     // translating original Point<dim> to gp point
 
     gp_Pnt P0 = point(origin);
-    gp_Ax1 gpaxis(P0, gp_Dir(direction[0], direction[1], direction[2]));
+    gp_Ax1 gpaxis(P0, gp_Dir(direction[0], dim > 1 ? direction[1] : 0, dim>2 ? direction[2] : 0));
     gp_Lin line(gpaxis);
 
     // destination point
@@ -413,7 +438,7 @@ namespace OpenCASCADE
     Assert(Inters.IsDone(), ExcMessage("Could not project point."));
 
     double minDistance = 1e7;
-    Point<3> result;
+    Point<dim> result;
     for (int i=0; i<Inters.NbPnt(); ++i)
       {
         const double distance = point(origin).Distance(Inters.Pnt(i+1));
@@ -421,15 +446,16 @@ namespace OpenCASCADE
         if (distance < minDistance)
           {
             minDistance = distance;
-            result = point(Inters.Pnt(i+1));
+            result = point<dim>(Inters.Pnt(i+1));
           }
       }
 
     return result;
   }
 
-  TopoDS_Edge interpolation_curve(std::vector<Point<3> > &curve_points,
-                                  const Tensor<1,3> &direction,
+  template <int dim>
+  TopoDS_Edge interpolation_curve(std::vector<Point<dim> > &curve_points,
+                                  const Tensor<1,dim> &direction,
                                   const bool closed,
                                   const double tolerance)
   {
@@ -439,7 +465,7 @@ namespace OpenCASCADE
     if (direction*direction > 0)
       {
         std::sort(curve_points.begin(), curve_points.end(),
-                  [&](const Point<3> &p1, const Point<3> &p2)
+                  [&](const Point<dim> &p1, const Point<dim> &p2)
         {
           return OpenCASCADE::point_compare(p1, p2, direction, tolerance);
         });
@@ -527,15 +553,12 @@ namespace OpenCASCADE
         unsigned int point_index             = start_point_index;
 
         // point_index and face_index always run together
-        std::vector<Point<3> > pointlist;
+        std::vector<Point<spacedim> > pointlist;
         do
           {
             visited_faces[face_index]            = true;
             auto current_point = vert_to_point[point_index];
-            if (spacedim==2)
-              pointlist.push_back(Point<3>(current_point[0],current_point[1],0));
-            else
-              pointlist.push_back(Point<3>(current_point[0],current_point[1],current_point[2]));
+            pointlist.push_back(current_point);
 
             // Get next point
             if (face_to_verts[face_index].first != point_index )
@@ -551,7 +574,7 @@ namespace OpenCASCADE
           }
         while (point_index != start_point_index);
 
-        interpolation_curves.push_back(interpolation_curve(pointlist, Tensor<1,3>(), true));
+        interpolation_curves.push_back(interpolation_curve(pointlist, Tensor<1,spacedim>(), true));
 
         finished = true;
         for (const auto &f : visited_faces)
@@ -566,11 +589,10 @@ namespace OpenCASCADE
   }
 
 
-
-  std::tuple<Point<3>, TopoDS_Shape, double, double>
-  project_point_and_pull_back(const TopoDS_Shape &in_shape,
-                              const Point<3> &origin,
-                              const double tolerance)
+  template<int dim>
+  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)
   {
     TopExp_Explorer exp;
     gp_Pnt Pproj = point(origin);
@@ -598,7 +620,7 @@ namespace OpenCASCADE
 
         SurfToProj->D0(proj_params.X(),proj_params.Y(),tmp_proj);
 
-        double distance = point(tmp_proj).distance(origin);
+        double distance = point<dim>(tmp_proj).distance(origin);
         if (distance < minDistance)
           {
             minDistance = distance;
@@ -643,16 +665,17 @@ namespace OpenCASCADE
         }
 
     Assert(counter > 0, ExcMessage("Could not find projection points."));
-    return std::tuple<Point<3>, TopoDS_Shape, double, double>
-           (point(Pproj),out_shape, u, v);
+    return std::tuple<Point<dim>, TopoDS_Shape, double, double>
+           (point<dim>(Pproj),out_shape, u, v);
   }
 
 
-  Point<3> closest_point(const TopoDS_Shape &in_shape,
-                         const Point<3> &origin,
-                         const double tolerance)
+  template<int dim>
+  Point<dim> closest_point(const TopoDS_Shape &in_shape,
+                           const Point<dim> &origin,
+                           const double tolerance)
   {
-    std::tuple<Point<3>, TopoDS_Shape, double, double>
+    std::tuple<Point<dim>, TopoDS_Shape, double, double>
     ref = project_point_and_pull_back(in_shape, origin, tolerance);
     return std::get<0>(ref);
   }
@@ -690,26 +713,27 @@ namespace OpenCASCADE
     return push_forward_and_differential_forms(face, u, v, tolerance);
   }
 
-  Point<3> push_forward(const TopoDS_Shape &in_shape,
-                        const double u,
-                        const double v)
+  template<int dim>
+  Point<dim> push_forward(const TopoDS_Shape &in_shape,
+                          const double u,
+                          const double v)
   {
     switch (in_shape.ShapeType())
       {
       case TopAbs_FACE:
       {
         BRepAdaptor_Surface surf(TopoDS::Face(in_shape));
-        return point(surf.Value(u,v));
+        return point<dim>(surf.Value(u,v));
       }
       case TopAbs_EDGE:
       {
         BRepAdaptor_Curve curve(TopoDS::Edge(in_shape));
-        return point(curve.Value(u));
+        return point<dim>(curve.Value(u));
       }
       default:
         Assert(false, ExcUnsupportedShape());
       }
-    return Point<3>();
+    return Point<dim>();
   }
 
   std::tuple<Point<3>,  Tensor<1,3>, double, double>
@@ -739,13 +763,14 @@ namespace OpenCASCADE
         Max_Curvature *= -1;
       }
 
-    return std::tuple<Point<3>, Tensor<1,3>, double, double>(point(Value), normal, Min_Curvature, Max_Curvature);
+    return std::tuple<Point<3>, Tensor<1,3>, double, double>(point<3>(Value), normal, Min_Curvature, Max_Curvature);
   }
 
 
 
+  template<int spacedim>
   void create_triangulation(const TopoDS_Face &face,
-                            Triangulation<2,3> &tria)
+                            Triangulation<2,spacedim> &tria)
   {
     BRepAdaptor_Surface surf(face);
     const double u0 = surf.FirstUParameter();
@@ -754,13 +779,13 @@ namespace OpenCASCADE
     const double v1 = surf.LastVParameter();
 
     std::vector<CellData<2> > cells;
-    std::vector<Point<3> > vertices;
+    std::vector<Point<spacedim> > vertices;
     SubCellData t;
 
-    vertices.push_back(point(surf.Value(u0,v0)));
-    vertices.push_back(point(surf.Value(u1,v0)));
-    vertices.push_back(point(surf.Value(u0,v1)));
-    vertices.push_back(point(surf.Value(u1,v1)));
+    vertices.push_back(point<spacedim>(surf.Value(u0,v0)));
+    vertices.push_back(point<spacedim>(surf.Value(u1,v0)));
+    vertices.push_back(point<spacedim>(surf.Value(u0,v1)));
+    vertices.push_back(point<spacedim>(surf.Value(u1,v1)));
 
     CellData<2> cell;
     for (unsigned int i=0; i<4; ++i)
@@ -770,22 +795,7 @@ namespace OpenCASCADE
     tria.create_triangulation(vertices, cells, t);
   }
 
-
-
-  // Explicit instantiations
-  template
-  std::vector<TopoDS_Edge> create_curves_from_triangulation_boundary
-  (const Triangulation<2, 2> &triangulation,
-   const Mapping<2, 2> &mapping);
-
-
-
-  template
-  std::vector<TopoDS_Edge> create_curves_from_triangulation_boundary
-  (const Triangulation<2, 3> &triangulation,
-   const Mapping<2, 3> &mapping);
-
-
+#include "utilities.inst"
 
 } // end namespace
 
diff --git a/source/opencascade/utilities.inst.in b/source/opencascade/utilities.inst.in
new file mode 100644 (file)
index 0000000..7fae0b9
--- /dev/null
@@ -0,0 +1,62 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 1998 - 2017 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE at
+// the top level of the deal.II distribution.
+//
+// ---------------------------------------------------------------------
+
+for (deal_II_dimension : DIMENSIONS)
+{
+#if deal_II_dimension > 1
+    // Explicit instantiations for dim = 2 and 3
+    template
+    std::vector<TopoDS_Edge> create_curves_from_triangulation_boundary
+    (const Triangulation<2, deal_II_dimension> &triangulation,
+     const Mapping<2, deal_II_dimension> &mapping);
+
+    template bool point_compare(const Point<deal_II_dimension>& p1,
+                                const Point<deal_II_dimension>& p2,
+                                const Tensor<1,deal_II_dimension>& direction,
+                                const double       tolerance );
+
+    template Point<deal_II_dimension> point(const gp_Pnt &p, const double &tolerance);
+
+    template  gp_Pnt point(const Point<deal_II_dimension> &p);
+
+    template
+    TopoDS_Edge interpolation_curve(std::vector<Point<deal_II_dimension> >& curve_points,
+                                    const Tensor<1, deal_II_dimension>& direction=Tensor<1,deal_II_dimension>(),
+                                    const bool closed,
+                                    const double tolerance);
+
+    template Point<deal_II_dimension> push_forward(const TopoDS_Shape &in_shape,
+            const double u,
+            const double v);
+
+    template Point<deal_II_dimension> line_intersection(const TopoDS_Shape &in_shape,
+            const Point<deal_II_dimension> &origin,
+            const Tensor<1,deal_II_dimension> &direction,
+            const double tolerance);
+
+    template void create_triangulation(const TopoDS_Face &face,
+                                       Triangulation<2,deal_II_dimension> &tria);
+
+    template std::tuple<Point<deal_II_dimension>, TopoDS_Shape, double, double>
+    project_point_and_pull_back(const TopoDS_Shape &in_shape,
+                                const Point<deal_II_dimension>& origin,
+                                const double tolerance);
+
+    template Point<deal_II_dimension> closest_point(const TopoDS_Shape &in_shape,
+            const Point<deal_II_dimension> &origin,
+            const double tolerance);
+
+#endif
+}
diff --git a/tests/opencascade/arclength_boundary_02.cc b/tests/opencascade/arclength_boundary_02.cc
new file mode 100644 (file)
index 0000000..47ffb03
--- /dev/null
@@ -0,0 +1,56 @@
+//-----------------------------------------------------------
+//
+//    Copyright (C) 2014 - 2015 by the deal.II authors
+//
+//    This file is subject to LGPL and may not be distributed
+//    without copyright and license information. Please refer
+//    to the file deal.II/doc/license.html for the  text  and
+//    further information on this license.
+//
+//-----------------------------------------------------------
+
+#include "../tests.h"
+
+#include <deal.II/opencascade/utilities.h>
+#include <deal.II/opencascade/boundary_lib.h>
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/grid_out.h>
+#include <deal.II/grid/grid_generator.h>
+
+#include <TopTools.hxx>
+#include <TopoDS_Shape.hxx>
+#include <Standard_Stream.hxx>
+
+// Create a Triangulation, interpolate its boundary points to a smooth
+// BSpline, and use that as an arlength Boundary Descriptor.
+
+using namespace OpenCASCADE;
+
+int main ()
+{
+  std::ofstream logfile("output");
+
+  // Create a bspline passign through the points
+  std::vector<Point<2> > pts;
+  pts.push_back(Point<2>(0,0));
+  pts.push_back(Point<2>(0,1));
+  pts.push_back(Point<2>(1,1));
+  pts.push_back(Point<2>(1,0));
+
+  TopoDS_Edge edge = interpolation_curve(pts, Tensor<1,2>(), true);
+  ArclengthProjectionLineManifold<2,2> manifold(edge);
+
+  Triangulation<2> tria;
+  GridGenerator::hyper_cube(tria);
+
+  tria.set_all_manifold_ids_on_boundary(1);
+  tria.begin_active()->face(2)->set_manifold_id(numbers::invalid_manifold_id);
+  tria.set_manifold(1, manifold);
+
+  tria.refine_global(3);
+  GridOut gridout;
+  gridout.write_msh(tria, logfile);
+
+  return 0;
+}
+
diff --git a/tests/opencascade/arclength_boundary_02.output b/tests/opencascade/arclength_boundary_02.output
new file mode 100644 (file)
index 0000000..fea54f2
--- /dev/null
@@ -0,0 +1,151 @@
+$NOD
+81
+1  0 0 0
+2  1 0 0
+3  0 1 0
+4  1 1 0
+5  0.5 0 0
+6  -0.251337 0.515857 0
+7  1.16406 0.505545 0
+8  0.504074 1.14672 0
+9  0.4896 0.521016 0
+10  0.25 0 0
+11  0.75 0 0
+12  -0.158882 0.242372 0
+13  -0.20277 0.796991 0
+14  1.11495 0.242077 0
+15  1.13535 0.7715 0
+16  0.242844 1.10455 0
+17  0.767092 1.12113 0
+18  0.4948 0.260508 0
+19  0.496837 0.833869 0
+20  0.119131 0.518436 0
+21  0.826831 0.51328 0
+22  0.180414 0.257274 0
+23  0.79253 0.255303 0
+24  0.174797 0.80468 0
+25  0.79798 0.801633 0
+26  0.125 0 0
+27  0.375 0 0
+28  0.625 0 0
+29  0.875 0 0
+30  -0.08483 0.11767 0
+31  -0.217128 0.375124 0
+32  -0.249712 0.660431 0
+33  -0.114709 0.911629 0
+34  1.06551 0.117217 0
+35  1.1481 0.372213 0
+36  1.16089 0.639761 0
+37  1.08298 0.894921 0
+38  0.117841 1.06055 0
+39  0.372227 1.13328 0
+40  0.636548 1.14359 0
+41  0.891132 1.07494 0
+42  0.4974 0.130254 0
+43  0.4922 0.390762 0
+44  0.493218 0.677443 0
+45  0.500456 0.990296 0
+46  -0.0661032 0.517147 0
+47  0.304365 0.519726 0
+48  0.658215 0.517148 0
+49  0.995446 0.509412 0
+50  0.215207 0.128637 0
+51  0.149773 0.387855 0
+52  0.0107661 0.249823 0
+53  0.337607 0.258891 0
+54  0.771265 0.127652 0
+55  0.80968 0.384292 0
+56  0.643665 0.257905 0
+57  0.953738 0.24869 0
+58  0.146964 0.661558 0
+59  0.208821 0.954615 0
+60  -0.0139862 0.800836 0
+61  0.335817 0.819275 0
+62  0.812405 0.657456 0
+63  0.782536 0.961381 0
+64  0.647409 0.817751 0
+65  0.966663 0.786567 0
+66  0.0672094 0.124472 0
+67  0.356303 0.129445 0
+68  -0.033678 0.381489 0
+69  0.320986 0.389308 0
+70  0.634332 0.128953 0
+71  0.915373 0.123867 0
+72  0.65094 0.387527 0
+73  0.978892 0.378252 0
+74  -0.0513737 0.660995 0
+75  0.320091 0.6695 0
+76  0.0516048 0.929231 0
+77  0.354022 0.976278 0
+78  0.652812 0.667449 0
+79  0.986646 0.648609 0
+80  0.641978 0.980671 0
+81  0.927966 0.926509 0
+$ENDNOD
+$ELM
+64
+1 3 0 0 4 1 26 66 30 
+2 3 0 0 4 26 10 50 66 
+3 3 0 0 4 30 66 52 12 
+4 3 0 0 4 66 50 22 52 
+5 3 0 0 4 10 27 67 50 
+6 3 0 0 4 27 5 42 67 
+7 3 0 0 4 50 67 53 22 
+8 3 0 0 4 67 42 18 53 
+9 3 0 0 4 12 52 68 31 
+10 3 0 0 4 52 22 51 68 
+11 3 0 0 4 31 68 46 6 
+12 3 0 0 4 68 51 20 46 
+13 3 0 0 4 22 53 69 51 
+14 3 0 0 4 53 18 43 69 
+15 3 0 0 4 51 69 47 20 
+16 3 0 0 4 69 43 9 47 
+17 3 0 0 4 5 28 70 42 
+18 3 0 0 4 28 11 54 70 
+19 3 0 0 4 42 70 56 18 
+20 3 0 0 4 70 54 23 56 
+21 3 0 0 4 11 29 71 54 
+22 3 0 0 4 29 2 34 71 
+23 3 0 0 4 54 71 57 23 
+24 3 0 0 4 71 34 14 57 
+25 3 0 0 4 18 56 72 43 
+26 3 0 0 4 56 23 55 72 
+27 3 0 0 4 43 72 48 9 
+28 3 0 0 4 72 55 21 48 
+29 3 0 0 4 23 57 73 55 
+30 3 0 0 4 57 14 35 73 
+31 3 0 0 4 55 73 49 21 
+32 3 0 0 4 73 35 7 49 
+33 3 0 0 4 6 46 74 32 
+34 3 0 0 4 46 20 58 74 
+35 3 0 0 4 32 74 60 13 
+36 3 0 0 4 74 58 24 60 
+37 3 0 0 4 20 47 75 58 
+38 3 0 0 4 47 9 44 75 
+39 3 0 0 4 58 75 61 24 
+40 3 0 0 4 75 44 19 61 
+41 3 0 0 4 13 60 76 33 
+42 3 0 0 4 60 24 59 76 
+43 3 0 0 4 33 76 38 3 
+44 3 0 0 4 76 59 16 38 
+45 3 0 0 4 24 61 77 59 
+46 3 0 0 4 61 19 45 77 
+47 3 0 0 4 59 77 39 16 
+48 3 0 0 4 77 45 8 39 
+49 3 0 0 4 9 48 78 44 
+50 3 0 0 4 48 21 62 78 
+51 3 0 0 4 44 78 64 19 
+52 3 0 0 4 78 62 25 64 
+53 3 0 0 4 21 49 79 62 
+54 3 0 0 4 49 7 36 79 
+55 3 0 0 4 62 79 65 25 
+56 3 0 0 4 79 36 15 65 
+57 3 0 0 4 19 64 80 45 
+58 3 0 0 4 64 25 63 80 
+59 3 0 0 4 45 80 40 8 
+60 3 0 0 4 80 63 17 40 
+61 3 0 0 4 25 65 81 63 
+62 3 0 0 4 65 15 37 81 
+63 3 0 0 4 63 81 41 17 
+64 3 0 0 4 81 37 4 41 
+$ENDELM

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.