]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Move computation of normal vector into anonymous namespace 4483/head
authorDaniel Arndt <daniel.arndt@iwr.uni-heidelberg.de>
Tue, 6 Jun 2017 13:14:14 +0000 (15:14 +0200)
committerDaniel Arndt <daniel.arndt@iwr.uni-heidelberg.de>
Tue, 6 Jun 2017 13:14:14 +0000 (15:14 +0200)
include/deal.II/grid/manifold_lib.h
source/grid/manifold_lib.cc
tests/mappings/mapping_q_mixed_manifolds_01.cc
tests/mappings/mapping_q_mixed_manifolds_02.cc

index 53d59fed5d87f6477a88aa8f29019635d949dc0e..60483e3af22d858b9400d6aa772f0ee0d042708a 100644 (file)
@@ -283,17 +283,17 @@ public:
                        const double tolerance = 1e-10);
 
   /**
-   * Compute the cartesian coordinates for a point given in cylindrical
+   * Compute the Cartesian coordinates for a point given in cylindrical
    * coordinates.
    */
   virtual Point<3>
   pull_back(const Point<spacedim> &space_point) const override;
 
   /**
-   * Compute the cylindrical coordinates \f$(r, \phi, \lambda)\f$ for the given
-   * point where \f$r\f$ denotes the distance from the axis,
-   * \f$\phi\f$ the angle between the given point and the computed normal
-   * direction and \f$\lambda\f$ the axial position.
+   * Compute the cylindrical coordinates $(r, \phi, \lambda)$ for the given
+   * point where $r$ denotes the distance from the axis,
+   * $\phi$ the angle between the given point and the computed normal
+   * direction and $\lambda$ the axial position.
    */
   virtual Point<spacedim>
   push_forward(const Point<3> &chart_point) const override;
@@ -306,13 +306,7 @@ public:
   get_new_point(const std::vector<Point<spacedim> > &surrounding_points,
                 const std::vector<double>           &weights) const override;
 
-private:
-  /**
-   * Compute a vector that is orthogonal to the given one.
-   * Used for initializing the member variable normal_direction.
-   */
-  Point<spacedim> compute_normal(const Tensor<1,spacedim> &vector) const;
-
+protected:
   /**
    * A vector orthogonal to direcetion.
    */
@@ -328,10 +322,19 @@ private:
    */
   const Point<spacedim> point_on_axis;
 
+private:
   /**
    * Relative tolerance to measure zero distances.
    */
   double tolerance;
+
+  // explicitly check for sensible template arguments, but not on windows
+  // because MSVC creates bogus warnings during normal compilation
+#ifndef DEAL_II_MSVC
+  static_assert (spacedim==3,
+                 "CylindricalManifold can only be used for spacedim==3!");
+#endif
+
 };
 
 
index f4d6e35663b623d692036e3b10893cc6bb60764a..3e8dba04a1fb17d7ca04bcd007fa8e26712cb269 100644 (file)
@@ -300,18 +300,49 @@ get_new_point (const std::vector<Point<spacedim> > &vertices,
 // CylindricalManifold
 // ============================================================
 
+namespace
+{
+  // helper function to compute a vector orthogonal to a given one.
+  template <int spacedim>
+  Point<spacedim>
+  compute_normal(const Tensor<1,spacedim> &vector)
+  {
+    AssertThrow(vector.norm() != 0.,
+                ExcMessage("The direction parameter must not be zero!"));
+    Point<3> normal;
+    if (std::abs(vector[0]) >= std::abs(vector[1])
+        && std::abs(vector[0]) >= std::abs(vector[2]))
+      {
+        normal[1]=-1.;
+        normal[2]=-1.;
+        normal[0]=(vector[1]+vector[2])/vector[0];
+      }
+    else if (std::abs(vector[1]) >= std::abs(vector[0])
+             && std::abs(vector[1]) >= std::abs(vector[2]))
+      {
+        normal[0]=-1.;
+        normal[2]=-1.;
+        normal[1]=(vector[0]+vector[2])/vector[1];
+      }
+    else
+      {
+        normal[0]=-1.;
+        normal[1]=-1.;
+        normal[2]=(vector[0]+vector[1])/vector[2];
+      }
+    return normal;
+  }
+}
+
+
+
 template <int dim, int spacedim>
 CylindricalManifold<dim, spacedim>::CylindricalManifold(const unsigned int axis,
                                                         const double tolerance) :
-  ChartManifold<dim,spacedim,3>(Tensor<1,3>({0,2.*numbers::PI,0})),
-              normal_direction(compute_normal(Point<3>::unit_vector(axis))),
-              direction (Point<3>::unit_vector(axis)),
-              point_on_axis (Point<3>()),
-              tolerance(tolerance)
-{
-  static_assert(spacedim==3,
-                "CylindricalManifold can only be used for spacedim==3!");
-}
+  CylindricalManifold<dim, spacedim>(Point<3>::unit_vector(axis),
+                                     Point<3>(),
+                                     tolerance)
+{}
 
 
 
@@ -324,42 +355,7 @@ CylindricalManifold<dim, spacedim>::CylindricalManifold(const Point<spacedim> &d
               direction (direction_/direction_.norm()),
               point_on_axis (point_on_axis_),
               tolerance(tolerance)
-{
-  static_assert(spacedim==3,
-                "CylindricalManifold can only be used for spacedim==3!");
-}
-
-
-
-template <int dim, int spacedim>
-Point<spacedim>
-CylindricalManifold<dim, spacedim>::compute_normal(const Tensor<1,spacedim> &vector) const
-{
-  AssertThrow(vector.norm() != 0.,
-              ExcMessage("The direction parameter must not be zero!"));
-  Point<3> normal;
-  if (std::abs(vector[0]) >= std::abs(vector[1])
-      && std::abs(vector[0]) >= std::abs(vector[2]))
-    {
-      normal[1]=-1.;
-      normal[2]=-1.;
-      normal[0]=(vector[1]+vector[2])/vector[0];
-    }
-  else if (std::abs(vector[1]) >= std::abs(vector[0])
-           && std::abs(vector[1]) >= std::abs(vector[2]))
-    {
-      normal[0]=-1.;
-      normal[2]=-1.;
-      normal[1]=(vector[0]+vector[2])/vector[1];
-    }
-  else
-    {
-      normal[0]=-1.;
-      normal[1]=-1.;
-      normal[2]=(vector[0]+vector[1])/vector[2];
-    }
-  return normal;
-}
+{}
 
 
 
@@ -369,7 +365,7 @@ CylindricalManifold<dim,spacedim>::
 get_new_point (const std::vector<Point<spacedim> > &surrounding_points,
                const std::vector<double>           &weights) const
 {
-  // First check if the the average in space lies on the axis.
+  // First check if the average in space lies on the axis.
   Point<spacedim> middle;
   double average_length = 0.;
   for (unsigned int i=0; i<surrounding_points.size(); ++i)
index dd49a5a2b37a6de7b522f3de5964dfbc37d0f4d7..cf55491a6d29ac027a8eca68ca4c2e09440c51ef 100644 (file)
 #include <string>
 #include <sstream>
 
+template <int dim>
+struct ManifoldWrapper
+{
+  Manifold<dim> *operator()(const Point<dim> &direction,
+                            const Point<dim> &center ) const;
+};
 
+template <>
+Manifold<2> *
+ManifoldWrapper<2>::operator()(const Point<2> &/*direction*/,
+                               const Point<2> &center ) const
+{
+  return new SphericalManifold<2>(center);
+}
+
+template <>
+Manifold<3> *
+ManifoldWrapper<3>::operator()(const Point<3> &direction,
+                               const Point<3> &center ) const
+{
+  return new CylindricalManifold<3>(direction, center);
+}
 
 template <int dim>
 void test()
@@ -45,8 +66,7 @@ void test()
   direction[dim-1] = 1.;
 
   std::shared_ptr<Manifold<dim> > cylinder_manifold
-  (dim == 2 ? static_cast<Manifold<dim>*>(new SphericalManifold<dim>(Point<dim>())) :
-   static_cast<Manifold<dim>*>(new CylindricalManifold<dim>(direction, Point<dim>())));
+  (ManifoldWrapper<dim>()(direction, Point<dim>()));
 
   // create cube and shift it to position 0.7
   Triangulation<dim> tria;
index ec551b0559edee5dec073768252dd2b505951d51..9f525c88bd120cd529cd1c72f34e0f3565104333 100644 (file)
@@ -150,7 +150,28 @@ void create_triangulation(Triangulation<3> &tria)
     }
 }
 
+template <int dim>
+struct ManifoldWrapper
+{
+  Manifold<dim> *operator()(const Point<dim> &direction,
+                            const Point<dim> &center ) const;
+};
+
+template <>
+Manifold<2> *
+ManifoldWrapper<2>::operator()(const Point<2> &/*direction*/,
+                               const Point<2> &center ) const
+{
+  return new SphericalManifold<2>(center);
+}
 
+template <>
+Manifold<3> *
+ManifoldWrapper<3>::operator()(const Point<3> &direction,
+                               const Point<3> &center ) const
+{
+  return new CylindricalManifold<3>(direction, center);
+}
 
 template <int dim>
 void test()
@@ -162,8 +183,7 @@ void test()
   direction[dim-1] = 1.;
 
   std::shared_ptr<Manifold<dim> > cylinder_manifold
-  (dim == 2 ? static_cast<Manifold<dim>*>(new SphericalManifold<dim>(center)) :
-   static_cast<Manifold<dim>*>(new CylindricalManifold<dim>(direction, center)));
+  (ManifoldWrapper<dim>()(direction, center));
   Triangulation<dim> tria;
   create_triangulation(tria);
   tria.set_manifold(MANIFOLD_ID, *cylinder_manifold);

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.