From: heltai Date: Tue, 3 Dec 2013 23:07:36 +0000 (+0000) Subject: Removed double inheritance. Now all Boundaries are derived from FlatManifold. X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=6e196bc2b124a0eb256e6867aba14c32575cfdd2;p=dealii-svn.git Removed double inheritance. Now all Boundaries are derived from FlatManifold. git-svn-id: https://svn.dealii.org/branches/branch_manifold_id@31865 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/include/deal.II/grid/manifold.h b/deal.II/include/deal.II/grid/manifold.h index 36b118cceb..36e73a0515 100644 --- a/deal.II/include/deal.II/grid/manifold.h +++ b/deal.II/include/deal.II/grid/manifold.h @@ -122,7 +122,30 @@ public: * exception of type ExcPureFunctionCalled. */ virtual - void project_to_manifold (Point &candidate) const; + Point project_to_manifold (const Point 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 > &normals, + const std::vector > &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 normal_vector(const Point 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 get_new_point(const std::vector > &surrounding_points, const std::vector &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 &candidate) const; + Point project_to_manifold (const Point candidate) const; }; diff --git a/deal.II/include/deal.II/grid/tria.h b/deal.II/include/deal.II/grid/tria.h index aad2764a57..e5dc7dc20e 100644 --- a/deal.II/include/deal.II/grid/tria.h +++ b/deal.II/include/deal.II/grid/tria.h @@ -40,8 +40,6 @@ DEAL_II_NAMESPACE_OPEN -template class Boundary; -template class StraightBoundary; template class Manifold; template class FlatManifold; @@ -1234,14 +1232,6 @@ private: typedef dealii::internal::Triangulation::Iterators 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 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 &boundary_object) DEAL_II_DEPRECATED; + const Manifold &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 &get_boundary (const types::manifold_id number) const; + const Manifold &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 , Triangulation > > boundary; + std::map , Triangulation > > boundary; /** diff --git a/deal.II/include/deal.II/grid/tria_accessor.h b/deal.II/include/deal.II/grid/tria_accessor.h index 6c2dcc55b9..ec8971380a 100644 --- a/deal.II/include/deal.II/grid/tria_accessor.h +++ b/deal.II/include/deal.II/grid/tria_accessor.h @@ -36,7 +36,6 @@ template class TriaRawIterator; template class TriaIterator; template class TriaActiveIterator; -template class Boundary; template class Manifold; @@ -1161,7 +1160,7 @@ public: * Triangulation::get_boundary() function * for the boundary object. */ - const Boundary &get_boundary () const; + const Manifold &get_boundary () const; /** * Return a constant reference to a diff --git a/deal.II/include/deal.II/grid/tria_accessor.templates.h b/deal.II/include/deal.II/grid/tria_accessor.templates.h index b0e5a447d4..df31c0d49c 100644 --- a/deal.II/include/deal.II/grid/tria_accessor.templates.h +++ b/deal.II/include/deal.II/grid/tria_accessor.templates.h @@ -27,8 +27,6 @@ #include #include #include -#include -#include #include #include @@ -1938,7 +1936,7 @@ TriaAccessor::at_boundary () const template -const Boundary & +const Manifold & TriaAccessor::get_boundary () const { Assert (structdim class Triangulation; * @author Wolfgang Bangerth, 1999, 2001, 2009, Ralf Hartmann, 2001, 2008 */ template - class Boundary : public Manifold + class Boundary : public FlatManifold { public: diff --git a/deal.II/include/deal.II/grid/tria_boundary_lib.h b/deal.II/include/deal.II/grid/tria_boundary_lib.h index 2f683b964e..a6f63a3ebb 100644 --- a/deal.II/include/deal.II/grid/tria_boundary_lib.h +++ b/deal.II/include/deal.II/grid/tria_boundary_lib.h @@ -83,6 +83,15 @@ public: virtual Point get_new_point_on_line (const typename Triangulation::line_iterator &line) const; + /** + * Refer to the general documentation of + * this class and the documentation of the + * base class. + */ + virtual Point + get_new_point (const std::vector > &surrounding_points, + const std::vector &weights) const; + /** * Refer to the general documentation of * this class and the documentation of the @@ -260,6 +269,16 @@ public: Point get_new_point_on_line (const typename Triangulation::line_iterator &line) const; + /** + * Refer to the general documentation of + * this class and the documentation of the + * base class. + */ + virtual Point + get_new_point (const std::vector > &surrounding_points, + const std::vector &weights) const; + + /** * Refer to the general * documentation of this class @@ -394,6 +413,16 @@ public: Point get_new_point_on_line (const typename Triangulation::line_iterator &line) const; + /** + * Refer to the general documentation of + * this class and the documentation of the + * base class. + */ + virtual Point + get_new_point (const std::vector > &surrounding_points, + const std::vector &weights) const; + + /** * Refer to the general documentation of * this class and the documentation of the @@ -570,6 +599,16 @@ public: virtual Point get_new_point_on_quad (const typename Triangulation::quad_iterator &quad) const; + /** + * Refer to the general documentation of + * this class and the documentation of the + * base class. + */ + virtual Point + get_new_point (const std::vector > &surrounding_points, + const std::vector &weights) const; + + /** * Refer to the general * documentation of this class @@ -680,6 +719,15 @@ public: */ virtual Point get_new_point_on_line (const typename Triangulation::line_iterator &line) const; + /** + * Refer to the general documentation of + * this class and the documentation of the + * base class. + */ + virtual Point + get_new_point (const std::vector > &surrounding_points, + const std::vector &weights) const; + /** * Construct a new point on a quad. @@ -767,6 +815,16 @@ public: virtual Point get_new_point_on_quad (const typename Triangulation::quad_iterator &quad) const; + /** + * Refer to the general documentation of + * this class and the documentation of the + * base class. + */ + virtual Point + get_new_point (const std::vector > &surrounding_points, + const std::vector &weights) const; + + /** * Construct a new points on a line. */ diff --git a/deal.II/include/deal.II/numerics/vector_tools.templates.h b/deal.II/include/deal.II/numerics/vector_tools.templates.h index 77cc717739..e814f047db 100644 --- a/deal.II/include/deal.II/numerics/vector_tools.templates.h +++ b/deal.II/include/deal.II/numerics/vector_tools.templates.h @@ -4373,8 +4373,7 @@ namespace VectorTools // case in tests/deal.II/no_flux_11. Point 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, diff --git a/deal.II/source/fe/mapping_c1.cc b/deal.II/source/fe/mapping_c1.cc index 5a6115220a..caa97ba586 100644 --- a/deal.II/source/fe/mapping_c1.cc +++ b/deal.II/source/fe/mapping_c1.cc @@ -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 &boundary - = line->get_triangulation().get_boundary(line->boundary_indicator()); - - Boundary::FaceVertexNormals face_vertex_normals; - boundary.get_normals_at_vertices (line, face_vertex_normals); + const Manifold &boundary + = line->get_boundary(); + + std::vector > points(2); + std::vector > 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 diff --git a/deal.II/source/grid/grid_tools.cc b/deal.II/source/grid/grid_tools.cc index 669dccb71b..7cd43ae0fb 100644 --- a/deal.II/source/grid/grid_tools.cc +++ b/deal.II/source/grid/grid_tools.cc @@ -1728,8 +1728,7 @@ next_cell: fix_up_object (const Iterator &object, const bool respect_manifold) { - const Boundary + const Manifold *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 diff --git a/deal.II/source/grid/manifold.cc b/deal.II/source/grid/manifold.cc index d5d82db9d5..7a6de0639f 100644 --- a/deal.II/source/grid/manifold.cc +++ b/deal.II/source/grid/manifold.cc @@ -32,11 +32,30 @@ Manifold::~Manifold () +template +Point +Manifold::project_to_manifold (const Point candidate) const +{ + Assert (false, ExcPureFunctionCalled()); + return candidate; +} + template void -Manifold::project_to_manifold (Point &candidate) const +Manifold::get_normals_at_points(std::vector > &normals, + const std::vector > &points) const +{ + AssertDimension(normals.size(), points.size()); + for(unsigned int i=0; i +Point +Manifold::normal_vector(const Point point) const { Assert (false, ExcPureFunctionCalled()); + return point; } template @@ -108,8 +127,7 @@ FlatManifold:: get_new_point (const std::vector > &surrounding_points, const std::vector &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 > &surrounding_points, Point p; for(unsigned int i=0; i -void -FlatManifold::project_to_manifold (Point &candidate) const -{} +Point +FlatManifold::project_to_manifold (const Point candidate) const +{ + return candidate; +} /* -------------------------- ManifoldChart --------------------- */ @@ -144,8 +164,7 @@ ManifoldChart:: get_new_point (const std::vector > &surrounding_points, const std::vector &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 &triangulation) { - // If the codimension is - // one, we cannot perform - // this check yet. - if (spacedim>dim) return; - + std::vector wp; + std::vector > sp; for (typename Triangulation::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 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 -const StraightBoundary -Triangulation::straight_boundary = StraightBoundary(); - - template const FlatManifold Triangulation::flat_manifold = FlatManifold(); @@ -9653,10 +9646,11 @@ Triangulation::set_mesh_smoothing(const MeshSmoothing mesh_smooth template void Triangulation::set_boundary (const types::manifold_id m_number, - const Boundary &boundary_object) + const Manifold &boundary_object) { types::boundary_id number = static_cast(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::set_manifold (const types::manifold_id m_number) template -const Boundary & +const Manifold & Triangulation::get_boundary (const types::manifold_id m_number) const { //look, if there is a boundary stored at //boundary_id number. - typename std::map, Triangulation > >::const_iterator it + typename std::map, Triangulation > >::const_iterator it = boundary.find(m_number); if (it != boundary.end()) @@ -9719,7 +9714,7 @@ Triangulation::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 &old_tria) faces = new internal::Triangulation::TriaFaces(*old_tria.faces); - typename std::map , Triangulation > >::const_iterator + typename std::map , Triangulation > >::const_iterator bdry_iterator = old_tria.boundary.begin(); for (; bdry_iterator != old_tria.boundary.end() ; bdry_iterator++) boundary[bdry_iterator->first] = bdry_iterator->second; diff --git a/deal.II/source/grid/tria_boundary_lib.cc b/deal.II/source/grid/tria_boundary_lib.cc index 1cbd41c77a..86d512782c 100644 --- a/deal.II/source/grid/tria_boundary_lib.cc +++ b/deal.II/source/grid/tria_boundary_lib.cc @@ -487,6 +487,23 @@ HyperBallBoundary::HyperBallBoundary (const Point p, compute_radius_automatically(false) {} +template +Point +HyperBallBoundary::get_new_point (const std::vector > &surrounding_points, + const std::vector &weights) const +{ + Assert(surrounding_points.size() == weights.size(), + ExcDimensionMismatch(surrounding_points.size(), weights.size())); + + double radius = 0; + for(unsigned int i=0; i p = FlatManifold::get_new_point(surrounding_points, weights); + p = p-center; + p = p/p.norm()*radius+center; + return p; +} template