]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Dummy implementation for CylindricalManifold for spacedim!=3 5695/head
authorDaniel Arndt <daniel.arndt@iwr.uni-heidelberg.de>
Tue, 2 Jan 2018 23:38:34 +0000 (00:38 +0100)
committerDaniel Arndt <daniel.arndt@iwr.uni-heidelberg.de>
Sun, 7 Jan 2018 12:58:19 +0000 (13:58 +0100)
include/deal.II/grid/manifold_lib.h
source/grid/manifold_lib.cc
source/grid/manifold_lib.inst.in

index 27f44282f15b226375d23a98a336137bbbbd2984..7686a0d7d0c554b2d985976498514fee43d2e69e 100644 (file)
@@ -920,6 +920,92 @@ private:
   FlatManifold<dim> chart_manifold;
 };
 
+#ifndef DOXYGEN
+/**
+ * Specialization for the only properly implemented spacedim parameter.
+ */
+template <int dim>
+class CylindricalManifold<dim,3> : public ChartManifold<dim,3,3>
+{
+public:
+  /**
+   * Constructor. Using default values for the constructor arguments yields a
+   * cylinder along the x-axis (<tt>axis=0</tt>). Choose <tt>axis=1</tt> or
+   * <tt>axis=2</tt> for a tube along the y- or z-axis, respectively. The
+   * tolerance value is used to determine if a point is on the axis.
+   */
+  CylindricalManifold (const unsigned int axis = 0,
+                       const double       tolerance = 1e-10);
+
+  /**
+   * Constructor. If constructed with this constructor, the manifold described
+   * is a cylinder with an axis that points in direction #direction and goes
+   * through the given #point_on_axis. The direction may be arbitrarily
+   * scaled, and the given point may be any point on the axis. The tolerance
+   * value is used to determine if a point is on the axis.
+   */
+  CylindricalManifold (const Point<3> &direction,
+                       const Point<3> &point_on_axis,
+                       const double    tolerance = 1e-10);
+
+  /**
+   * Compute the Cartesian coordinates for a point given in cylindrical
+   * coordinates.
+   */
+  virtual Point<3>
+  pull_back(const Point<3> &space_point) const override;
+
+  /**
+   * 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<3>
+  push_forward(const Point<3> &chart_point) const override;
+
+  /**
+   * Compute the derivatives of the mapping from cylindrical coordinates
+   * $(r, \phi, \lambda)$ to cartesian coordinates 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 DerivativeForm<1, 3, 3>
+  push_forward_gradient(const Point<3> &chart_point) const override;
+
+  /**
+   * Compute new points on the CylindricalManifold. See the documentation of
+   * the base class for a detailed description of what this function does.
+   */
+  virtual Point<3>
+  get_new_point (const ArrayView<const Point<3>> &surrounding_points,
+                 const ArrayView<const double>   &weights) const override;
+
+protected:
+  /**
+   * A vector orthogonal to the normal direction.
+   */
+  const Tensor<1,3> normal_direction;
+
+  /**
+   * The direction vector of the axis.
+   */
+
+  const Tensor<1,3> direction;
+  /**
+   * An arbitrary point on the axis.
+   */
+  const Point<3> point_on_axis;
+
+private:
+  /**
+   * Relative tolerance to measure zero distances.
+   */
+  double tolerance;
+
+};
+#endif //DOXYGEN
+
 DEAL_II_NAMESPACE_CLOSE
 
 #endif
index 52e836103812ad3558e8bd3db987e5d17b446f7a..6840286dcd819a84f2f898970b7df58aed6bf73f 100644 (file)
@@ -813,19 +813,29 @@ get_new_point (const ArrayView<const Tensor<1,3>> &directions,
 template <int dim, int spacedim>
 CylindricalManifold<dim, spacedim>::CylindricalManifold(const unsigned int axis,
                                                         const double tolerance) :
-  CylindricalManifold<dim, spacedim>(Point<3>::unit_vector(axis),
-                                     Point<3>(),
+  CylindricalManifold<dim, spacedim>(Point<spacedim>::unit_vector(axis),
+                                     Point<spacedim>(),
                                      tolerance)
 {}
 
 
 
+template <int dim>
+CylindricalManifold<dim, 3>::CylindricalManifold(const unsigned int axis,
+                                                 const double tolerance) :
+  CylindricalManifold<dim, 3>(Point<3>::unit_vector(axis),
+                              Point<3>(),
+                              tolerance)
+{}
+
+
+
 template <int dim, int spacedim>
 CylindricalManifold<dim, spacedim>::CylindricalManifold(const Point<spacedim> &direction_,
                                                         const Point<spacedim> &point_on_axis_,
                                                         const double tolerance) :
-  ChartManifold<dim,3,3>(Tensor<1,3>({0,2.*numbers::PI,0})),
-              normal_direction(internal::compute_normal(direction_, true)),
+  ChartManifold<dim,spacedim,3>(Tensor<1,3>({0,2.*numbers::PI,0})),
+              normal_direction(Tensor<1, spacedim>()),
               direction (direction_/direction_.norm()),
               point_on_axis (point_on_axis_),
               tolerance(tolerance)
@@ -838,14 +848,40 @@ CylindricalManifold<dim, spacedim>::CylindricalManifold(const Point<spacedim> &d
 
 
 
+template <int dim>
+CylindricalManifold<dim, 3>::CylindricalManifold(const Point<3> &direction_,
+                                                 const Point<3> &point_on_axis_,
+                                                 const double tolerance) :
+  ChartManifold<dim,3,3>(Tensor<1,3>({0,2.*numbers::PI,0})),
+              normal_direction(internal::compute_normal(direction_, true)),
+              direction (direction_/direction_.norm()),
+              point_on_axis (point_on_axis_),
+              tolerance(tolerance)
+{}
+
+
+
 template <int dim, int spacedim>
 Point<spacedim>
 CylindricalManifold<dim,spacedim>::
 get_new_point (const ArrayView<const Point<spacedim>> &surrounding_points,
                const ArrayView<const double>          &weights) const
+{
+  Assert (spacedim==3,
+          ExcMessage("CylindricalManifold can only be used for spacedim==3!"));
+  return Point<spacedim>();
+}
+
+
+
+template <int dim>
+Point<3>
+CylindricalManifold<dim,3>::
+get_new_point (const ArrayView<const Point<3>> &surrounding_points,
+               const ArrayView<const double>          &weights) const
 {
   // First check if the average in space lies on the axis.
-  Point<spacedim> middle;
+  Point<3> middle;
   double average_length = 0.;
   for (unsigned int i=0; i<surrounding_points.size(); ++i)
     {
@@ -856,10 +892,10 @@ get_new_point (const ArrayView<const Point<spacedim>> &surrounding_points,
   const double lambda = middle*direction;
 
   if ((middle-direction*lambda).square() < tolerance*average_length)
-    return Point<spacedim>()+direction*lambda;
+    return Point<3>()+direction*lambda;
   else // If not, using the ChartManifold should yield valid results.
-    return ChartManifold<dim, spacedim, 3>::get_new_point(surrounding_points,
-                                                          weights);
+    return ChartManifold<dim, 3, 3>::get_new_point(surrounding_points,
+                                                   weights);
 }
 
 
@@ -867,6 +903,17 @@ get_new_point (const ArrayView<const Point<spacedim>> &surrounding_points,
 template <int dim, int spacedim>
 Point<3>
 CylindricalManifold<dim, spacedim>::pull_back(const Point<spacedim> &space_point) const
+{
+  Assert (spacedim==3,
+          ExcMessage("CylindricalManifold can only be used for spacedim==3!"));
+  return Point<3>();
+}
+
+
+
+template <int dim>
+Point<3>
+CylindricalManifold<dim, 3>::pull_back(const Point<3> &space_point) const
 {
   // First find the projection of the given point to the axis.
   const Tensor<1,3> normalized_point = space_point-point_on_axis;
@@ -889,6 +936,17 @@ CylindricalManifold<dim, spacedim>::pull_back(const Point<spacedim> &space_point
 template <int dim, int spacedim>
 Point<spacedim>
 CylindricalManifold<dim, spacedim>::push_forward(const Point<3> &chart_point) const
+{
+  Assert (spacedim==3,
+          ExcMessage("CylindricalManifold can only be used for spacedim==3!"));
+  return Point<spacedim>();
+}
+
+
+
+template<int dim>
+Point<3>
+CylindricalManifold<dim, 3>::push_forward(const Point<3> &chart_point) const
 {
   // Rotate the orthogonal direction by the given angle.
   // Formula from Section 5.2 in
@@ -909,6 +967,18 @@ CylindricalManifold<dim, spacedim>::push_forward(const Point<3> &chart_point) co
 template<int dim, int spacedim>
 DerivativeForm<1, 3, spacedim>
 CylindricalManifold<dim, spacedim>::push_forward_gradient(const Point<3> &chart_point) const
+{
+  Assert (spacedim==3,
+          ExcMessage("CylindricalManifold can only be used for spacedim==3!"));
+
+  return DerivativeForm<1,3,spacedim>();
+}
+
+
+
+template<int dim>
+DerivativeForm<1, 3, 3>
+CylindricalManifold<dim, 3>::push_forward_gradient(const Point<3> &chart_point) const
 {
   Tensor<2, 3> derivatives;
 
index d9d00a7c7b77045286cb9e80cc7b4265e682f451..7e181074c76440c8015057da1dfe9443b15f2cd6 100644 (file)
 //
 // ---------------------------------------------------------------------
 
-for (deal_II_dimension : DIMENSIONS)
-{
-    template class CylindricalManifold<deal_II_dimension, 3>;
-}
-
 for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension :  SPACE_DIMENSIONS)
 {
 #if deal_II_dimension <= deal_II_space_dimension
     template class PolarManifold<deal_II_dimension, deal_II_space_dimension>;
     template class SphericalManifold<deal_II_dimension, deal_II_space_dimension>;
     template class TransfiniteInterpolationManifold<deal_II_dimension, deal_II_space_dimension>;
+    template class CylindricalManifold<deal_II_dimension, deal_II_space_dimension>;
 #endif
 #if deal_II_dimension == deal_II_space_dimension
     template class TorusManifold<deal_II_dimension>;

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.