// $Id$
// Version: $Name$
//
-// Copyright (C) 1998, 1999, 2000 by the deal.II authors
+// Copyright (C) 1998, 1999, 2000, 2001 by the deal.II authors
//
// This file is subject to QPL and may not be distributed
// without copyright and license information. Please refer
* @ref{HyperBallBoundary} creating a hyperball with given radius
* around a given center point.
*
- * @author Wolfgang Bangerth, 1999
+ * @author Wolfgang Bangerth, 1999, Ralf Hartmann, 2001
*/
template <int dim>
class Boundary : public Subscriptor
virtual Point<dim>
get_new_point_on_quad (const typename Triangulation<dim>::quad_iterator &quad) const;
+ /**
+ * Return intermediate points
+ * on the boundary line.
+ *
+ * The number of points requested
+ * is given by the size of the
+ * vector @p{points}. It is the
+ * task of the derived classes to
+ * arrange the points in
+ * approximately equal distances.
+ *
+ * This function is needed by the
+ * @p{MappingQ} class. As this
+ * function is not needed for Q1
+ * mappings, it is not made pure
+ * virtual, to avoid the need to
+ * overload it. The default
+ * implementation throws an error
+ * in any case, however.
+ */
+ virtual void
+ get_intermediate_points_on_line (const typename Triangulation<dim>::line_iterator &line,
+ vector<Point<dim> > &points) const;
+
/**
* Exception.
*/
* placing new points in the middle of old ones, it rather assumes that the
* boundary of the domain is given by the polygon/polyhedron defined by the
* boundary of the initial coarse triangulation.
+ *
+ * @author Wolfgang Bangerth, 1998, Ralf Hartmann, 2001
*/
template <int dim>
class StraightBoundary : public Boundary<dim> {
*/
virtual Point<dim>
get_new_point_on_quad (const typename Triangulation<dim>::quad_iterator &quad) const;
+
+ /**
+ * Gives @p{n=points.size()}
+ * points that splits the
+ * p{StraightBoundary} line into
+ * p{n+1} partitions of equal
+ * lengths.
+ *
+ * Refer to the general documentation of
+ * this class and the documentation of the
+ * base class.
+ */
+ virtual void
+ get_intermediate_points_on_line (const typename Triangulation<dim>::line_iterator &line,
+ vector<Point<dim> > &points) const;
+
};
// $Id$
// Version: $Name$
//
-// Copyright (C) 1998, 1999, 2000 by the deal.II authors
+// Copyright (C) 1998, 1999, 2000, 2001 by the deal.II authors
//
// This file is subject to QPL and may not be distributed
// without copyright and license information. Please refer
* @ref{Boundary}, which would seem natural, since this way we can use the
* @ref{StraightBoundary}@p{<dim>::in_between(neighbors)} function.
*
- * @author Wolfgang Bangerth, 1998
+ * @author Wolfgang Bangerth, 1998, Ralf Hartmann, 2001
*/
template <int dim>
class HyperBallBoundary : public StraightBoundary<dim>
virtual Point<dim>
get_new_point_on_quad (const typename Triangulation<dim>::quad_iterator &quad) const;
+ /**
+ * Refer to the general documentation of
+ * this class and the documentation of the
+ * base class.
+ *
+ * Only implemented for @p{dim=2}.
+ */
+ virtual void
+ get_intermediate_points_on_line (const typename Triangulation<dim>::line_iterator &line,
+ vector<Point<dim> > &points) const;
+
/**
* Return the center of the ball.
*/
// $Id$
// Version: $Name$
//
-// Copyright (C) 1998, 1999, 2000 by the deal.II authors
+// Copyright (C) 1998, 1999, 2000, 2001 by the deal.II authors
//
// This file is subject to QPL and may not be distributed
// without copyright and license information. Please refer
};
+template <int dim>
+void
+Boundary<dim>::get_intermediate_points_on_line (
+ const typename Triangulation<dim>::line_iterator &,
+ vector<Point<dim> > &) const
+{
+ Assert (false, ExcPureVirtualFunctionCalled());
+};
+
+
template <int dim>
Point<dim>
StraightBoundary<dim>::get_new_point_on_line (const typename Triangulation<dim>::line_iterator &line) const
#endif
+template <int dim>
+void
+StraightBoundary<dim>::get_intermediate_points_on_line (
+ const typename Triangulation<dim>::line_iterator &line,
+ vector<Point<dim> > &points) const
+{
+ const unsigned int n=points.size();
+ Assert(n>0, ExcInternalError());
+ const double part=1./(n+1);
+ double position=part;
+ for (unsigned int i=0; i<n; ++i, position+=part)
+ points[i]=(1-position)*line->vertex(0) + position*line->vertex(1);
+};
+
+
// explicit instantiations
template class Boundary<deal_II_dimension>;
template class StraightBoundary<deal_II_dimension>;
// $Id$
// Version: $Name$
//
-// Copyright (C) 1998, 1999, 2000 by the deal.II authors
+// Copyright (C) 1998, 1999, 2000, 2001 by the deal.II authors
//
// This file is subject to QPL and may not be distributed
// without copyright and license information. Please refer
};
+template <>
+void
+HyperBallBoundary<2>::get_intermediate_points_on_line (
+ const Triangulation<2>::line_iterator &line,
+ vector<Point<2> > &points) const
+{
+ const unsigned int n=points.size();
+ Assert(n>0, ExcInternalError());
+ if (n==1)
+ points[0]=get_new_point_on_line(line);
+ else
+ {
+ Point<2> v0=line->vertex(0)-center,
+ v1=line->vertex(1)-center;
+ double alpha0=atan(v0(1)/v0(0)),
+ alpha1=atan(v1(1)/v1(0));
+ double dalpha=(alpha1-alpha0)/(n+1);
+ double alpha=alpha0+dalpha;
+
+ for (unsigned int i=0; i<n; ++i, alpha+=dalpha)
+ {
+ points[i](0)=radius*cos(alpha)+center(0);
+ points[i](1)=radius*sin(alpha)+center(1);
+ }
+ }
+}
+
+
+
+template <int dim>
+void
+HyperBallBoundary<dim>::get_intermediate_points_on_line (
+ const typename Triangulation<dim>::line_iterator &,
+ vector<Point<dim> > &) const
+{
+ Assert(false, ExcNotImplemented());
+}
+
template <int dim>
Point<dim>