From 271717a05f087b08371b1d9041c1eca2e491fbee Mon Sep 17 00:00:00 2001 From: heltai Date: Wed, 29 Jan 2014 10:44:06 +0000 Subject: [PATCH] Simplified interface of Manifold. Removed single argument versions of normal_vector and project_to_manifold git-svn-id: https://svn.dealii.org/branches/branch_manifold_id@32341 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/include/deal.II/grid/manifold.h | 59 +++++++------------ deal.II/include/deal.II/grid/tria_boundary.h | 15 ++--- .../include/deal.II/grid/tria_boundary_lib.h | 31 ++++++---- deal.II/source/grid/manifold.cc | 27 ++------- deal.II/source/grid/tria_boundary_lib.cc | 26 +++++--- 5 files changed, 71 insertions(+), 87 deletions(-) diff --git a/deal.II/include/deal.II/grid/manifold.h b/deal.II/include/deal.II/grid/manifold.h index 5743e42d6c..1b0f147d5c 100644 --- a/deal.II/include/deal.II/grid/manifold.h +++ b/deal.II/include/deal.II/grid/manifold.h @@ -92,32 +92,12 @@ public: * * Derived classes do not need to implement this function unless * mesh smoothing algorithms are used with a particular boundary - * object. The default implementation of this function calls the - * method with the same name, but without surrounding_points. - * - * This version of the projection function can be overwritten when - * only an approximation of the manifold is available, and this - * approximation requires surrounding points to be defined. + * object. */ virtual Point project_to_manifold (const std::vector > &surrounding_points, const Point &candidate) const; - /** - * Given a point which lies close to the given manifold, it modifies - * it and projects it to manifold itself. This function is used in - * some mesh smoothing algorithms that try to move around points in - * order to improve the mesh quality but need to ensure that points - * that were on the manfold remain on the manifold. - * - * Derived classes do not need to implement this function unless - * mesh smoothing algorithms are used with a particular boundary - * object. The default implementation of this function throws an - * exception of type ExcPureFunctionCalled. - */ - virtual - 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 @@ -130,20 +110,6 @@ public: void get_normals_at_points(const std::vector > &points, std::vector > &normals) const; - /** - * Compute the normal to the manifold at the given point. The - * default implementation of this function just calls - * normal_vector. The difference from this function and that one, - * is that one can use this function when only an approximation of - * the normal is required, for example when computing normals at - * quadrature points. Note that this function only makes sense - * when the manifold is a codimension one manifold. - */ - virtual - Point - normal_vector(const std::vector > &vertices, - const Point &point) const; - /** * Compute the normal to the manifold at the given point. The * default implementation of this function throws an @@ -151,7 +117,9 @@ public: * manifold is a codimension one manifold. */ virtual - Point normal_vector(const Point &point) const; + Point + normal_vector(const std::vector > &vertices, + const Point &point) const; }; @@ -205,8 +173,6 @@ public: const std::vector &weights) const; - using Manifold::project_to_manifold; - /** * Project to FlatManifold: do nothing. Note however that this * function can be overloaded by derived classes, which will then @@ -215,9 +181,24 @@ public: * this class. */ virtual - Point project_to_manifold (const Point &candidate) const; + Point project_to_manifold (const std::vector > &points, + const Point &candidate) const; private: + /** + * The periodicity of this Manifold. Periodicity affects the way a + * middle point is computed. It is assumed that if two points are + * more than half period distant, then the distance should be + * computed by crossing the periodicity boundary, i.e., the average + * is computed by adding a full period to the sum of the two. For + * example, if along direction 0 we have 2*pi periodicity, then the + * average of (2*pi-eps) and (eps) is not pi, but 2*pi (or zero), + * since, on a periodic manifold, these two points are at distance + * 2*eps and not (2*pi-eps). + * + * A periodicity 0 along one direction means no periodicity. This is + * the default value for all directions. + */ const Point periodicity; }; diff --git a/deal.II/include/deal.II/grid/tria_boundary.h b/deal.II/include/deal.II/grid/tria_boundary.h index 72b31d4db4..882304fa53 100644 --- a/deal.II/include/deal.II/grid/tria_boundary.h +++ b/deal.II/include/deal.II/grid/tria_boundary.h @@ -155,14 +155,15 @@ public: virtual void get_normals_at_points(const std::vector > &points, std::vector > &normals) const; + /** - * Given a point and a vector of vertices, compute an approximation - * of the normal to the manifold. - */ - virtual - Point - normal_vector(const std::vector > &vertices, - const Point &point) const; + * Given a point and a vector of vertices, compute an approximation + * of the normal to the manifold. + */ + virtual + Point + normal_vector(const std::vector > &vertices, + const Point &point) const; }; 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 a83b3d89a6..8daf006889 100644 --- a/deal.II/include/deal.II/grid/tria_boundary_lib.h +++ b/deal.II/include/deal.II/grid/tria_boundary_lib.h @@ -71,21 +71,24 @@ public: CylinderBoundary (const double radius, const Point &direction, const Point &point_on_axis); - + /** * Since this class is derived from FlatManifold, we need to * overload only the project_to_manifold function in order for the * boundary to function properly. */ virtual Point - project_to_manifold (const Point &candidate) const; + project_to_manifold (const std::vector > & vertices, + const Point &candidate) const; + /** * Compute the normal to the surface. */ - virtual Point - normal_vector (const Point &point) const; - + virtual Point + normal_vector (const std::vector > & vertices, + const Point &point) const; + /** * Return the radius of the cylinder. */ @@ -186,7 +189,8 @@ public: * boundary to function properly. */ virtual Point - project_to_manifold (const Point &candidate) const; + project_to_manifold (const std::vector > & vertices, + const Point &candidate) const; /** @@ -197,7 +201,8 @@ public: */ virtual Point - normal_vector (const Point &point) const; + normal_vector (const std::vector > & vertices, + const Point &point) const; protected: /** @@ -266,7 +271,8 @@ public: */ virtual Point - normal_vector (const Point &p) const; + normal_vector (const std::vector > & vertices, + const Point &p) const; /** @@ -370,7 +376,8 @@ public: */ virtual Point - normal_vector (const Point &p) const; + normal_vector (const std::vector > & vertices, + const Point &p) const; }; @@ -442,7 +449,8 @@ public: * Compute the normal to the surface. */ virtual Point - normal_vector (const Point &point) const; + normal_vector (const std::vector > & vertices, + const Point &point) const; private: /** @@ -495,7 +503,8 @@ public: * Get the normal from cartesian coordinates. This normal does not * have unit length. */ - Point normal_vector(const Point &p) const; + Point normal_vector(const std::vector > & vertices, + const Point &p) const; private: double get_correct_angle(const double angle,const double x,const double y) const; diff --git a/deal.II/source/grid/manifold.cc b/deal.II/source/grid/manifold.cc index dc59e03bf1..fe6243d883 100644 --- a/deal.II/source/grid/manifold.cc +++ b/deal.II/source/grid/manifold.cc @@ -30,22 +30,13 @@ Manifold::~Manifold () {} - - -template -Point -Manifold::project_to_manifold (const Point &candidate) const -{ - Assert (false, ExcPureFunctionCalled()); - return candidate; -} - template Point Manifold::project_to_manifold (const std::vector > &, const Point &candidate) const -{ - return project_to_manifold(candidate); +{ + Assert (false, ExcPureFunctionCalled()); + return candidate; } @@ -56,20 +47,13 @@ Manifold::get_normals_at_points(const std::vector > &p { AssertDimension(normals.size(), points.size()); for(unsigned int i=0; i Point Manifold::normal_vector(const std::vector > &, const Point &point) const -{ - return normal_vector(point); -} - -template -Point -Manifold::normal_vector(const Point &point) const { Assert (false, ExcPureFunctionCalled()); return point; @@ -122,7 +106,8 @@ get_new_point (const std::vector > &surrounding_points, template Point -FlatManifold::project_to_manifold (const Point &candidate) const +FlatManifold::project_to_manifold (const std::vector > & vertices, + const Point &candidate) const { return candidate; } diff --git a/deal.II/source/grid/tria_boundary_lib.cc b/deal.II/source/grid/tria_boundary_lib.cc index 831775a029..77044ede78 100644 --- a/deal.II/source/grid/tria_boundary_lib.cc +++ b/deal.II/source/grid/tria_boundary_lib.cc @@ -61,7 +61,8 @@ CylinderBoundary::get_axis_vector (const unsigned int axis) template Point CylinderBoundary:: -project_to_manifold (const Point &middle) const +project_to_manifold (const std::vector > & vertices, + const Point &middle) const { // The temptative point is computed from the get_new_point function // of the FlatManifold class, which in turns calls this class @@ -85,7 +86,8 @@ project_to_manifold (const Point &middle) const template Point CylinderBoundary:: -normal_vector (const Point &vertex) const +normal_vector (const std::vector > & vertices, + const Point &vertex) const { const Point vector_from_axis = (vertex-point_on_axis) - ((vertex-point_on_axis) * direction) * direction; @@ -134,7 +136,8 @@ double ConeBoundary::get_radius (Point x) const template Point ConeBoundary:: -project_to_manifold (const Point &middle) const +project_to_manifold (const std::vector > & vertices, + const Point &middle) const { const Point axis = x_1 - x_0; // To project it on the boundary of the cone we first compute the orthogonal @@ -149,7 +152,8 @@ project_to_manifold (const Point &middle) const template Point ConeBoundary:: -normal_vector (const Point &vertex) const +normal_vector (const std::vector > & vertices, + const Point &vertex) const { const Point axis = x_1 - x_0; @@ -270,7 +274,8 @@ HyperBallBoundary::get_intermediate_point(const Point &p template Point HyperBallBoundary:: -normal_vector (const Point &vertex) const +normal_vector (const std::vector > & vertices, + const Point &vertex) const { return vertex-center; } @@ -327,7 +332,8 @@ HalfHyperBallBoundary::get_new_point(const std::vector > &surrou template Point HalfHyperBallBoundary:: -normal_vector (const Point &vertex) const +normal_vector (const std::vector > & vertices, + const Point &vertex) const { // Check whether the point is inside the flat part, or outside. If // it is both at distance R from the center and with zero first @@ -397,7 +403,8 @@ HalfHyperShellBoundary::HalfHyperShellBoundary (const Point ¢er, template Point HalfHyperShellBoundary:: -normal_vector (const Point &vertex) const +normal_vector (const std::vector > & vertices, + const Point &vertex) const { if (vertex(0) == this->center(0)) { @@ -406,7 +413,7 @@ normal_vector (const Point &vertex) const return n; } else - return HyperShellBoundary::normal_vector (vertex); + return HyperShellBoundary::normal_vector (vertices, vertex); } @@ -523,7 +530,8 @@ TorusBoundary<2,3>:: get_surf_norm_from_sp(const Point<2> &surfP) const //Normal field without unit length template <> Point<3> -TorusBoundary<2,3>::normal_vector(const Point<3> &p) const +TorusBoundary<2,3>::normal_vector(const std::vector > & vertices, + const Point<3> &p) const { Point<2> surfP=pull_back(p); -- 2.39.5