]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Added Axis Projection.
authorLuca Heltai <luca.heltai@sissa.it>
Tue, 12 Aug 2014 15:58:41 +0000 (17:58 +0200)
committerLuca Heltai <luca.heltai@sissa.it>
Fri, 31 Oct 2014 10:06:31 +0000 (11:06 +0100)
include/deal.II/grid/occ_boundary_lib.h
include/deal.II/grid/occ_utilities.h
source/grid/occ_boundary_lib.cc
source/grid/occ_boundary_lib.inst.in
source/grid/occ_utilities.cc

index e46be218febe4f4b6533a733cc12448c902cb8aa..eb53542afafb4660dc72334a9077156ff5906cc7 100644 (file)
@@ -73,13 +73,58 @@ namespace OpenCASCADE
     project_to_manifold (const std::vector<Point<spacedim> > &surrounding_points,
                         const Point<spacedim> &candidate) const;
 
+
+  private:
+    /**
+     * The topological shape which is used internally to project
+     * points. You can construct one such a shape by calling the
+     * OpenCASCADE::read_IGES() function, which will create a
+     * TopoDS_Shape with the geometry contained in the IGES file.
+     */ 
+    const TopoDS_Shape sh;
+
     /**
-     * Exception thrown when the point specified as argument does not
-     * lie between #tolerance from the given TopoDS_Shape.
+     * Relative tolerance used by this class to compute distances.
      */
-    DeclException1 (ExcPointNotOnManifold,
-                   Point<spacedim>,
-                   <<"The point [ "<<arg1<<" ] is not on the manifold.");
+    const double tolerance;
+  };
+
+  /**
+   * A Boundary object based on OpenCASCADE TopoDS_Shape where new
+   * points are first computed using the FlatManifold class, and then
+   * projected along the direction given at construction time, using
+   * OpenCASCADE utilities.
+   *
+   * This class makes no assumptions on the shape you pass to it, and
+   * the topological dimension of the Manifold is inferred from the
+   * TopoDS_Shape itself. 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.
+   * 
+   * @author Luca Heltai, Andrea Mola, 2011--2014.
+   */
+  template <int dim, int spacedim>
+  class AxisProjectionBoundary : public Boundary<dim,spacedim> {
+    public:
+    AxisProjectionBoundary(const TopoDS_Shape sh, 
+                          const Point<3> direction, 
+                          const double tolerance=1e-7);
+    
+    /**
+     * Perform the actual projection onto the manifold. This function,
+     * in debug mode, checks that each of the #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 normal
+     * projection algorithms.
+     */
+    virtual Point<spacedim>
+    project_to_manifold (const std::vector<Point<spacedim> > &surrounding_points,
+                        const Point<spacedim> &candidate) const;
 
   private:
     /**
@@ -90,6 +135,11 @@ namespace OpenCASCADE
      */ 
     const TopoDS_Shape sh;
 
+    /**
+     * Direction used to project new points on the shape. 
+     */
+    const Point<3> direction;
+
     /**
      * Relative tolerance used by this class to compute distances.
      */
index 45b1d52b0f2ba289d176b5cf183b465b55762b34..d483b58e11bdec13a5688428fe79de6b8c0d31f0 100644 (file)
@@ -164,6 +164,19 @@ namespace OpenCASCADE
                         double &v, 
                         const double tolerance=1e-7);
 
+  /**
+   * Intersect a line passing through the given #origin point along
+   * #direction and the given topological shape. If there is more than
+   * one intersection, it will return the closest one.
+   *
+   * The optional #tolerance parameter is used to compute distances.
+   */
+  Point<3> axis_intersection(const TopoDS_Shape in_shape, 
+                            const Point<3> origin, 
+                            const Point<3> direction, 
+                            const double tolerance=1e-7);
+  
+
   /**
    * Convert OpenCASCADE point into a Point<3>.
    */
@@ -185,8 +198,34 @@ namespace OpenCASCADE
   inline bool point_compare(const dealii::Point<3> &p1, const dealii::Point<3> &p2,
                            const dealii::Point<3> direction=Point<3>(),
                            const double tolerance=1e-10);
-                                     
-}
+
+
+  /**
+   * Exception thrown when the point specified as argument does not
+   * lie between #tolerance from the given TopoDS_Shape.
+   */
+  DeclException1 (ExcPointNotOnManifold,
+                 Point<3>,
+                 <<"The point [ "<<arg1<<" ] is not on the manifold.");
+
+  /**
+   * Exception thrown when the point specified as argument cannot be
+   * projected to the manifold.
+   */                        
+  DeclException1 (ExcProjectionFailed, 
+                 Point<3>,
+                 <<"Projection of point [ "<< arg1
+                 << " ] failed.");
+
+  /**
+   * Thrown when internal OpenCASCADE utilities fail to return the OK
+   * status.
+   */
+  DeclException1 (ExcOCCError, 
+                 IFSelect_ReturnStatus, 
+                 <<"An OpenCASCADE routine failed with return status "
+                 <<arg1);
+} 
 /*@}*/
 
 DEAL_II_NAMESPACE_CLOSE
index 2ef10faf2e67e051122b72bf55896760125d303b..64724e91480c9f5b28c3800f86bd756dff2f8e05 100644 (file)
@@ -6,6 +6,8 @@ DEAL_II_NAMESPACE_OPEN
 
 namespace OpenCASCADE
 {
+
+  /*============================== NormalProjectionBoundary ==============================*/
   template <int dim, int spacedim>
   NormalProjectionBoundary<dim,spacedim>::NormalProjectionBoundary(const TopoDS_Shape sh, 
                                                                   const double tolerance) :
@@ -27,10 +29,41 @@ namespace OpenCASCADE
             .distance(surrounding_points[i]) < (surrounding_points[i].norm()>0 ? 
                                                 tolerance*surrounding_points[i].norm() :
                                                 tolerance),
-            ExcPointNotOnManifold(surrounding_points[i]));
+            // ExcPointNotOnManifold(surrounding_points[i]));
+            ExcMessage("Points are not on manifold"));
     
     return closest_point(sh, candidate, out_shape, u, v);
   }
+
+  /*============================== AxisProjectionBoundary ==============================*/
+  template <int dim, int spacedim>
+  AxisProjectionBoundary<dim,spacedim>::AxisProjectionBoundary(const TopoDS_Shape sh, 
+                                                              const Point<3> direction,
+                                                              const double tolerance) :
+    sh(sh),
+    direction(direction),
+    tolerance(tolerance) 
+  {
+    Assert(spacedim == 3, ExcImpossibleInDim(spacedim));
+  }
+  
+  
+  template <int dim, int spacedim>
+  Point<spacedim>  AxisProjectionBoundary<dim,spacedim>::
+  project_to_manifold (const std::vector<Point<spacedim> > &surrounding_points,
+                      const Point<spacedim> &candidate) const {
+    TopoDS_Shape out_shape;
+    double u=0, v=0;
+    for(unsigned int i=0; i<surrounding_points.size(); ++i)
+      Assert(closest_point(sh, surrounding_points[i], out_shape, u, v)
+            .distance(surrounding_points[i]) < (surrounding_points[i].norm()>0 ? 
+                                                tolerance*surrounding_points[i].norm() :
+                                                tolerance),
+            ExcPointNotOnManifold(surrounding_points[i]));
+    
+    return axis_intersection(sh, candidate, direction, tolerance);
+  }
+
   
   // Explicit instantiations
 #include "occ_boundary_lib.inst"  
index 6e97bc7391f5246c7ff2b31eefe9a3d7c9bb60d0..3471c0aea80e0a5c1ed1b484f55491622cf14568 100644 (file)
@@ -19,6 +19,7 @@
 for (deal_II_dimension : DIMENSIONS)
   {
        template class NormalProjectionBoundary<deal_II_dimension, 3>;
+       template class AxisProjectionBoundary<deal_II_dimension, 3>;
   }
 
 
index 73b04d3ed129c9c0a2e8b79b3748ccb14cff4997..1bf4571e1b567a5fffb8e80aaca028c7e89387db 100644 (file)
 #include <TColStd_SequenceOfTransient.hxx>
 #include <TColStd_HSequenceOfTransient.hxx>
 #include <TopExp_Explorer.hxx>
-// #include <gp_Pnt.hxx>
-// #include <gp_Vec.hxx>
+#include <gp_Pnt.hxx>
+#include <gp_Lin.hxx>
+#include <gp_Vec.hxx>
 #include <GeomAPI_ProjectPointOnSurf.hxx>
 #include <GeomAPI_ProjectPointOnCurve.hxx>
+#include <IntCurvesFace_ShapeIntersector.hxx>
 // #include <Standard_Real.hxx>
 // #include <Standard_Integer.hxx>
 // #include <BRep_Tool.hxx>
@@ -243,6 +245,29 @@ namespace OpenCASCADE
    }
 
 
+  Point<3> axis_intersection(const TopoDS_Shape in_shape, 
+                            const Point<3> origin, 
+                            const Point<3> direction, 
+                            const double tolerance)
+  {
+    // translating original Point<dim> to gp point
+
+    gp_Pnt P0 = Pnt(origin);
+    gp_Ax1 gpaxis(P0, gp_Dir(direction(0), direction(1), direction(2)));
+    gp_Lin line(gpaxis);
+
+    // destination point
+    gp_Pnt Pproj(0.0,0.0,0.0);
+
+    // we prepare now the surface for the projection we get the whole
+    // shape from the iges model
+    IntCurvesFace_ShapeIntersector Inters;
+    Inters.Load(in_shape,tolerance);
+       
+    Inters.PerformNearest(line,-RealLast(),+RealLast());
+    Assert(Inters.IsDone(), ExcMessage("Could not project point."));
+    return Pnt(Inters.Pnt(1));
+  }
 
   TopoDS_Edge interpolation_curve(std::vector<Point<3> > &curve_points,
                                   const Point<3> direction,

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.