From 87d88ac4d59cb749d4055e65c032b7aafc5a5904 Mon Sep 17 00:00:00 2001 From: hartmann Date: Mon, 29 Jan 2001 09:07:40 +0000 Subject: [PATCH] Remove old version of HyperBallBoundary::get_intermediate_points_on_line that only work for special cases in 2d by a new implementation that (hopefully) works for all cases in 2d and 3d. git-svn-id: https://svn.dealii.org/trunk@3816 0785d39b-7218-0410-832d-ea1e28bc413d --- .../deal.II/source/grid/tria_boundary_lib.cc | 95 +++++++++++++------ 1 file changed, 66 insertions(+), 29 deletions(-) diff --git a/deal.II/deal.II/source/grid/tria_boundary_lib.cc b/deal.II/deal.II/source/grid/tria_boundary_lib.cc index 23012394ae..dc11bc36b3 100644 --- a/deal.II/deal.II/source/grid/tria_boundary_lib.cc +++ b/deal.II/deal.II/source/grid/tria_boundary_lib.cc @@ -75,13 +75,24 @@ get_new_point_on_quad (const typename Triangulation::quad_iterator &quad) c }; -#if deal_II_dimension == 2 +#if deal_II_dimension == 1 template <> void -HyperBallBoundary<2>::get_intermediate_points_on_line ( - const Triangulation<2>::line_iterator &line, - std::vector > &points) const +HyperBallBoundary<1>::get_intermediate_points_on_line ( + const Triangulation<1>::line_iterator &line, + std::vector > &points) const +{ + Assert(false, ExcInternalError()); +} + +#else + +template +void +HyperBallBoundary::get_intermediate_points_on_line ( + const typename Triangulation::line_iterator &line, + typename std::vector > &points) const { const unsigned int n=points.size(); Assert(n>0, ExcInternalError()); @@ -89,42 +100,68 @@ HyperBallBoundary<2>::get_intermediate_points_on_line ( points[0]=get_new_point_on_line(line); else { - Point<2> v0=line->vertex(0)-center, - v1=line->vertex(1)-center; + Point v0=line->vertex(0)-center, + v1=line->vertex(1)-center; + const double length=sqrt((v1-v0).square()); + double eps=1e-14; Assert(fabs(v0.square()-radius*radius) S; - S[0][0]=cos(alpha); - S[1][0]=sin(alpha); - S[0][1]=-S[1][0]; - S[1][1]=S[0][0]; - - contract(points[0], S, v0); - for (unsigned int i=1; i pm=StraightBoundary::get_new_point_on_line (line)-center; + + const double h=sqrt(pm.square()); + + double beta=0; + unsigned int left_index=0, right_index=0; + if ((n+1)%2==0) + { + // if the number of + // parts is even insert + // the midpoint + left_index=(n-1)/2; + right_index=left_index; + points[left_index]=pm; + --left_index; + ++right_index; + beta=d_alpha; + } + else + { + beta=(n/2+1)*d_alpha-alpha/2.; + right_index=n/2; + left_index=n/2-1; + } + + // n even: m=n/2, + // n odd: m=(n-1)/2 + const unsigned int m=n/2; + for (unsigned int i=0; i=1, ExcInternalError()); + + double d=h*tan(beta); + points[right_index]=pm+d/length*(v1-v0); + points[left_index]=pm-d/length*(v1-v0); + } + + // project the points from the + // straight line to the + // HyperBallBoundary for (unsigned int i=0; i -void -HyperBallBoundary::get_intermediate_points_on_line ( - const typename Triangulation::line_iterator &, - typename std::vector > &) const -{ - Assert(false, ExcNotImplemented()); -} - - #if deal_II_dimension == 3 template <> -- 2.39.5