]> https://gitweb.dealii.org/ - dealii.git/commitdiff
SphericalManifold: deprecate the public center member. 17187/head
authorDavid Wells <drwells@email.unc.edu>
Sat, 6 Jul 2024 14:15:06 +0000 (10:15 -0400)
committerDavid Wells <drwells@email.unc.edu>
Sat, 6 Jul 2024 14:46:19 +0000 (10:46 -0400)
include/deal.II/grid/manifold_lib.h
source/grid/manifold_lib.cc

index 11c5d6ba5d9cb8388b3f8fd1fde3a4b8e14c6ba0..710398d67135a996aa27de632f266e52410b61b9 100644 (file)
@@ -348,10 +348,21 @@ public:
 
   /**
    * The center of the spherical coordinate system.
+   *
+   * @deprecated Use get_center() instead.
    */
+  DEAL_II_DEPRECATED_EARLY_WITH_COMMENT(
+    "Access the center with get_center() instead.")
   const Point<spacedim> center;
 
 private:
+  /**
+   * The center of the spherical coordinate system.
+   *
+   * @note This exists to avoid warnings when using center internally.
+   */
+  const Point<spacedim> p_center;
+
   /**
    * Return a point on the spherical manifold which is intermediate
    * with respect to the surrounding points. This function uses a linear
index 889b0175ed445b73fd652481bcbee0b7437864eb..7238642b787a14eb8be4a2a0c6e97458e119bab2 100644 (file)
@@ -359,6 +359,7 @@ template <int dim, int spacedim>
 SphericalManifold<dim, spacedim>::SphericalManifold(
   const Point<spacedim> center)
   : center(center)
+  , p_center(center)
   , polar_manifold(center)
 {}
 
@@ -392,8 +393,8 @@ SphericalManifold<dim, spacedim>::get_intermediate_point(
   if (spacedim == 1)
     return Point<spacedim>(w * p2 + (1 - w) * p1);
 
-  const Tensor<1, spacedim> v1 = p1 - center;
-  const Tensor<1, spacedim> v2 = p2 - center;
+  const Tensor<1, spacedim> v1 = p1 - p_center;
+  const Tensor<1, spacedim> v2 = p2 - p_center;
   const double              r1 = v1.norm();
   const double              r2 = v2.norm();
 
@@ -408,11 +409,11 @@ SphericalManifold<dim, spacedim>::get_intermediate_point(
 
   // Points are collinear with the center (allow for 8*eps as a tolerance)
   if (cosgamma < -1 + 8. * std::numeric_limits<double>::epsilon())
-    return center;
+    return p_center;
 
   // Points are along a line, in which case e1 and e2 are essentially the same.
   if (cosgamma > 1 - 8. * std::numeric_limits<double>::epsilon())
-    return Point<spacedim>(center + w * v2 + (1 - w) * v1);
+    return Point<spacedim>(p_center + w * v2 + (1 - w) * v1);
 
   // Find the angle sigma that corresponds to arclength equal to w. acos
   // should never be undefined because we have ruled out the two special cases
@@ -434,7 +435,7 @@ SphericalManifold<dim, spacedim>::get_intermediate_point(
   const Tensor<1, spacedim> P = std::cos(sigma) * e1 + std::sin(sigma) * n;
 
   // Project this point on the manifold.
-  return Point<spacedim>(center + (w * r2 + (1.0 - w) * r1) * P);
+  return Point<spacedim>(p_center + (w * r2 + (1.0 - w) * r1) * P);
 }
 
 
@@ -450,8 +451,8 @@ SphericalManifold<dim, spacedim>::get_tangent_vector(
 
   Assert(p1 != p2, ExcMessage("p1 and p2 should not concide."));
 
-  const Tensor<1, spacedim> v1 = p1 - center;
-  const Tensor<1, spacedim> v2 = p2 - center;
+  const Tensor<1, spacedim> v1 = p1 - p_center;
+  const Tensor<1, spacedim> v2 = p2 - p_center;
   const double              r1 = v1.norm();
   const double              r2 = v2.norm();
 
@@ -500,12 +501,12 @@ SphericalManifold<dim, spacedim>::normal_vector(
   // (tangential to the sphere).  In this case, the normal vector is
   // easy to compute since it is proportional to the vector from the
   // center to the point 'p'.
-  if (spherical_face_is_horizontal<dim, spacedim>(face, center))
+  if (spherical_face_is_horizontal<dim, spacedim>(face, p_center))
     {
       // So, if this is a "horizontal" face, then just compute the normal
       // vector as the one from the center to the point 'p', adequately
       // scaled.
-      const Tensor<1, spacedim> unnormalized_spherical_normal = p - center;
+      const Tensor<1, spacedim> unnormalized_spherical_normal = p - p_center;
       const Tensor<1, spacedim> normalized_spherical_normal =
         unnormalized_spherical_normal / unnormalized_spherical_normal.norm();
       return normalized_spherical_normal;
@@ -553,7 +554,7 @@ SphericalManifold<dim, spacedim>::get_normals_at_vertices(
   // (tangential to the sphere).  In this case, the normal vector is
   // easy to compute since it is proportional to the vector from the
   // center to the point 'p'.
-  if (spherical_face_is_horizontal<dim, spacedim>(face, center))
+  if (spherical_face_is_horizontal<dim, spacedim>(face, p_center))
     {
       // So, if this is a "horizontal" face, then just compute the normal
       // vector as the one from the center to the point 'p', adequately
@@ -561,7 +562,7 @@ SphericalManifold<dim, spacedim>::get_normals_at_vertices(
       for (unsigned int vertex = 0;
            vertex < GeometryInfo<spacedim>::vertices_per_face;
            ++vertex)
-        face_vertex_normals[vertex] = face->vertex(vertex) - center;
+        face_vertex_normals[vertex] = face->vertex(vertex) - p_center;
     }
   else
     Manifold<dim, spacedim>::get_normals_at_vertices(face, face_vertex_normals);
@@ -783,7 +784,7 @@ SphericalManifold<dim, spacedim>::do_get_new_points(
   double max_distance = 0.;
   for (unsigned int i = 0; i < surrounding_points.size(); ++i)
     {
-      directions[i] = surrounding_points[i] - center;
+      directions[i] = surrounding_points[i] - p_center;
       distances[i]  = directions[i].norm();
 
       if (distances[i] != 0.)
@@ -824,7 +825,7 @@ SphericalManifold<dim, spacedim>::do_get_new_points(
       // the Newton iteration in step 2, which would crash.
       if (new_candidates[row].first == 0.0)
         {
-          new_points[row]               = center;
+          new_points[row]               = p_center;
           accurate_point_was_found[row] = true;
           continue;
         }
@@ -852,7 +853,7 @@ SphericalManifold<dim, spacedim>::do_get_new_points(
         {
           for (unsigned int row = 0; row < weight_rows; ++row)
             new_points[row] =
-              center + new_candidates[row].first * new_candidates[row].second;
+              p_center + new_candidates[row].first * new_candidates[row].second;
 
           return;
         }
@@ -960,7 +961,7 @@ SphericalManifold<dim, spacedim>::do_get_new_points(
                 new_candidates[merged_weights_index[row]].second;
 
             new_points[row] =
-              center + new_candidates[row].first * new_candidates[row].second;
+              p_center + new_candidates[row].first * new_candidates[row].second;
           }
     }
 }

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.