]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Remove deprecated function Manifold::get_new_point(Quadrature).
authorWolfgang Bangerth <bangerth@colostate.edu>
Fri, 5 May 2017 17:09:12 +0000 (11:09 -0600)
committerWolfgang Bangerth <bangerth@colostate.edu>
Fri, 5 May 2017 17:09:12 +0000 (11:09 -0600)
include/deal.II/grid/manifold.h
include/deal.II/grid/manifold_lib.h
source/grid/manifold.cc
source/grid/manifold_lib.cc

index a99688ee169188b3cc1f32f31f2b835f2e2d77ed..c7712ab24802aff2df3fea71c38fcf4f820c6962 100644 (file)
@@ -355,27 +355,6 @@ public:
                           const Point<spacedim> &p2,
                           const double w) const;
 
-  /**
-   * Return the point which shall become the new vertex surrounded by the
-   * given points which make up the quadrature. We use a quadrature object,
-   * which should be filled with the surrounding points together with
-   * appropriate weights.
-   *
-   * In its default implementation it uses a pair-wise reduction of
-   * the points in the quadrature formula by calling the function
-   * get_intermediate_point() on the first two points, then on the
-   * resulting point and the next, until all points in the quadrature
-   * have been taken into account. User classes can get away by simply
-   * implementing the get_intermediate_point() function. Notice that
-   * by default the get_intermediate_point() function calls the
-   * project_to_manifold() function with the convex combination of its
-   * arguments. For simple situations you may get away by implementing
-   * only the project_to_manifold() function.
-   */
-  virtual
-  Point<spacedim>
-  get_new_point (const Quadrature<spacedim> &quad) const DEAL_II_DEPRECATED;
-
   /**
    * Return the point which shall become the new vertex surrounded by the
    * given points @p surrounding_points. @p weights contains appropriate
@@ -702,31 +681,6 @@ public:
   FlatManifold (const Tensor<1,spacedim> &periodicity = Tensor<1,spacedim>(),
                 const double tolerance=1e-10);
 
-  /**
-   * Let the new point be the average sum of surrounding vertices.
-   *
-   * This particular implementation constructs the weighted average of the
-   * surrounding points, and then calls internally the function
-   * project_to_manifold(). The reason why we do it this way, is to allow lazy
-   * programmers to implement only the project_to_manifold() function for their
-   * own Manifold classes which are small (or trivial) perturbations of a flat
-   * manifold. This is the case whenever the coarse mesh is a decent
-   * approximation of the manifold geometry. In this case, the middle point of
-   * a cell is close to true middle point of the manifold, and a projection
-   * may suffice.
-   *
-   * For most simple geometries, it is possible to get reasonable results by
-   * deriving your own Manifold class from FlatManifold, and write a new
-   * interface only for the project_to_manifold function. You will have good
-   * approximations also with large deformations, as long as in the coarsest
-   * mesh size you are trying to refine, the middle point is not too far from
-   * the manifold mid point, i.e., as long as the coarse mesh size is small
-   * enough.
-   */
-  virtual
-  Point<spacedim>
-  get_new_point(const Quadrature<spacedim> &quad) const DEAL_II_DEPRECATED;
-
   /**
    * Let the new point be the average sum of surrounding vertices.
    *
@@ -962,15 +916,6 @@ public:
    */
   virtual ~ChartManifold ();
 
-
-  /**
-   * Refer to the general documentation of this class and the documentation of
-   * the base class for more information.
-   */
-  virtual
-  Point<spacedim>
-  get_new_point(const Quadrature<spacedim> &quad) const DEAL_II_DEPRECATED;
-
   /**
    * Refer to the general documentation of this class and the documentation of
    * the base class for more information.
index cde73c68d76f6efaecebe4179e1bb1a74bce9895..cc7d5e757c9b4efce6ea5484041599d6f602d41d 100644 (file)
@@ -227,16 +227,6 @@ public:
   get_tangent_vector (const Point<spacedim> &x1,
                       const Point<spacedim> &x2) const;
 
-  /**
-   * Return a point on the spherical manifold which is intermediate
-   * with respect to the surrounding points.
-   *
-   * @deprecated Use the other function that takes points and weights separately instead.
-   */
-  virtual
-  Point<spacedim>
-  get_new_point(const dealii::Quadrature<spacedim> &quadrature) const DEAL_II_DEPRECATED;
-
   /**
    * Return a point on the spherical manifold which is intermediate
    * with respect to the surrounding points.
@@ -292,13 +282,6 @@ public:
                        const Point<spacedim> &point_on_axis,
                        const double tolerance = 1e-10);
 
-  /**
-   * Compute new points on the CylindricalManifold. See the documentation of
-   * the base class for a detailed description of what this function does.
-   */
-  virtual Point<spacedim>
-  get_new_point(const Quadrature<spacedim> &quad) const DEAL_II_DEPRECATED;
-
   /**
    * Compute new points on the CylindricalManifold. See the documentation of
    * the base class for a detailed description of what this function does.
index b080fb8aa863437df40f545e28137c3ebe3c4d7b..1de1e07e8ca4370840f85ad528e2f993e0403495 100644 (file)
@@ -80,16 +80,6 @@ get_intermediate_point (const Point<spacedim> &p1,
 
 
 
-template <int dim, int spacedim>
-Point<spacedim>
-Manifold<dim, spacedim>::
-get_new_point (const Quadrature<spacedim> &quad) const
-{
-  return get_new_point(quad.get_points(),quad.get_weights());
-}
-
-
-
 template <int dim, int spacedim>
 Point<spacedim>
 Manifold<dim, spacedim>::
@@ -551,16 +541,6 @@ FlatManifold<dim,spacedim>::FlatManifold (const Tensor<1,spacedim> &periodicity,
 
 
 
-template <int dim, int spacedim>
-Point<spacedim>
-FlatManifold<dim, spacedim>::
-get_new_point (const Quadrature<spacedim> &quad) const
-{
-  return get_new_point(quad.get_points(),quad.get_weights());
-}
-
-
-
 template <int dim, int spacedim>
 Point<spacedim>
 FlatManifold<dim, spacedim>::
@@ -745,16 +725,6 @@ ChartManifold<dim,spacedim,chartdim>::ChartManifold (const Tensor<1,chartdim> &p
 
 
 
-template <int dim, int spacedim, int chartdim>
-Point<spacedim>
-ChartManifold<dim,spacedim,chartdim>::
-get_new_point (const Quadrature<spacedim> &quad) const
-{
-  return get_new_point(quad.get_points(),quad.get_weights());
-}
-
-
-
 template <int dim, int spacedim, int chartdim>
 Point<spacedim>
 ChartManifold<dim,spacedim,chartdim>::
index b95593149a92db181c1657f96411d61cb89a4b0f..46c46eb371d4aaab87ad9fddc3e2538b947123ab 100644 (file)
@@ -269,13 +269,7 @@ get_tangent_vector (const Point<spacedim> &p1,
   return (r1-r2)*e1 + r1*gamma*tg;
 }
 
-template <int dim, int spacedim>
-Point<spacedim>
-SphericalManifold<dim, spacedim>::
-get_new_point (const Quadrature<spacedim> &quad) const
-{
-  return get_new_point(quad.get_points(),quad.get_weights());
-}
+
 
 template <int dim, int spacedim>
 Point<spacedim>
@@ -329,16 +323,6 @@ CylindricalManifold<dim,spacedim>::CylindricalManifold(const Point<spacedim> &di
 
 
 
-template <int dim, int spacedim>
-Point<spacedim>
-CylindricalManifold<dim,spacedim>::
-get_new_point (const Quadrature<spacedim> &quad) const
-{
-  return get_new_point(quad.get_points(),quad.get_weights());
-}
-
-
-
 template <int dim, int spacedim>
 Point<spacedim>
 CylindricalManifold<dim,spacedim>::
@@ -378,7 +362,7 @@ get_new_point (const std::vector<Point<spacedim> > &surrounding_points,
 
 
 // ============================================================
-// FunctionChartManifold
+// FunctionManifold
 // ============================================================
 template <int dim, int spacedim, int chartdim>
 FunctionManifold<dim,spacedim,chartdim>::FunctionManifold
@@ -396,6 +380,8 @@ FunctionManifold<dim,spacedim,chartdim>::FunctionManifold
   AssertDimension(pull_back_function.n_components, chartdim);
 }
 
+
+
 template <int dim, int spacedim, int chartdim>
 FunctionManifold<dim,spacedim,chartdim>::FunctionManifold
 (const std::string push_forward_expression,
@@ -419,6 +405,8 @@ FunctionManifold<dim,spacedim,chartdim>::FunctionManifold
   pull_back_function = pb;
 }
 
+
+
 template <int dim, int spacedim, int chartdim>
 FunctionManifold<dim,spacedim,chartdim>::~FunctionManifold()
 {
@@ -434,6 +422,8 @@ FunctionManifold<dim,spacedim,chartdim>::~FunctionManifold()
     }
 }
 
+
+
 template <int dim, int spacedim, int chartdim>
 Point<spacedim>
 FunctionManifold<dim,spacedim,chartdim>::push_forward(const Point<chartdim> &chart_point) const
@@ -458,6 +448,7 @@ FunctionManifold<dim,spacedim,chartdim>::push_forward(const Point<chartdim> &cha
 }
 
 
+
 template <int dim, int spacedim, int chartdim>
 DerivativeForm<1,chartdim, spacedim>
 FunctionManifold<dim,spacedim,chartdim>::push_forward_gradient(const Point<chartdim> &chart_point) const
@@ -472,6 +463,7 @@ FunctionManifold<dim,spacedim,chartdim>::push_forward_gradient(const Point<chart
 }
 
 
+
 template <int dim, int spacedim, int chartdim>
 Point<chartdim>
 FunctionManifold<dim,spacedim,chartdim>::pull_back(const Point<spacedim> &space_point) const
@@ -499,6 +491,8 @@ TorusManifold<dim>::pull_back(const Point<3> &p) const
   return Point<3>(phi, theta, w);
 }
 
+
+
 template <int dim>
 Point<3>
 TorusManifold<dim>::push_forward(const Point<3> &chart_point) const
@@ -514,6 +508,7 @@ TorusManifold<dim>::push_forward(const Point<3> &chart_point) const
 }
 
 
+
 template <int dim>
 TorusManifold<dim>::TorusManifold (const double R, const double r)
   : ChartManifold<dim,3,3> (Point<3>(2*numbers::PI, 2*numbers::PI, 0.0)),
@@ -525,6 +520,8 @@ TorusManifold<dim>::TorusManifold (const double R, const double r)
   Assert (r>0.0, ExcMessage("inner radius must be positive."));
 }
 
+
+
 template <int dim>
 DerivativeForm<1,3,3>
 TorusManifold<dim>::push_forward_gradient(const Point<3> &chart_point) const

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.