From 677ca5aabbebd5e395a930391a16296476a77bf9 Mon Sep 17 00:00:00 2001 From: heltai Date: Wed, 15 Jan 2014 22:55:08 +0000 Subject: [PATCH] Added back default computation of normals for FlatManifolds git-svn-id: https://svn.dealii.org/branches/branch_manifold_id@32223 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/include/deal.II/grid/tria.h | 9 +++ deal.II/include/deal.II/grid/tria_boundary.h | 22 +++++-- deal.II/source/fe/mapping_c1.cc | 5 +- deal.II/source/grid/tria.cc | 14 +++-- deal.II/source/grid/tria_boundary.cc | 62 +++++++++++++++++++- 5 files changed, 96 insertions(+), 16 deletions(-) diff --git a/deal.II/include/deal.II/grid/tria.h b/deal.II/include/deal.II/grid/tria.h index b577b3cfd1..5eea08a8bc 100644 --- a/deal.II/include/deal.II/grid/tria.h +++ b/deal.II/include/deal.II/grid/tria.h @@ -42,6 +42,7 @@ DEAL_II_NAMESPACE_OPEN template class Manifold; template class FlatManifold; +template class StraightBoundary; template class TriaAccessor; template class TriaAccessor<0,1,spacedim>; @@ -1240,6 +1241,14 @@ public: */ static const FlatManifold flat_manifold; + /** + * Default boundary object. This is used + * for those objects for which no + * boundary description has been explicitly + * set using set_manifold(). + */ + + static const StraightBoundary straight_boundary; /** * Declare some symbolic names diff --git a/deal.II/include/deal.II/grid/tria_boundary.h b/deal.II/include/deal.II/grid/tria_boundary.h index b80a2bcd2f..a2374e777a 100644 --- a/deal.II/include/deal.II/grid/tria_boundary.h +++ b/deal.II/include/deal.II/grid/tria_boundary.h @@ -128,11 +128,11 @@ public: * polygon/polyhedron defined by the boundary of the initial coarse * triangulation. * - * @deprecated The functionality of this class is equivalent, but - * less general, to that of FlatManifold. Please use - * FlatManifold instead of this - * class. StraightBoundary will be removed in future - * releases. + * The functionality of this class is equivalent, in most aspects, + * to that of FlatManifold. The only difference between + * this class and the FlatManifold class is that for this + * class the normal is well defined, and therefore there exist the + * specialized method which computes it. * * @ingroup boundary * @@ -140,13 +140,23 @@ public: * Heltai 2013, 2014 */ template -class StraightBoundary : public Boundary +class StraightBoundary : public FlatManifold { public: /** * Default constructor. Some compilers require this for some reasons. */ StraightBoundary (); + + /** + * Given a vector of points, return the normals to the Manifold in + * those points. Input arguments must be of the same size. + */ + virtual + void get_normals_at_points(std::vector > &normals, + const std::vector > &points) const; + + } DEAL_II_DEPRECATED; DEAL_II_NAMESPACE_CLOSE diff --git a/deal.II/source/fe/mapping_c1.cc b/deal.II/source/fe/mapping_c1.cc index 3196ececa5..9097bb47fd 100644 --- a/deal.II/source/fe/mapping_c1.cc +++ b/deal.II/source/fe/mapping_c1.cc @@ -67,11 +67,8 @@ 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 Manifold &boundary - = line->get_boundary(); - std::vector > face_vertex_normals(2); - boundary.get_normals_at_points (face_vertex_normals, points); + line->get_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/tria.cc b/deal.II/source/grid/tria.cc index aaae9a3da4..f363e7c15f 100644 --- a/deal.II/source/grid/tria.cc +++ b/deal.II/source/grid/tria.cc @@ -27,6 +27,7 @@ #include #include #include +#include #include #include #include @@ -9547,6 +9548,9 @@ template const FlatManifold Triangulation::flat_manifold = FlatManifold(); +template +const StraightBoundary +Triangulation::straight_boundary = StraightBoundary(); template const unsigned int @@ -9714,7 +9718,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 flat_manifold; + return straight_boundary; } } @@ -9723,8 +9727,8 @@ template const Manifold & Triangulation::get_manifold (const types::manifold_id m_number) const { - //look, if there is a boundary stored at - //boundary_id number. + //look, if there is a manifold stored at + //manifold_id number. typename std::map, Triangulation > >::const_iterator it = manifold.find(m_number); @@ -9737,8 +9741,8 @@ Triangulation::get_manifold (const types::manifold_id m_number) c { //if we have not found an entry //connected with number, we return - //flat_manifold - return flat_manifold; + //straight_boundary + return straight_boundary; } } diff --git a/deal.II/source/grid/tria_boundary.cc b/deal.II/source/grid/tria_boundary.cc index 3ed27b97b2..d066860414 100644 --- a/deal.II/source/grid/tria_boundary.cc +++ b/deal.II/source/grid/tria_boundary.cc @@ -13,7 +13,8 @@ // the top level of the deal.II distribution. // // --------------------------------------------------------------------- - + +#include #include DEAL_II_NAMESPACE_OPEN @@ -31,6 +32,65 @@ template StraightBoundary::StraightBoundary () {} +template +void +StraightBoundary::get_normals_at_points(std::vector > &normals, + const std::vector > &points) const +{ + AssertDimension(normals.size(), points.size()); + switch(dim) + { + case 1: + Assert(false, ExcNotImplemented()); + break; + case 2: + { + Assert(points.size() >= 2, ExcMessage("At least two points!")); + for(unsigned int i=0; i t = points[1]-points[0]; + Point n; n[0] = t[1]; n[1] = -t[0]; + normals[i] = n; + } + } + break; + case 3: + { + AssertDimension(points.size(), 4); + + const unsigned int vertices_per_face = GeometryInfo<3>::vertices_per_face; + static const unsigned int neighboring_vertices[4][2]= + { {1,2},{3,0},{0,3},{2,1}}; + for (unsigned int vertex=0; vertex tangents[2] + = { points[neighboring_vertices[vertex][0]] + - points[vertex], + points[neighboring_vertices[vertex][1]] + - points[vertex] }; + + // then compute the normal by + // taking the cross + // product. since the normal is + // not required to be + // normalized, no problem here + cross_product (normals[vertex], + tangents[0], tangents[1]); + }; + } + break; + default: + Assert(false, ExcNotImplemented()); + break; + } +} + + // explicit instantiations #include "tria_boundary.inst" -- 2.39.5