]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Removed double inheritance. Now all Boundaries are derived from FlatManifold.
authorheltai <heltai@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 3 Dec 2013 23:07:36 +0000 (23:07 +0000)
committerheltai <heltai@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 3 Dec 2013 23:07:36 +0000 (23:07 +0000)
git-svn-id: https://svn.dealii.org/branches/branch_manifold_id@31865 0785d39b-7218-0410-832d-ea1e28bc413d

12 files changed:
deal.II/include/deal.II/grid/manifold.h
deal.II/include/deal.II/grid/tria.h
deal.II/include/deal.II/grid/tria_accessor.h
deal.II/include/deal.II/grid/tria_accessor.templates.h
deal.II/include/deal.II/grid/tria_boundary.h
deal.II/include/deal.II/grid/tria_boundary_lib.h
deal.II/include/deal.II/numerics/vector_tools.templates.h
deal.II/source/fe/mapping_c1.cc
deal.II/source/grid/grid_tools.cc
deal.II/source/grid/manifold.cc
deal.II/source/grid/tria.cc
deal.II/source/grid/tria_boundary_lib.cc

index 36b118ccebb3856a25d59bc218da44ea10cbcf30..36e73a0515ffd8835e36cde8a091edac6f9a8e70 100644 (file)
@@ -122,7 +122,30 @@ public:
    * exception of type ExcPureFunctionCalled.
    */
   virtual
-  void project_to_manifold (Point<spacedim> &candidate) const;
+  Point<spacedim> project_to_manifold (const Point<spacedim> candidate) const;
+
+  /**
+    * Given a vector of points, return the normals to the Manifold in
+    * those points. Notice that the Manifold may or may not be a
+    * codimension one manifold. If it is not, this function will throw
+    * an Error. Input arguments must be of the same size. The default
+    * implementation calls for each point the function
+    * normal_vector. Derived classes can overload that function.
+    */
+    virtual
+    void get_normals_at_points(std::vector<Point<spacedim> > &normals,
+                              const std::vector<Point<spacedim> > &points) const;
+
+    /** 
+     * Compute the normal to the manifold at the given point. The
+     * default implementation of this function trhows an
+     * Exception. Notice that this function only makes sense when the
+     * manifold is a codimension one manifold.
+     */
+    virtual
+    Point<spacedim> normal_vector(const Point<spacedim> point) const;
+    
+    
 };
 
 
@@ -148,18 +171,29 @@ public:
   /**
    * Let the new point be the average sum of surrounding vertices.
    *
-   * Refer to the general documentation of this class and the
-   * documentation of the base class for more information.
+   * This particular implementation constructs the weighted average of
+   * the surrounding points, and then calss internally the function
+   * project_to_manifold. The reason why we do it this way, is to
+   * allow lazy programmers to implement only the project_to_manifold
+   * function for their own Manifold classes which are small (or
+   * trivial) perturbations of a flat manifold. For most simple
+   * geometries, it is possible to get reasonable results by deriving
+   * your own Manifold class from FlatManifold, and write a new
+   * interface only for the project_to_manifold function.
    */
     virtual Point<spacedim>
     get_new_point(const std::vector<Point<spacedim> > &surrounding_points,
                  const std::vector<double> &weights) const;
 
   /**
-   *  Project to FlatManifold: do nothing.
+   *  Project to FlatManifold: do nothing. Note however that this
+   *  function can be overloaded by derived classes, which will then
+   *  benefit from the logic behind the get_new_point class which are
+   *  often very similar (if not identical) to the one implemented in
+   *  this class.
    */
   virtual
-  void project_to_manifold (Point<spacedim> &candidate) const;
+  Point<spacedim> project_to_manifold (const Point<spacedim> candidate) const;
 };
 
 
index aad2764a5794494230993e2a3c9436c33e482b55..e5dc7dc20e214e98d3946e5991653ca274c26bbe 100644 (file)
@@ -40,8 +40,6 @@
 
 DEAL_II_NAMESPACE_OPEN
 
-template <int dim, int spacedim> class Boundary;
-template <int dim, int spacedim> class StraightBoundary;
 template <int spacedim> class Manifold;
 template <int spacedim> class FlatManifold;
 
@@ -1234,14 +1232,6 @@ private:
   typedef dealii::internal::Triangulation::Iterators<dim, spacedim> IteratorSelector;
 
 public:
-  /**
-   * Default boundary object. This is used
-   * for those boundaries for which no
-   * boundary object has been explicitly
-   * set using set_boundary().
-   */
-  static const StraightBoundary<dim,spacedim> straight_boundary;
-
   /**
    * Default manifold object. This is used
    * for those objects for which no
@@ -1815,7 +1805,7 @@ public:
    * @see @ref GlossBoundaryIndicator "Glossary entry on boundary indicators"
    */
   void set_boundary (const types::manifold_id   number,
-                     const Boundary<dim,spacedim> &boundary_object) DEAL_II_DEPRECATED;
+                     const Manifold<spacedim> &boundary_object) DEAL_II_DEPRECATED;
   /**
    * Assign a Manifold object to a certain part of the
    * triangulation. The Manifold object is used to find the location
@@ -1891,7 +1881,7 @@ public:
    *
    * @see @ref GlossBoundaryIndicator "Glossary entry on boundary indicators"
    */
-  const Boundary<dim,spacedim> &get_boundary (const types::manifold_id number) const;
+  const Manifold<spacedim> &get_boundary (const types::manifold_id number) const;
 
   /**
    * Return a constant reference to a Manifold object used for this
@@ -3214,7 +3204,7 @@ private:
    *  Collection of boundary objects. We store only objects, which are
    *  not of type StraightBoundary.
    */
-  std::map<types::manifold_id, SmartPointer<const Boundary<dim, spacedim> , Triangulation<dim, spacedim> > >  boundary;
+  std::map<types::manifold_id, SmartPointer<const Manifold<spacedim> , Triangulation<dim, spacedim> > >  boundary;
 
     
   /**
index 6c2dcc55b9ae2dabeefac6898e5d2f82f56f7925..ec8971380ac671fb47f68f365be4a9863a6838da 100644 (file)
@@ -36,7 +36,6 @@ template <typename Accessor> class TriaRawIterator;
 template <typename Accessor> class TriaIterator;
 template <typename Accessor> class TriaActiveIterator;
 
-template <int dim, int spacedim> class Boundary;
 template <int spacedim> class Manifold;
 
 
@@ -1161,7 +1160,7 @@ public:
    * Triangulation::get_boundary() function
    * for the boundary object.
    */
-  const Boundary<dim,spacedim> &get_boundary () const;
+  const Manifold<spacedim> &get_boundary () const;
     
   /**
    * Return a constant reference to a
index b0e5a447d4cd9b650edc80c15a53df1eeef30faf..df31c0d49c1ac6f78f93bbb42906e6593db7d2ca 100644 (file)
@@ -27,8 +27,6 @@
 #include <deal.II/grid/tria_iterator.h>
 #include <deal.II/grid/tria_accessor.h>
 #include <deal.II/grid/tria_iterator.templates.h>
-#include <deal.II/grid/manifold.h>
-#include <deal.II/grid/tria_boundary.h>
 #include <deal.II/distributed/tria.h>
 
 #include <cmath>
@@ -1938,7 +1936,7 @@ TriaAccessor<structdim, dim, spacedim>::at_boundary () const
 
 
 template <int structdim, int dim, int spacedim>
-const Boundary<dim,spacedim> &
+const Manifold<spacedim> &
 TriaAccessor<structdim, dim, spacedim>::get_boundary () const
 {
   Assert (structdim<spacedim, ExcImpossibleInDim(dim));
index c7894c5bc7638d93281001e7ffd183495b90375a..e2a58c8a8e02007cce74106779c2057e12e75ca1 100644 (file)
@@ -82,7 +82,7 @@ template <int dim, int space_dim> class Triangulation;
  * @author Wolfgang Bangerth, 1999, 2001, 2009, Ralf Hartmann, 2001, 2008
  */
 template <int dim, int spacedim=dim>
-  class Boundary : public Manifold<spacedim>
+  class Boundary : public FlatManifold<spacedim>
 {
 public:
 
index 2f683b964ed3ac41c39c16db5429d05fc9ee9e06..a6f63a3ebb643902804dd8cceb983805b737dedb 100644 (file)
@@ -83,6 +83,15 @@ public:
   virtual Point<spacedim>
   get_new_point_on_line (const typename Triangulation<dim,spacedim>::line_iterator &line) const;
 
+  /**
+   * Refer to the general documentation of
+   * this class and the documentation of the
+   * base class.
+   */
+  virtual Point<spacedim>
+  get_new_point (const std::vector<Point<spacedim> > &surrounding_points,
+                const std::vector<double> &weights) const;
+
   /**
    * Refer to the general documentation of
    * this class and the documentation of the
@@ -260,6 +269,16 @@ public:
   Point<dim>
   get_new_point_on_line (const typename Triangulation<dim>::line_iterator &line) const;
 
+  /**
+   * Refer to the general documentation of
+   * this class and the documentation of the
+   * base class.
+   */
+  virtual Point<dim>
+  get_new_point (const std::vector<Point<dim> > &surrounding_points,
+                const std::vector<double> &weights) const;
+
+
   /**
    * Refer to the general
    * documentation of this class
@@ -394,6 +413,16 @@ public:
   Point<spacedim>
   get_new_point_on_line (const typename Triangulation<dim,spacedim>::line_iterator &line) const;
 
+  /**
+   * Refer to the general documentation of
+   * this class and the documentation of the
+   * base class.
+   */
+  virtual Point<spacedim>
+  get_new_point (const std::vector<Point<spacedim> > &surrounding_points,
+                const std::vector<double> &weights) const;
+
+
   /**
    * Refer to the general documentation of
    * this class and the documentation of the
@@ -570,6 +599,16 @@ public:
   virtual Point<dim>
   get_new_point_on_quad (const typename Triangulation<dim>::quad_iterator &quad) const;
 
+  /**
+   * Refer to the general documentation of
+   * this class and the documentation of the
+   * base class.
+   */
+  virtual Point<dim>
+  get_new_point (const std::vector<Point<dim> > &surrounding_points,
+                const std::vector<double> &weights) const;
+
+    
   /**
    * Refer to the general
    * documentation of this class
@@ -680,6 +719,15 @@ public:
    */
   virtual Point<dim>
   get_new_point_on_line (const typename Triangulation<dim>::line_iterator &line) const;
+  /**
+   * Refer to the general documentation of
+   * this class and the documentation of the
+   * base class.
+   */
+  virtual Point<dim>
+  get_new_point (const std::vector<Point<dim> > &surrounding_points,
+                const std::vector<double> &weights) const;
+
 
   /**
    * Construct a new point on a quad.
@@ -767,6 +815,16 @@ public:
   virtual Point<spacedim>
   get_new_point_on_quad (const typename Triangulation<dim,spacedim>::quad_iterator &quad) const;
 
+  /**
+   * Refer to the general documentation of
+   * this class and the documentation of the
+   * base class.
+   */
+  virtual Point<spacedim>
+  get_new_point (const std::vector<Point<spacedim> > &surrounding_points,
+                const std::vector<double> &weights) const;
+
+
   /**
    * Construct a new points on a line.
    */
index 77cc717739a52e512e661f0f3a1a74888d1f9569..e814f047db2e0ce8d870620e72ce74375e1a732b 100644 (file)
@@ -4373,8 +4373,7 @@ namespace VectorTools
                    // case in tests/deal.II/no_flux_11.
                     Point<dim> normal_vector
                       = (cell->face(face_no)->get_boundary()
-                         .normal_vector (cell->face(face_no),
-                                         fe_values.quadrature_point(i)));
+                         .normal_vector (fe_values.quadrature_point(i)));
                    if (normal_vector * fe_values.normal_vector(i) < 0)
                      normal_vector *= -1;
                     Assert (std::fabs(normal_vector.norm() - 1) < 1e-14,
index 5a6115220a5fadf5960ce6d4fa60c5c5422872e6..caa97ba5869438cde8cf34a708714017c7df7c90 100644 (file)
@@ -64,11 +64,14 @@ MappingC1<2>::add_line_support_points (const Triangulation<2>::cell_iterator &ce
         {
           // first get the normal vectors at the two vertices of this line
           // from the boundary description
-          const Boundary<dim> &boundary
-            = line->get_triangulation().get_boundary(line->boundary_indicator());
-
-          Boundary<dim>::FaceVertexNormals face_vertex_normals;
-          boundary.get_normals_at_vertices (line, face_vertex_normals);
+          const Manifold<dim> &boundary
+            = line->get_boundary();
+
+         std::vector<Point<dim> > points(2);
+         std::vector<Point<dim> > face_vertex_normals(2);
+         points[0] = line->vertex(0);
+         points[1] = line->vertex(1);
+          boundary.get_normals_at_points (face_vertex_normals, points);
 
           // then transform them into interpolation points for a cubic
           // polynomial
index 669dccb71ba7d2e12d19b20dae822c151ea39211..7cd43ae0fb54ff9627f9d260f14461a04eb693f3 100644 (file)
@@ -1728,8 +1728,7 @@ next_cell:
       fix_up_object (const Iterator &object,
                      const bool respect_manifold)
       {
-        const Boundary<Iterator::AccessorType::dimension,
-              Iterator::AccessorType::space_dimension>
+        const Manifold<Iterator::AccessorType::space_dimension>
               *manifold = (respect_manifold ?
                            &object->get_boundary() :
                            0);
@@ -1799,16 +1798,14 @@ next_cell:
                        /
                        eps);
                 else
-                  gradient[d]
-                    = ((objective_function (object,
-                                            manifold->project_to_surface(object,
-                                                                         object_mid_point + h))
-                        -
-                        objective_function (object,
-                                            manifold->project_to_surface(object,
-                                                                         object_mid_point - h)))
-                       /
-                       eps);
+                   gradient[d]
+                     = ((objective_function (object,
+                                             manifold->project_to_manifold(object_mid_point + h))
+                         -
+                         objective_function (object,
+                                             manifold->project_to_manifold(object_mid_point - h)))
+                        /
+                        eps);
               }
 
             // sometimes, the
@@ -1840,8 +1837,7 @@ next_cell:
                                 gradient;
 
             if (respect_manifold == true)
-              object_mid_point = manifold->project_to_surface(object,
-                                                              object_mid_point);
+              object_mid_point = manifold->project_to_manifold(object_mid_point);
 
             // compute current value of the
             // objective function
index d5d82db9d52efb448d2bea8bcc8f4018a696d5d8..7a6de0639f7225c25484e1ef8d4307a299ae98ff 100644 (file)
@@ -32,11 +32,30 @@ Manifold<spacedim>::~Manifold ()
 
 
 
+template <int spacedim>
+Point<spacedim>
+Manifold<spacedim>::project_to_manifold (const Point<spacedim> candidate) const 
+{
+  Assert (false, ExcPureFunctionCalled());
+  return candidate;
+}
+
 template <int spacedim>
 void
-Manifold<spacedim>::project_to_manifold (Point<spacedim> &candidate) const 
+Manifold<spacedim>::get_normals_at_points(std::vector<Point<spacedim> > &normals,
+                                         const std::vector<Point<spacedim> > &points) const
+{
+  AssertDimension(normals.size(), points.size());
+  for(unsigned int i=0; i<normals.size(); ++i)
+    normals[i] = normal_vector(points[i]);
+}
+
+template <int spacedim>
+Point<spacedim>
+Manifold<spacedim>::normal_vector(const Point<spacedim> point) const
 {
   Assert (false, ExcPureFunctionCalled());
+  return point;
 }
 
 template <int spacedim>
@@ -108,8 +127,7 @@ FlatManifold<spacedim>::
 get_new_point (const std::vector<Point<spacedim> > &surrounding_points,
               const std::vector<double> &weights) const
 {
-  Assert(surrounding_points.size() == weights.size(),
-        ExcDimensionMismatch(surrounding_points.size(), weights.size()));
+  AssertDimension(surrounding_points.size(), weights.size());
 
 #ifdef DEBUG
   double sum=0;
@@ -122,13 +140,15 @@ get_new_point (const std::vector<Point<spacedim> > &surrounding_points,
   Point<spacedim> p;
   for(unsigned int i=0; i<surrounding_points.size(); ++i)
     p += surrounding_points[i]*weights[i];
-  return p;
+  return project_to_manifold(p);
 }
 
 template <int spacedim>
-void
-FlatManifold<spacedim>::project_to_manifold (Point<spacedim> &candidate) const 
-{}
+Point<spacedim> 
+FlatManifold<spacedim>::project_to_manifold (const Point<spacedim> candidate) const 
+{
+  return candidate;
+}
 
 
 /* -------------------------- ManifoldChart --------------------- */
@@ -144,8 +164,7 @@ ManifoldChart<dim, spacedim>::
 get_new_point (const std::vector<Point<spacedim> > &surrounding_points,
               const std::vector<double> &weights) const
 {
-  Assert(surrounding_points.size() == weights.size(),
-        ExcDimensionMismatch(surrounding_points.size(), weights.size()));
+  AssertDimension(surrounding_points.size(), weights.size());
   chart_points.resize(surrounding_points.size());
   
   for(unsigned int i=0; i<surrounding_points.size(); ++i) 
index 4d944f8ee8e21d1b4bb30f6fba00da6efcfff59e..aaae9a3da47cdb55078ef3e26d28a7155436f74d 100644 (file)
@@ -9136,11 +9136,8 @@ namespace internal
       void
       prevent_distorted_boundary_cells (Triangulation<dim,spacedim> &triangulation)
       {
-        // If the codimension is
-        // one, we cannot perform
-        // this check yet.
-        if (spacedim>dim) return;
-
+       std::vector<double> wp;
+       std::vector<Point<spacedim> > sp;
         for (typename Triangulation<dim,spacedim>::cell_iterator
              cell=triangulation.begin(); cell!=triangulation.end(); ++cell)
           if (cell->at_boundary() &&
@@ -9173,8 +9170,9 @@ namespace internal
                         // the new point on the
                         // boundary would be
                         // this one.
+                       get_default_points_and_weights(sp,wp,face);
                         const Point<spacedim> new_bound
-                          = face->get_boundary().get_new_point_on_face(face);
+                          = face->get_boundary().get_new_point(sp,wp);
 //                        triangulation.get_boundary(face->boundary_indicator())
 //                             .get_new_point_on_face (face);
                         // to check it,
@@ -9545,11 +9543,6 @@ namespace internal
 }
 
 
-template <int dim, int spacedim>
-const StraightBoundary<dim,spacedim>
-Triangulation<dim, spacedim>::straight_boundary = StraightBoundary<dim,spacedim>();
-
-
 template <int dim, int spacedim>
 const FlatManifold<spacedim>
 Triangulation<dim, spacedim>::flat_manifold = FlatManifold<spacedim>();
@@ -9653,10 +9646,11 @@ Triangulation<dim, spacedim>::set_mesh_smoothing(const MeshSmoothing mesh_smooth
 template <int dim, int spacedim>
 void
 Triangulation<dim, spacedim>::set_boundary (const types::manifold_id m_number,
-                                            const Boundary<dim, spacedim> &boundary_object)
+                                            const Manifold<spacedim> &boundary_object)
 {
   types::boundary_id number = static_cast<types::boundary_id>(m_number);
-  Assert(number < numbers::internal_face_boundary_id, ExcIndexRange(number,0,numbers::internal_face_boundary_id));
+  Assert(number < numbers::internal_face_boundary_id,
+        ExcIndexRange(number,0,numbers::internal_face_boundary_id));
 
   boundary[m_number] = &boundary_object;
 }
@@ -9701,12 +9695,13 @@ Triangulation<dim, spacedim>::set_manifold (const types::manifold_id m_number)
 
 
 template <int dim, int spacedim>
-const Boundary<dim, spacedim> &
+const Manifold<spacedim> &
 Triangulation<dim, spacedim>::get_boundary (const types::manifold_id m_number) const
 {
   //look, if there is a boundary stored at
   //boundary_id number.
-  typename std::map<types::manifold_id, SmartPointer<const Boundary<dim, spacedim>, Triangulation<dim, spacedim> > >::const_iterator it
+  typename std::map<types::manifold_id,
+    SmartPointer<const Manifold<spacedim>, Triangulation<dim, spacedim> > >::const_iterator it
     = boundary.find(m_number);
 
   if (it != boundary.end())
@@ -9719,7 +9714,7 @@ Triangulation<dim, spacedim>::get_boundary (const types::manifold_id m_number) c
       //if we have not found an entry
       //connected with number, we return
       //straight_boundary
-      return straight_boundary;
+      return flat_manifold;
     }
 }
 
@@ -9847,7 +9842,8 @@ copy_triangulation (const Triangulation<dim, spacedim> &old_tria)
 
   faces         = new internal::Triangulation::TriaFaces<dim>(*old_tria.faces);
 
-  typename std::map<types::boundary_id, SmartPointer<const Boundary<dim, spacedim> , Triangulation<dim, spacedim> > >::const_iterator
+  typename std::map<types::boundary_id,
+    SmartPointer<const Manifold<spacedim> , Triangulation<dim, spacedim> > >::const_iterator
   bdry_iterator = old_tria.boundary.begin();
   for (; bdry_iterator != old_tria.boundary.end() ; bdry_iterator++)
     boundary[bdry_iterator->first] = bdry_iterator->second;
index 1cbd41c77a3848531d5927b8cebe859081c39d2b..86d512782c86a57ae760ee28c2433d313624efdb 100644 (file)
@@ -487,6 +487,23 @@ HyperBallBoundary<dim,spacedim>::HyperBallBoundary (const Point<spacedim> p,
   compute_radius_automatically(false)
 {}
 
+template <int dim, int spacedim>
+Point<spacedim>
+HyperBallBoundary<dim,spacedim>::get_new_point (const std::vector<Point<spacedim> > &surrounding_points,
+              const std::vector<double> &weights) const
+{
+  Assert(surrounding_points.size() == weights.size(),
+        ExcDimensionMismatch(surrounding_points.size(), weights.size()));
+
+  double radius = 0;
+  for(unsigned int i=0; i<surrounding_points.size(); ++i)
+    radius += weights[i]*(surrounding_points[i]-center).norm();
+
+  Point<spacedim> p = FlatManifold<spacedim>::get_new_point(surrounding_points, weights);
+  p = p-center;
+  p = p/p.norm()*radius+center;  
+  return p; 
+}
 
 
 template <int dim, int spacedim>

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.