From d5acb5a70adeb15d5044589ada6566c0358ea380 Mon Sep 17 00:00:00 2001 From: Luca Heltai Date: Sun, 8 Apr 2018 19:27:24 +0200 Subject: [PATCH] Added clone members to (deprecated) Boundary objects. --- include/deal.II/grid/tria_boundary.h | 27 ++++-- include/deal.II/grid/tria_boundary_lib.h | 99 +++++++++++++------- source/grid/tria_boundary.cc | 18 ++++ source/grid/tria_boundary_lib.cc | 65 +++++++++++++ tests/grid/oned_vertex_manifold.cc | 9 -- tests/manifold/chart_manifold_03_embedded.cc | 11 ++- 6 files changed, 176 insertions(+), 53 deletions(-) diff --git a/include/deal.II/grid/tria_boundary.h b/include/deal.II/grid/tria_boundary.h index 8aa25aba87..d50da52db7 100644 --- a/include/deal.II/grid/tria_boundary.h +++ b/include/deal.II/grid/tria_boundary.h @@ -91,6 +91,10 @@ public: */ virtual ~Boundary () = default; + /** + * Clone this Boundary. + */ + virtual std::unique_ptr > clone() const override; /** * Return intermediate points on a line spaced according to the interior @@ -249,6 +253,11 @@ public: */ StraightBoundary (); + /** + * Clone this Boundary. + */ + virtual std::unique_ptr > clone() const override; + /** * Let the new point be the arithmetic mean of the two vertices of the line. * @@ -256,7 +265,7 @@ public: * the base class for more information. */ virtual Point - get_new_point_on_line (const typename Triangulation::line_iterator &line) const; + get_new_point_on_line (const typename Triangulation::line_iterator &line) const override; /** * Let the new point be the arithmetic mean of the four vertices of this @@ -268,7 +277,7 @@ public: */ virtual Point - get_new_point_on_quad (const typename Triangulation::quad_iterator &quad) const; + get_new_point_on_quad (const typename Triangulation::quad_iterator &quad) const override; /** * Return n=points.size() points that split the StraightBoundary @@ -280,7 +289,7 @@ public: virtual void get_intermediate_points_on_line (const typename Triangulation::line_iterator &line, - std::vector > &points) const; + std::vector > &points) const override; /** * Return n=points.size()=m*m points that splits the @@ -292,7 +301,7 @@ public: virtual void get_intermediate_points_on_quad (const typename Triangulation::quad_iterator &quad, - std::vector > &points) const; + std::vector > &points) const override; /** * Implementation of the function declared in the base class. @@ -303,7 +312,7 @@ public: virtual Tensor<1,spacedim> normal_vector (const typename Triangulation::face_iterator &face, - const Point &p) const; + const Point &p) const override; /** * Compute the normals to the boundary at the vertices of the given face. @@ -314,7 +323,7 @@ public: virtual void get_normals_at_vertices (const typename Triangulation::face_iterator &face, - typename Boundary::FaceVertexNormals &face_vertex_normals) const; + typename Boundary::FaceVertexNormals &face_vertex_normals) const override; /** * Given a candidate point and a line segment characterized by the iterator, @@ -333,7 +342,7 @@ public: virtual Point project_to_surface (const typename Triangulation::line_iterator &line, - const Point &candidate) const; + const Point &candidate) const override; /** * Same function as above but for a point that is to be projected onto the @@ -349,7 +358,7 @@ public: virtual Point project_to_surface (const typename Triangulation::quad_iterator &quad, - const Point &candidate) const; + const Point &candidate) const override; /** * Same function as above but for a point that is to be projected onto the @@ -366,7 +375,7 @@ public: virtual Point project_to_surface (const typename Triangulation::hex_iterator &hex, - const Point &candidate) const; + const Point &candidate) const override; }; diff --git a/include/deal.II/grid/tria_boundary_lib.h b/include/deal.II/grid/tria_boundary_lib.h index 321f6d0a85..ff256667a2 100644 --- a/include/deal.II/grid/tria_boundary_lib.h +++ b/include/deal.II/grid/tria_boundary_lib.h @@ -72,19 +72,24 @@ public: const Point &direction, const Point &point_on_axis); + /** + * Clone this Boundary object. + */ + virtual std::unique_ptr > clone() const override; + /** * Refer to the general documentation of this class and the documentation of * the base class. */ virtual Point - get_new_point_on_line (const typename Triangulation::line_iterator &line) const; + get_new_point_on_line (const typename Triangulation::line_iterator &line) const override; /** * Refer to the general documentation of this class and the documentation of * the base class. */ virtual Point - get_new_point_on_quad (const typename Triangulation::quad_iterator &quad) const; + get_new_point_on_quad (const typename Triangulation::quad_iterator &quad) const override; /** * Refer to the general documentation of this class and the documentation of @@ -94,7 +99,7 @@ public: */ virtual void get_intermediate_points_on_line (const typename Triangulation::line_iterator &line, - std::vector > &points) const; + std::vector > &points) const override; /** * Refer to the general documentation of this class and the documentation of @@ -104,7 +109,7 @@ public: */ virtual void get_intermediate_points_on_quad (const typename Triangulation::quad_iterator &quad, - std::vector > &points) const; + std::vector > &points) const override; /** * Compute the normals to the boundary at the vertices of the given face. @@ -114,7 +119,7 @@ public: */ virtual void get_normals_at_vertices (const typename Triangulation::face_iterator &face, - typename Boundary::FaceVertexNormals &face_vertex_normals) const; + typename Boundary::FaceVertexNormals &face_vertex_normals) const override; /** * Return the radius of the cylinder. @@ -216,6 +221,11 @@ public: const Point x_0, const Point x_1); + /** + * Clone this Boundary object. + */ + virtual std::unique_ptr > clone() const override; + /** * Return the radius of the (truncated) cone at given point x on * the axis. @@ -228,7 +238,7 @@ public: */ virtual Point - get_new_point_on_line (const typename Triangulation::line_iterator &line) const; + get_new_point_on_line (const typename Triangulation::line_iterator &line) const override; /** * Refer to the general documentation of this class and the documentation of @@ -236,7 +246,7 @@ public: */ virtual Point - get_new_point_on_quad (const typename Triangulation::quad_iterator &quad) const; + get_new_point_on_quad (const typename Triangulation::quad_iterator &quad) const override; /** * Refer to the general documentation of this class and the documentation of @@ -247,7 +257,7 @@ public: virtual void get_intermediate_points_on_line (const typename Triangulation::line_iterator &line, - std::vector > &points) const; + std::vector > &points) const override; /** * Refer to the general documentation of this class and the documentation of @@ -258,7 +268,7 @@ public: virtual void get_intermediate_points_on_quad (const typename Triangulation::quad_iterator &quad, - std::vector > &points) const; + std::vector > &points) const override; /** * Compute the normals to the boundary at the vertices of the given face. @@ -269,7 +279,7 @@ public: virtual void get_normals_at_vertices (const typename Triangulation::face_iterator &face, - typename Boundary::FaceVertexNormals &face_vertex_normals) const; + typename Boundary::FaceVertexNormals &face_vertex_normals) const override; /** * Implementation of the function declared in the base class. @@ -280,7 +290,7 @@ public: virtual Tensor<1,dim> normal_vector (const typename Triangulation::face_iterator &face, - const Point &p) const; + const Point &p) const override; protected: @@ -351,13 +361,18 @@ public: HyperBallBoundary (const Point p = Point(), const double radius = 1.0); + /** + * Clone this Boundary object. + */ + virtual std::unique_ptr > clone() const override; + /** * Refer to the general documentation of this class and the documentation of * the base class. */ virtual Point - get_new_point_on_line (const typename Triangulation::line_iterator &line) const; + get_new_point_on_line (const typename Triangulation::line_iterator &line) const override; /** * Refer to the general documentation of this class and the documentation of @@ -365,7 +380,7 @@ public: */ virtual Point - get_new_point_on_quad (const typename Triangulation::quad_iterator &quad) const; + get_new_point_on_quad (const typename Triangulation::quad_iterator &quad) const override; /** * Refer to the general documentation of this class and the documentation of @@ -376,7 +391,7 @@ public: virtual void get_intermediate_points_on_line (const typename Triangulation::line_iterator &line, - std::vector > &points) const; + std::vector > &points) const override; /** * Refer to the general documentation of this class and the documentation of @@ -387,7 +402,7 @@ public: virtual void get_intermediate_points_on_quad (const typename Triangulation::quad_iterator &quad, - std::vector > &points) const; + std::vector > &points) const override; /** * Implementation of the function declared in the base class. @@ -398,7 +413,7 @@ public: virtual Tensor<1,spacedim> normal_vector (const typename Triangulation::face_iterator &face, - const Point &p) const; + const Point &p) const override; /** * Compute the normals to the boundary at the vertices of the given face. @@ -409,7 +424,7 @@ public: virtual void get_normals_at_vertices (const typename Triangulation::face_iterator &face, - typename Boundary::FaceVertexNormals &face_vertex_normals) const; + typename Boundary::FaceVertexNormals &face_vertex_normals) const override; /** * Return the center of the ball. @@ -493,17 +508,22 @@ public: HalfHyperBallBoundary (const Point p = Point(), const double radius = 1.0); + /** + * Clone this Boundary object. + */ + virtual std::unique_ptr > clone() const override; + /** * Check if on the line x==0, otherwise pass to the base class. */ virtual Point - get_new_point_on_line (const typename Triangulation::line_iterator &line) const; + get_new_point_on_line (const typename Triangulation::line_iterator &line) const override; /** * Check if on the line x==0, otherwise pass to the base class. */ virtual Point - get_new_point_on_quad (const typename Triangulation::quad_iterator &quad) const; + get_new_point_on_quad (const typename Triangulation::quad_iterator &quad) const override; /** * Refer to the general documentation of this class and the documentation of @@ -513,7 +533,7 @@ public: */ virtual void get_intermediate_points_on_line (const typename Triangulation::line_iterator &line, - std::vector > &points) const; + std::vector > &points) const override; /** * Refer to the general documentation of this class and the documentation of @@ -523,7 +543,7 @@ public: */ virtual void get_intermediate_points_on_quad (const typename Triangulation::quad_iterator &quad, - std::vector > &points) const; + std::vector > &points) const override; /** * Compute the normals to the boundary at the vertices of the given face. @@ -533,7 +553,7 @@ public: */ virtual void get_normals_at_vertices (const typename Triangulation::face_iterator &face, - typename Boundary::FaceVertexNormals &face_vertex_normals) const; + typename Boundary::FaceVertexNormals &face_vertex_normals) const override; }; @@ -562,6 +582,11 @@ public: * radius as argument. This radius will be ignored */ HyperShellBoundary (const Point ¢er = Point()); + + /** + * Clone this Boundary object. + */ + virtual std::unique_ptr > clone() const override; }; @@ -597,17 +622,22 @@ public: const double inner_radius = -1, const double outer_radius = -1); + /** + * Clone this Boundary object. + */ + virtual std::unique_ptr > clone() const override; + /** * Construct a new point on a line. */ virtual Point - get_new_point_on_line (const typename Triangulation::line_iterator &line) const; + get_new_point_on_line (const typename Triangulation::line_iterator &line) const override; /** * Construct a new point on a quad. */ virtual Point - get_new_point_on_quad (const typename Triangulation::quad_iterator &quad) const; + get_new_point_on_quad (const typename Triangulation::quad_iterator &quad) const override; /** * Refer to the general documentation of this class and the documentation of @@ -617,7 +647,7 @@ public: */ virtual void get_intermediate_points_on_line (const typename Triangulation::line_iterator &line, - std::vector > &points) const; + std::vector > &points) const override; /** * Refer to the general documentation of this class and the documentation of @@ -627,7 +657,7 @@ public: */ virtual void get_intermediate_points_on_quad (const typename Triangulation::quad_iterator &quad, - std::vector > &points) const; + std::vector > &points) const override; /** * Compute the normals to the boundary at the vertices of the given face. @@ -637,7 +667,7 @@ public: */ virtual void get_normals_at_vertices (const typename Triangulation::face_iterator &face, - typename Boundary::FaceVertexNormals &face_vertex_normals) const; + typename Boundary::FaceVertexNormals &face_vertex_normals) const override; private: /** @@ -668,32 +698,37 @@ public: */ TorusBoundary (const double R, const double r); + /** + * Clone this Boundary object. + */ + virtual std::unique_ptr > clone() const override; + //Boundary Refinement Functions /** * Construct a new point on a line. */ virtual Point - get_new_point_on_line (const typename Triangulation::line_iterator &line) const; + get_new_point_on_line (const typename Triangulation::line_iterator &line) const override; /** * Construct a new point on a quad. */ virtual Point - get_new_point_on_quad (const typename Triangulation::quad_iterator &quad) const; + get_new_point_on_quad (const typename Triangulation::quad_iterator &quad) const override; /** * Construct a new points on a line. */ virtual void get_intermediate_points_on_line ( const typename Triangulation< dim, spacedim >::line_iterator &line, - std::vector< Point< spacedim > > &points) const; + std::vector< Point< spacedim > > &points) const override; /** * Construct a new points on a quad. */ virtual void get_intermediate_points_on_quad ( const typename Triangulation< dim, spacedim >::quad_iterator &quad, - std::vector< Point< spacedim > > &points ) const; + std::vector< Point< spacedim > > &points ) const override; /** * Get the normal from cartesian coordinates. This normal does not have unit @@ -701,7 +736,7 @@ public: */ virtual void get_normals_at_vertices ( const typename Triangulation< dim, spacedim >::face_iterator &face, - typename Boundary::FaceVertexNormals &face_vertex_normals) const; + typename Boundary::FaceVertexNormals &face_vertex_normals) const override; private: //Handy functions diff --git a/source/grid/tria_boundary.cc b/source/grid/tria_boundary.cc index 82a64a2829..6cbab1f677 100644 --- a/source/grid/tria_boundary.cc +++ b/source/grid/tria_boundary.cc @@ -29,6 +29,16 @@ DEAL_II_NAMESPACE_OPEN /* -------------------------- Boundary --------------------- */ +template +std::unique_ptr > +Boundary::clone() const +{ + AssertThrow(false, ExcPureFunctionCalled()); + return std::unique_ptr >(); +} + + + template void Boundary:: @@ -194,6 +204,14 @@ get_new_point_on_line (const typename Triangulation::line_iterato } +template +std::unique_ptr > +StraightBoundary::clone() const +{ + return std_cxx14::make_unique >(); +} + + namespace { // compute the new midpoint of a quad -- diff --git a/source/grid/tria_boundary_lib.cc b/source/grid/tria_boundary_lib.cc index 62b44b7783..f0864857cf 100644 --- a/source/grid/tria_boundary_lib.cc +++ b/source/grid/tria_boundary_lib.cc @@ -18,6 +18,7 @@ #include #include #include +#include #include @@ -47,6 +48,16 @@ CylinderBoundary::CylinderBoundary (const double radius, {} + +template +std::unique_ptr > +CylinderBoundary::clone() const +{ + return std_cxx14::make_unique >(radius, direction, point_on_axis); +} + + + template Point CylinderBoundary::get_axis_vector (const unsigned int axis) @@ -285,6 +296,14 @@ ConeBoundary::ConeBoundary (const double radius_0, +template +std::unique_ptr > ConeBoundary::clone() const +{ + return std_cxx14::make_unique >(radius_0,radius_1,x_0,x_1); +} + + + template double ConeBoundary::get_radius (Point x) const { @@ -515,6 +534,15 @@ HyperBallBoundary::HyperBallBoundary (const Point p, +template +std::unique_ptr > +HyperBallBoundary::clone() const +{ + return std_cxx14::make_unique >(center, radius); +} + + + template Point HyperBallBoundary::get_new_point_on_line (const typename Triangulation::line_iterator &line) const @@ -820,6 +848,15 @@ HalfHyperBallBoundary::HalfHyperBallBoundary (const Point center, +template +std::unique_ptr > HalfHyperBallBoundary::clone() const +{ + return std_cxx14::make_unique >(this->get_center(), + this->get_radius()); +} + + + template Point HalfHyperBallBoundary:: @@ -960,6 +997,16 @@ HyperShellBoundary::HyperShellBoundary (const Point ¢er) } + +template +std::unique_ptr > +HyperShellBoundary::clone() const +{ + return std_cxx14::make_unique >(this->get_center()); +} + + + /* ---------------------------------------------------------------------- */ @@ -983,6 +1030,17 @@ HalfHyperShellBoundary::HalfHyperShellBoundary (const Point ¢er, +template +std::unique_ptr > +HalfHyperShellBoundary::clone() const +{ + return std_cxx14::make_unique >(this->get_center(), + inner_radius, + outer_radius); +} + + + template Point HalfHyperShellBoundary:: @@ -1240,6 +1298,13 @@ TorusBoundary<2,3>::TorusBoundary (const double R_, Assert (R>r, ExcMessage("Outer radius must be greater than inner radius.")); } +template +std::unique_ptr > +TorusBoundary::clone() const +{ + return std_cxx14::make_unique >(R,r); +} + template diff --git a/tests/grid/oned_vertex_manifold.cc b/tests/grid/oned_vertex_manifold.cc index 79838a2fcc..c14d2c521d 100644 --- a/tests/grid/oned_vertex_manifold.cc +++ b/tests/grid/oned_vertex_manifold.cc @@ -51,15 +51,6 @@ int main() deallog << "current vertex manifold id: " << cell->face(vertex_n)->manifold_id() << std::endl; - - if (cell->face(vertex_n)->manifold_id() == 42) - { - // since the reference is const, the pointer we compare it to - // must also be a to a const object - const auto &spherical_reference = spherical_manifold; - Assert(&spherical_reference == &cell->face(vertex_n)->get_manifold(), - ExcMessage("manifolds should be the same.")); - } } } } diff --git a/tests/manifold/chart_manifold_03_embedded.cc b/tests/manifold/chart_manifold_03_embedded.cc index f73486a431..c2cff33f69 100644 --- a/tests/manifold/chart_manifold_03_embedded.cc +++ b/tests/manifold/chart_manifold_03_embedded.cc @@ -36,10 +36,15 @@ public: ChartManifold (periodicity) {} + virtual std::unique_ptr > + clone () const override + { + return std::unique_ptr >(new MyFlatManifold(this->get_periodicity())); + } virtual Point - pull_back(const Point &space_point) const + pull_back(const Point &space_point) const override { Point p; for (unsigned int d=0; d - push_forward(const Point &chart_point) const + push_forward(const Point &chart_point) const override { Point p; for (unsigned int d=0; d - push_forward_gradient(const Point &chart_point) const + push_forward_gradient(const Point &chart_point) const override { DerivativeForm<1,spacedim+1,spacedim> x; for (unsigned int d=0; d