]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Fixed also last comments by @bangerth, using std::tuple where appropriate. We split...
authorLuca Heltai <luca.heltai@sissa.it>
Thu, 16 Oct 2014 19:30:30 +0000 (21:30 +0200)
committerLuca Heltai <luca.heltai@sissa.it>
Fri, 31 Oct 2014 10:07:12 +0000 (11:07 +0100)
Added documentation directory to step-54.

Replaced std::tuple with std_cxx11::tuple

Fixed std_cxx11 also in the tests.

15 files changed:
examples/step-54/doc/builds-on [new file with mode: 0644]
examples/step-54/doc/intro.dox [new file with mode: 0644]
examples/step-54/doc/kind [new file with mode: 0644]
examples/step-54/doc/results.dox [new file with mode: 0644]
examples/step-54/doc/tooltip [new file with mode: 0644]
include/deal.II/opencascade/boundary_lib.h
include/deal.II/opencascade/utilities.h
source/opencascade/boundary_lib.cc
source/opencascade/utilities.cc
tests/opencascade/closest_point_to_shape.cc
tests/opencascade/create_tria_01.cc
tests/opencascade/create_tria_01.output
tests/opencascade/iges_describe.cc
tests/opencascade/normal_projection_01.cc
tests/opencascade/plane_intersection_01.cc

diff --git a/examples/step-54/doc/builds-on b/examples/step-54/doc/builds-on
new file mode 100644 (file)
index 0000000..2302bdb
--- /dev/null
@@ -0,0 +1 @@
+step-49
diff --git a/examples/step-54/doc/intro.dox b/examples/step-54/doc/intro.dox
new file mode 100644 (file)
index 0000000..e8b6366
--- /dev/null
@@ -0,0 +1,25 @@
+<br>
+
+<i>This program was contributed by Andrea Mola and Luca Heltai.</i>
+
+@note This program elaborates on concepts of industrial geometry, using 
+tools that interface with OpenCASCADE (http://www.opencascade.org) that 
+allow the specification of arbitrary IGES files as boundaries for your
+geometries.
+
+
+<a name="Intro"></a>
+<h1>Introduction</h1>
+
+In some of the previous tutorial programs (step-1 and step-5 among others)
+we have seen that the deal.II mesh refinement is able to 
+
+-mesh refinement -> rule to refine boundary cells
+-set_boundary
+-so far only elementary shapes were considered
+
+-this tutorial shows how to refine a mesh on top of a complex surface
+generated with a CAD modeler and saved in an IGES file.
+
+
+<h2>Where geometry and meshes intersect</h2>
diff --git a/examples/step-54/doc/kind b/examples/step-54/doc/kind
new file mode 100644 (file)
index 0000000..c1d9154
--- /dev/null
@@ -0,0 +1 @@
+techniques
diff --git a/examples/step-54/doc/results.dox b/examples/step-54/doc/results.dox
new file mode 100644 (file)
index 0000000..ef0e365
--- /dev/null
@@ -0,0 +1,11 @@
+<h1>Results</h1>
+
+Running the program produces a series of mesh files <code>mesh.vtu</code> 
+that we can
+visualize with any of the usual visualization programs that can read the VTK
+file format. 
+
+<img src="http://www.dealii.org/images/steps/developer/step-54.mesh.png" alt="">
+<img src="http://www.dealii.org/images/steps/developer/step-54.mesh.png" alt="">
+<img src="http://www.dealii.org/images/steps/developer/step-54.mesh.png" alt="">
+<img src="http://www.dealii.org/images/steps/developer/step-54.mesh.png" alt="">
\ No newline at end of file
diff --git a/examples/step-54/doc/tooltip b/examples/step-54/doc/tooltip
new file mode 100644 (file)
index 0000000..9afee0b
--- /dev/null
@@ -0,0 +1 @@
+Geometry: Using industry standard IGES files as boundary descriptors.
index 4a0553d91f81ca8f7740839b002b54c4c0b88845..cffff0fe79236362309ffd446e0bd8c696e19c15 100644 (file)
@@ -124,7 +124,7 @@ namespace OpenCASCADE
    * results if the triangulation to be refined is close to the
    * boundary of the given TopoDS_Shape, or when the direction you use
    * at construction time does not intersect the shape. An exception
-   * is thrown when this appens.
+   * is thrown when this happens.
    *
    * @author Luca Heltai, Andrea Mola, 2011--2014.
    */
@@ -177,9 +177,11 @@ namespace OpenCASCADE
   /**
    * A Boundary object based on OpenCASCADE TopoDS_Shape where new
    * points are first computed by averaging the surrounding points in
-   * the same way as FlatManifold does, and then projecting it using
+   * the same way as FlatManifold does, and then projecting them using
    * OpenCASCADE utilities onto the manifold along a direction which
-   * is an estimation of the surrounding points (hence mesh cell) normal.
+   * is an estimation of the surrounding points (hence mesh cell)
+   * normal.
+   *
    * The direction normal to the mesh is particularly useful because
    * it is the direction in which the mesh is missing nodes. For
    * instance, during the refinement of a cell a new node is initially
@@ -194,6 +196,7 @@ namespace OpenCASCADE
    * triangles are identified with the points assigned, and the normals
    * of such triangles are averaged to obtain the approximation of
    * the normal to the cell.
+   *
    * The case in which 2 surrounding points are present (i.e.:a cell
    * edge is being refined) is of course more tricky. The average of
    * the CAD surface normals at the 2 surrounding points is first
@@ -202,23 +205,21 @@ namespace OpenCASCADE
    * to have the new point with equal distance with respect to the
    * surrounding points
    *
-   * This class only operates with CAD faces and makes the
-   * assumption that the shape you pass to it
-   * contains at least a face. If that is not the case, an Exeption
-   * is thrown. In debug mode there is a sanity check to
-   * make sure that the surrounding points (the ones used in
-   * project_to_manifold()) actually live on the Manifold, i.e.,
-   * calling OpenCASCADE::closest_point() on those points leaves them
-   * untouched. If this is not the case, an ExcPointNotOnManifold is
-   * thrown.
+   * This class only operates with CAD faces and makes the assumption
+   * that the shape you pass to it contains at least one face. If that
+   * is not the case, an Exception is thrown. In debug mode there is a
+   * sanity check to make sure that the surrounding points (the ones
+   * used in project_to_manifold()) actually live on the Manifold,
+   * i.e., calling OpenCASCADE::closest_point() on those points leaves
+   * them untouched. If this is not the case, an ExcPointNotOnManifold
+   * is thrown.
    *
    *
    * Notice that this type of Boundary descriptor may fail to give
    * results if the triangulation to be refined is close to the
    * boundary of the given TopoDS_Shape, or when the normal direction
-   * estimated from the surrounding points does not intersect the shape.
-   * An exception
-   * is thrown when this appens.
+   * estimated from the surrounding points does not intersect the
+   * shape.  An exception is thrown when this happens.
    *
    * @author Luca Heltai, Andrea Mola, 2011--2014.
    */
@@ -228,7 +229,8 @@ namespace OpenCASCADE
   public:
     /**
      * Construct a Boundary object which will project points on the
-     * TopoDS_Shape @p sh, along the given @p direction.
+     * TopoDS_Shape @p sh, along a direction which is approximately
+     * normal to the mesh cell.
      */
     NormalToMeshProjectionBoundary(const TopoDS_Shape &sh,
                                    const double tolerance=1e-7);
@@ -238,9 +240,6 @@ namespace OpenCASCADE
      * in debug mode, checks that each of the @p surrounding_points is
      * within tolerance from the given TopoDS_Shape. If this is not
      * the case, an exception is thrown.
-     *
-     * The projected point is computed using OpenCASCADE directional
-     * projection algorithms.
      */
     virtual Point<spacedim>
     project_to_manifold (const std::vector<Point<spacedim> > &surrounding_points,
index 8b6c27c3d0923a838ef33cc7db8f0b4006470c40..e93de35ae83c43f63be7a71a65e663db1721e28d 100644 (file)
@@ -95,14 +95,12 @@ namespace OpenCASCADE
   /**
    * Count the subobjects of a shape. This function is useful to
    * gather information about the TopoDS_Shape passed as argument. It
-   * counts the number of faces, edges and vertices (the only
+   * returns the number of faces, edges and vertices (the only
    * topological entities associated with actual geometries) which are
    * contained in the given shape.
    */
-  void count_elements(const TopoDS_Shape &shape,
-                      unsigned int &n_faces,
-                      unsigned int &n_edges,
-                      unsigned int &n_vertices);
+  std_cxx11::tuple<unsigned int, unsigned int, unsigned int>
+  count_elements(const TopoDS_Shape &shape);
 
   /**
    * Read IGES files and translate their content into openCascade
@@ -123,16 +121,17 @@ namespace OpenCASCADE
 
   /**
     * This function returns the tolerance associated with the shape.
-    * Each CAD geometrical object is defined along with a tolerance, which indicates
-    * possible inaccuracy of its placement. For instance, the tolerance tol of a vertex
-    * indicates that it can be located in any point contained in a sphere centered
-    * in the nominal position and having radius tol. While carrying out an operation
-    * such as projecting a point onto a
-    * surface (which will in turn have its tolerance) we must keep in mind that the
-    * precision of the projection will be limited by the tolerance with which the
-    * surface is built.
-    * The tolerance is computed taking the maximum tolerance among the subshapes
-    * composing the shape.
+    * Each CAD geometrical object is defined along with a tolerance,
+    * which indicates possible inaccuracy of its placement. For
+    * instance, the tolerance of a vertex indicates that it can be
+    * located in any point contained in a sphere centered in the
+    * nominal position and having radius tol. While carrying out an
+    * operation such as projecting a point onto a surface (which will
+    * in turn have its tolerance) we must keep in mind that the
+    * precision of the projection will be limited by the tolerance
+    * with which the surface is built.  The tolerance is computed
+    * taking the maximum tolerance among the subshapes composing the
+    * shape.
     */
   double get_shape_tolerance(const TopoDS_Shape &shape);
 
@@ -219,38 +218,75 @@ namespace OpenCASCADE
                                std::vector<TopoDS_Shell> &shells,
                                std::vector<TopoDS_Wire> &wires);
 
+  /**
+   * Project the point @p origin on the topological shape given by @p
+   * in_shape, and returns the projected point, the subshape which
+   * contains the point and the parametric u and v coordinates of the
+   * point within the resulting shape. If the shape is not elementary,
+   * all its subshapes are iterated, faces first, then edges, and the
+   * returned shape is the closest one to the point @p origin. If the
+   * returned shape is an edge, then only the u coordinate is filled
+   * with sensible information, and the v coordinate is set to zero.
+   *
+   * This function returns a tuple containing the projected point, the
+   * shape, the u coordinate and the v coordinate (which is different
+   * from zero only if the resulting shape is a face).
+   */
+  std_cxx11::tuple<Point<3>, TopoDS_Shape, double, double>
+  project_point_and_pull_back(const TopoDS_Shape &in_shape,
+                              const Point<3> &origin,
+                              const double tolerance=1e-7);
 
   /**
-   * Get the closest point to the given topological shape. If the
-   * shape is not elementary, all its subshapes are iterated, faces
-   * first, then edges, and the closest point is returned together
-   * with the shape which contains it, and the u v coordinates of the
-   * point. If the returned shape is an edge, then only the u
-   * coordinate is filled with sensible information, and the v
-   * coordinate is set to zero.
+   * Return the projection of the point @p origin on the topological
+   * shape given by @p in_shape. If the shape is not elementary, all
+   * its subshapes 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,
-                         TopoDS_Shape &out_shape,
-                         double &u,
-                         double &v,
+  Point<3> closest_point(const TopoDS_Shape &in_shape,
+                         const Point<3> &origin,
                          const double tolerance=1e-7);
 
+  /**
+   * Given an elementary shape @p in_shape and the reference
+   * coordinates within the shape, returns the corresponding point in
+   * real space. If the shape is a TopoDS_Edge, the @p v coordinate is
+   * ignored. Only edges or faces, as returned by the function
+   * project_point_and_pull_back(), can be 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);
+
 
   /**
-   * Get the closest point to the given topological shape. If the
-   * shape is not elementary, all its subshapes are iterated, faces
-   * first, then edges, and the closest point is returned together
-   * with the shape which contains it, and the u v coordinates of the
-   * point. If the returned shape is an edge, then only the u
-   * coordinate is filled with sensible information, and the v
-   * coordinate is set to zero.
+   * Given a TopoDS_Face @p face and the reference coordinates within
+   * this face, returns the corresponding point in real space, the
+   * normal to the surface at that point and the mean curvature as a
+   * tuple.
    */
-  Point<3>  closest_point_and_differential_forms(const TopoDS_Shape in_shape,
-                                                 const Point<3> origin,
-                                                 Point<3> &surface_normal,
-                                                 double &mean_curvature,
-                                                 const double tolerance=1e-7);
+  std_cxx11::tuple<Point<3>, Point<3>, double >
+  push_forward_and_differential_forms(const TopoDS_Face &face,
+                                      const double u,
+                                      const double v,
+                                      const double tolerance=1e-7);
+
+
+  /**
+   * Get the closest point to the given topological shape, together
+   * with the normal and the mean curvature at that point. If the
+   * shape is not elementary, all its sub-faces (only the faces) are
+   * iterated, faces first, and only the closest point is
+   * returned. This function will throw an exception if the @p
+   * in_shape does not contain at least one face.
+   */
+  std_cxx11::tuple<Point<3>, Point<3>, double>
+  closest_point_and_differential_forms(const TopoDS_Shape &in_shape,
+                                       const Point<3> &origin,
+                                       const double tolerance=1e-7);
+
 
   /**
    * Intersect a line passing through the given @p origin point along
@@ -323,7 +359,6 @@ namespace OpenCASCADE
    * Trying to make operations on the wrong type of shapes.
    */
   DeclException0(ExcUnsupportedShape);
-
 }
 
 
index 5a1864753092bcba29d84a51d804486d5b968017..236c4229f5be5336d4beb9da9ca40749b3bc7f32 100644 (file)
@@ -68,15 +68,14 @@ namespace OpenCASCADE
   project_to_manifold (const std::vector<Point<spacedim> > &surrounding_points,
                        const Point<spacedim> &candidate) const
   {
-    TopoDS_Shape out_shape;
-    double u=0, v=0;
+#ifdef DEBUG
     for (unsigned int i=0; i<surrounding_points.size(); ++i)
-      Assert(closest_point(sh, surrounding_points[i], out_shape, u, v)
+      Assert(closest_point(sh, surrounding_points[i], tolerance)
              .distance(surrounding_points[i]) <
              std::max(tolerance*surrounding_points[i].norm(), tolerance),
              ExcPointNotOnManifold(surrounding_points[i]));
-
-    return closest_point(sh, candidate, out_shape, u, v);
+#endif
+    return closest_point(sh, candidate,tolerance);
   }
 
 
@@ -98,16 +97,14 @@ namespace OpenCASCADE
   project_to_manifold (const std::vector<Point<spacedim> > &surrounding_points,
                        const Point<spacedim> &candidate) const
   {
-    TopoDS_Shape out_shape;
-    double u=0, v=0;
+#ifdef DEBUG
     for (unsigned int i=0; i<surrounding_points.size(); ++i)
-      Assert(closest_point(sh, surrounding_points[i], out_shape, u, v)
+      Assert(closest_point(sh, surrounding_points[i],tolerance)
              .distance(surrounding_points[i]) <
              std::max(tolerance*surrounding_points[i].norm(), tolerance),
              ExcPointNotOnManifold(surrounding_points[i]));
-
+#endif
     return line_intersection(sh, candidate, direction, tolerance);
-
   }
 
 
@@ -121,15 +118,10 @@ namespace OpenCASCADE
   {
     Assert(spacedim == 3, ExcNotImplemented());
 
-    unsigned int n_faces;
-    unsigned int n_edges;
-    unsigned int n_vertices;
-    count_elements(sh,
-                   n_faces,
-                   n_edges,
-                   n_vertices);
+    std_cxx11::tuple<unsigned int, unsigned int, unsigned int>
+    counts = count_elements(sh);
 
-    Assert(n_faces > 0, ExcMessage("NormalToMeshProjectionBoundary needs a shape containing faces to operate."));
+    Assert(std_cxx11::get<0>(counts) > 0, ExcMessage("NormalToMeshProjectionBoundary needs a shape containing faces to operate."));
   }
 
 
@@ -141,13 +133,15 @@ namespace OpenCASCADE
     TopoDS_Shape out_shape;
     double u=0, v=0;
     Point<3> average_normal(0.0,0.0,0.0);
+#ifdef DEBUG
     for (unsigned int i=0; i<surrounding_points.size(); ++i)
       {
-        Assert(closest_point(sh, surrounding_points[i], out_shape, u, v)
+        Assert(closest_point(sh, surrounding_points[i], tolerance)
                .distance(surrounding_points[i]) <
                std::max(tolerance*surrounding_points[i].norm(), tolerance),
                ExcPointNotOnManifold(surrounding_points[i]));
       }
+#endif
 
     switch (surrounding_points.size())
       {
@@ -155,10 +149,12 @@ namespace OpenCASCADE
       {
         for (unsigned int i=0; i<surrounding_points.size(); ++i)
           {
-            Point<3> surface_normal;
-            double mean_curvature;
-            closest_point_and_differential_forms(sh, surrounding_points[i], surface_normal, mean_curvature);
-            average_normal += surface_normal;
+            std_cxx11::tuple<Point<3>, Point<3>, double>
+            p_and_diff_forms =
+              closest_point_and_differential_forms(sh,
+                                                   surrounding_points[i],
+                                                   tolerance);
+            average_normal += std_cxx11::get<1>(p_and_diff_forms);
           }
 
         average_normal/=2.0;
@@ -166,11 +162,10 @@ namespace OpenCASCADE
         Assert(average_normal.norm() > 1e-4,
                ExcMessage("Failed to refine cell: the average of the surface normals at the surrounding edge turns out to be a null vector, making the projection direction undetermined."));
 
-        Point<3> P = (surrounding_points[0]+surrounding_points[1])/2;
-        Point<3> N = surrounding_points[0]-surrounding_points[1];
-        N = N/sqrt(N.square());
-        average_normal = average_normal-(average_normal*N)*N;
-        average_normal = average_normal/average_normal.norm();
+        Point<3> T = surrounding_points[0]-surrounding_points[1];
+        T /= T.norm();
+        average_normal = average_normal-(average_normal*T)*T;
+        average_normal /= average_normal.norm();
         break;
       }
       case 8:
@@ -194,10 +189,10 @@ namespace OpenCASCADE
 
         average_normal = (n1+n2+n3+n4)/4.0;
 
-        Assert(average_normal.norm() > 1e-4,
+        Assert(average_normal.norm() > tolerance,
                ExcMessage("Failed to refine cell: the average of the surface normals at the surrounding edge turns out to be a null vector, making the projection direction undetermined."));
 
-        average_normal = average_normal/average_normal.norm();
+        average_normal /= average_normal.norm();
         break;
       }
       default:
index c44b63ade471976c9aa9b5243869ce722a8fd791..51eb8d7ab6b96b5cb3520adb62ce42f4f3928f1a 100644 (file)
@@ -62,13 +62,11 @@ DEAL_II_NAMESPACE_OPEN
 
 namespace OpenCASCADE
 {
-  void count_elements(const TopoDS_Shape &shape,
-                      unsigned int &n_faces,
-                      unsigned int &n_edges,
-                      unsigned int &n_vertices)
+  std_cxx11::tuple<unsigned int, unsigned int, unsigned int>
+  count_elements(const TopoDS_Shape &shape)
   {
     TopExp_Explorer exp;
-    n_faces=0, n_edges=0, n_vertices=0;
+    unsigned int n_faces=0, n_edges=0, n_vertices=0;
     for (exp.Init(shape, TopAbs_FACE);
          exp.More(); exp.Next(), ++n_faces)
       {}
@@ -78,6 +76,7 @@ namespace OpenCASCADE
     for (exp.Init(shape, TopAbs_VERTEX);
          exp.More(); exp.Next(), ++n_vertices)
       {}
+    return std_cxx11::tuple<unsigned int, unsigned int, unsigned int>(n_faces, n_edges, n_vertices);
   }
 
   void extract_geometrical_shapes(const TopoDS_Shape &shape,
@@ -347,7 +346,6 @@ namespace OpenCASCADE
 
     double minDistance = 1e7;
     double distance;
-    int lowest_dist_int = 0;
     Point<3> result;
     for (int i=0; i<Inters.NbPnt(); ++i)
       {
@@ -357,7 +355,6 @@ namespace OpenCASCADE
           {
             minDistance = distance;
             result = point(Inters.Pnt(i+1));
-            lowest_dist_int = i+1;
           }
       }
 
@@ -395,15 +392,11 @@ namespace OpenCASCADE
     return out_shape;
   }
 
-
-  Point<3> closest_point(const TopoDS_Shape in_shape,
-                         const Point<3> origin,
-                         TopoDS_Shape &out_shape,
-                         double &u,
-                         double &v,
-                         const double tolerance)
+  std_cxx11::tuple<Point<3>, TopoDS_Shape, double, double>
+  project_point_and_pull_back(const TopoDS_Shape &in_shape,
+                              const Point<3> &origin,
+                              const double tolerance)
   {
-
     TopExp_Explorer exp;
     gp_Pnt Pproj = point(origin);
 
@@ -412,8 +405,10 @@ namespace OpenCASCADE
 
     unsigned int counter = 0;
     unsigned int face_counter = 0;
-    u=0;
-    v=0;
+
+    TopoDS_Shape out_shape;
+    double u=0;
+    double v=0;
 
     for (exp.Init(in_shape, TopAbs_FACE); exp.More(); exp.Next())
       {
@@ -473,53 +468,89 @@ namespace OpenCASCADE
         }
 
     Assert(counter > 0, ExcMessage("Could not find projection points."));
-    return point(Pproj);
+    return std_cxx11::tuple<Point<3>, TopoDS_Shape, double, double>
+           (point(Pproj),out_shape, u, v);
   }
 
 
-  Point<3> closest_point_and_differential_forms(const TopoDS_Shape in_shape,
-                                                const Point<3> origin,
-                                                Point<3> &surface_normal,
-                                                double &mean_curvature,
-                                                const double tolerance)
-
+  Point<3> closest_point(const TopoDS_Shape &in_shape,
+                         const Point<3> &origin,
+                         const double tolerance)
   {
+    std_cxx11::tuple<Point<3>, TopoDS_Shape, double, double>
+    ref = project_point_and_pull_back(in_shape, origin, tolerance);
+    return std_cxx11::get<0>(ref);
+  }
 
-    TopoDS_Shape out_shape;
-    double u, v;
+  std_cxx11::tuple<Point<3>, Point<3>, double>
+  closest_point_and_differential_forms(const TopoDS_Shape &in_shape,
+                                       const Point<3> &origin,
+                                       const double tolerance)
+
+  {
+    std_cxx11::tuple<Point<3>, TopoDS_Shape, double, double>
+    shape_and_params = project_point_and_pull_back(in_shape,
+                                                   origin,
+                                                   tolerance);
 
-    Point<3> result = closest_point(in_shape,
-                                    origin,
-                                    out_shape,
-                                    u,
-                                    v,
-                                    tolerance);
+    TopoDS_Shape &out_shape = std_cxx11::get<1>(shape_and_params);
+    double &u = std_cxx11::get<2>(shape_and_params);
+    double &v = std_cxx11::get<3>(shape_and_params);
 
     // just a check here: the number of faces in out_shape must be 1, otherwise
     // something is wrong
-    unsigned int n_faces, n_edges, n_vertices;
-    count_elements(out_shape,
-                   n_faces,
-                   n_edges,
-                   n_vertices);
+    std_cxx11::tuple<unsigned int, unsigned int, unsigned int> numbers =
+      count_elements(out_shape);
 
-    Assert(n_faces > 0, ExcMessage("Could not find normal: the shape containing the closest point has 0 faces."));
-    Assert(n_faces < 2, ExcMessage("Could not find normal: the shape containing the closest point has more than 1 face."));
+    Assert(std_cxx11::get<0>(numbers) > 0,
+           ExcMessage("Could not find normal: the shape containing the closest point has 0 faces."));
+    Assert(std_cxx11::get<0>(numbers) < 2,
+           ExcMessage("Could not find normal: the shape containing the closest point has more than 1 face."));
 
 
     TopExp_Explorer exp;
-    exp.Init(in_shape, TopAbs_FACE);
+    exp.Init(out_shape, TopAbs_FACE);
     TopoDS_Face face = TopoDS::Face(exp.Current());
+    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)
+  {
+    switch (in_shape.ShapeType())
+      {
+      case TopAbs_FACE:
+      {
+        BRepAdaptor_Surface surf(TopoDS::Face(in_shape));
+        return point(surf.Value(u,v));
+      }
+      case TopAbs_EDGE:
+      {
+        BRepAdaptor_Curve curve(TopoDS::Edge(in_shape));
+        return point(curve.Value(u));
+      }
+      default:
+        Assert(false, ExcUnsupportedShape());
+      }
+    return Point<3>();
+  }
+
+  std_cxx11::tuple<Point<3>, Point<3>, double >
+  push_forward_and_differential_forms(const TopoDS_Face &face,
+                                      const double u,
+                                      const double v,
+                                      const double tolerance)
+  {
     Handle(Geom_Surface) SurfToProj = BRep_Tool::Surface(face);
-    GeomLProp_SLProps props(SurfToProj, u, v, 1, tolerance);
+    GeomLProp_SLProps props(SurfToProj, u, v, 1, 1e-7);
+    gp_Pnt Value = props.Value();
+    Assert(props.IsNormalDefined(), ExcMessage("Normal is not well defined!"));
     gp_Dir Normal = props.Normal();
+    Assert(props.IsCurvatureDefined(), ExcMessage("Curvature is not well defined!"));
     Standard_Real Mean_Curvature = props.MeanCurvature();
-
-    surface_normal = Point<3>(Normal.X(),Normal.Y(),Normal.Z());
-    mean_curvature = double(Mean_Curvature);
-
-    return result;
+    Point<3> normal = Point<3>(Normal.X(),Normal.Y(),Normal.Z());
+    return std_cxx11::tuple<Point<3>, Point<3>, double>(point(Value), normal, Mean_Curvature);
   }
 
 
index c648b19aef47021103702240250ca51a27f6b320..f8162944789d828e9f6a31659043b301cf90dc14 100644 (file)
@@ -70,7 +70,13 @@ int main ()
   TopoDS_Shape sh;
   for(unsigned int i=0; i<points.size(); ++i) 
     {
-      Point<3> pp = closest_point(edge, points[i], sh, u, v);
+      std_cxx11::tuple<Point<3>, TopoDS_Shape, double, double> ref = 
+       project_point_and_pull_back(edge, points[i]);
+      
+      Point<3> pp = std_cxx11::get<0>(ref);
+      sh = std_cxx11::get<1>(ref);
+      u = std_cxx11::get<2>(ref);
+      v = std_cxx11::get<3>(ref);
       
       deallog << "Origin: " << points[i]
              << ", on unit circle: " << pp
index 226f116c46fc4e9f6851420f574f8f7e4339b1f1..907c270f71eacb2b59e3bbbb4312906f9289a2ce 100644 (file)
@@ -9,7 +9,7 @@
 //
 //-----------------------------------------------------------
 
-// Read the first face of goteborg.iges, attach it to a projector,
+// Read the file side.iges, attach it to a projector,
 // create a single-cell Triangulation, and refine it with the new
 // projector.
 
@@ -44,7 +44,7 @@ int main ()
 {
   std::ofstream logfile("output");
   
-  TopoDS_Shape sh = read_IGES(SOURCE_DIR "/iges_files/goteborg.iges", 1);
+  TopoDS_Shape sh = read_IGES(SOURCE_DIR "/iges_files/side.iges");
   std::vector<TopoDS_Face> faces;
   std::vector<TopoDS_Edge> edges;
   std::vector<TopoDS_Vertex> vertices;
index 89325f344a423ac4adffad2187ecd1a5e7ad7268..7b880ba9927affeb66002ce91bf23893eb1b7891 100644 (file)
@@ -1,30 +1,30 @@
 $NOD
 25
-1  -1986.99 -2.3276e-17 -149.998
-2  -2637.58 3.55271e-15 -249.603
-3  -2008.61 314.696 428.208
-4  -3059.55 1.42109e-14 501.418
-5  -2304.65 -7.10543e-15 -245.601
-6  -2000.53 167.884 134.23
-7  -2725.34 0.0425119 204.053
-8  -2545.92 193.951 440.032
-9  -2404.02 65.5321 139.195
-10  -2142.3 -5.17421e-15 -216.101
-11  -2469.32 0 -267.915
-12  -1995.92 96.6019 -14.6159
-13  -2004.03 238.906 282.434
-14  -2632.21 0.000934116 -15.8842
-15  -2877.56 -5.67115e-09 370.719
-16  -2279.15 262.151 429.296
-17  -2808.9 111.45 459.573
-18  -2354.15 31.1791 -52.9914
-19  -2472.24 120.934 295.122
-20  -2202.2 116.44 136.819
-21  -2565.7 37.2635 170.327
-22  -2173.99 60.762 -34.0318
-23  -2511.39 14.632 -39.6776
-24  -2239.11 183.341 286.741
-25  -2682.97 67.745 321.647
+1  0 0 0.001
+2  -2.11667e-14 0.001 0.001
+3  0 0 0
+4  -2.11667e-14 0.001 0
+5  -0.000234467 0.000444022 0.001
+6  3.40745e-16 -4.38101e-16 0.0005
+7  -2.11537e-14 0.001 0.0005
+8  -0.000234467 0.000444022 0
+9  -0.000234467 0.000444022 0.0005
+10  -0.000139999 0.00020943 0.001
+11  -0.00021601 0.000764021 0.001
+12  5.11117e-16 -6.57151e-16 0.00075
+13  5.11117e-16 -6.57151e-16 0.00025
+14  -2.11472e-14 0.001 0.00075
+15  -2.11472e-14 0.001 0.00025
+16  -0.000139999 0.00020943 0
+17  -0.00021601 0.000764021 0
+18  -0.000234467 0.000444022 0.00075
+19  -0.000234467 0.000444022 0.00025
+20  -0.000139999 0.00020943 0.0005
+21  -0.00021601 0.000764021 0.0005
+22  -0.000139999 0.00020943 0.00075
+23  -0.00021601 0.000764021 0.00075
+24  -0.000139999 0.00020943 0.00025
+25  -0.00021601 0.000764021 0.00025
 $ENDNOD
 $ELM
 16
index 3a4633c450634316e6a079529d79a169a3927d67..f2f181d878ae36c0f8ce11fa6fed77054c397bdd 100644 (file)
@@ -36,11 +36,15 @@ int main ()
 
   extract_geometrical_shapes(sh, faces, edges, vertices);
   
-  unsigned int nf=0, ne=0, nv=0;
-  count_elements(sh, nf, ne, nv);
-  deallog << "Shape contains " << nf << " faces, "
-         << ne << " edges, and " 
-         << nv << " vertices." << std::endl;
+  std_cxx11::tuple<unsigned int, unsigned int, unsigned int>
+    n = count_elements(sh);
+  unsigned int nf=std_cxx11::get<0>(n);
+  unsigned int ne=std_cxx11::get<1>(n);
+  unsigned int nv=std_cxx11::get<2>(n);
+
+  deallog << "Shape contains " << std_cxx11::get<0>(n) << " faces, "
+         << std_cxx11::get<1>(n) << " edges, and " 
+         << std_cxx11::get<2>(n) << " vertices." << std::endl;
   
   if(nf != faces.size())
     deallog << "Error!" << std::endl;
index 635a1439e259ab43c9c7edb6fc2ad1fc3ebf9683..bdc5e73e04906602786d511cded8ea677b1b061a 100644 (file)
@@ -44,7 +44,7 @@ int main ()
   std::ofstream logfile("output");
 
   gp_Pnt center(.5,.5,.5);
-  Standard_Real radius(Point<3>().distance(Pnt(center)));
+  Standard_Real radius(Point<3>().distance(point(center)));
 
   
   TopoDS_Face face = BRepPrimAPI_MakeSphere(center, radius);
index e215376ba9a7b1ccd04dd66eeaf9f80b78db225b..3aa521047797fc869bb7eef6d9890f36b001bf5d 100644 (file)
@@ -49,7 +49,7 @@ int main()
 
   // Create a sphere
   gp_Pnt center(.5,.5,.5);
-  Standard_Real radius(Point<3>().distance(Pnt(center)));
+  Standard_Real radius(Point<3>().distance(point(center)));
 
   TopoDS_Face face = BRepPrimAPI_MakeSphere(center, radius);
   

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.