]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Transfinite interpolation manifold.
authorMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Thu, 18 May 2017 16:49:10 +0000 (18:49 +0200)
committerMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Thu, 1 Jun 2017 14:13:39 +0000 (16:13 +0200)
Added tests.
Adjusted division by zero in MappingQGeneric.
Disable using additional points on TransfiniteInterpolationManifold

14 files changed:
doc/doxygen/images/circular_mesh_only_boundary_manifold.png [new file with mode: 0644]
doc/doxygen/images/circular_mesh_transfinite_interpolation.png [new file with mode: 0644]
include/deal.II/grid/manifold_lib.h
source/fe/mapping_q_generic.cc
source/grid/manifold_lib.cc
source/grid/manifold_lib.inst.in
tests/manifold/transfinite_manifold_01.cc [new file with mode: 0644]
tests/manifold/transfinite_manifold_01.output [new file with mode: 0644]
tests/manifold/transfinite_manifold_02.cc [new file with mode: 0644]
tests/manifold/transfinite_manifold_02.output [new file with mode: 0644]
tests/manifold/transfinite_manifold_03.cc [new file with mode: 0644]
tests/manifold/transfinite_manifold_03.output [new file with mode: 0644]
tests/manifold/transfinite_manifold_04.cc [new file with mode: 0644]
tests/manifold/transfinite_manifold_04.output [new file with mode: 0644]

diff --git a/doc/doxygen/images/circular_mesh_only_boundary_manifold.png b/doc/doxygen/images/circular_mesh_only_boundary_manifold.png
new file mode 100644 (file)
index 0000000..fb502d5
Binary files /dev/null and b/doc/doxygen/images/circular_mesh_only_boundary_manifold.png differ
diff --git a/doc/doxygen/images/circular_mesh_transfinite_interpolation.png b/doc/doxygen/images/circular_mesh_transfinite_interpolation.png
new file mode 100644 (file)
index 0000000..c0ace94
Binary files /dev/null and b/doc/doxygen/images/circular_mesh_transfinite_interpolation.png differ
index cc7d5e757c9b4efce6ea5484041599d6f602d41d..f6ab6b14932227b276972dcb04d992422180d8cd 100644 (file)
@@ -503,6 +503,282 @@ private:
   double r, R;
 };
 
+
+
+/**
+ * A mapping class that extends curved boundary descriptions into the interior
+ * of the computational domain. The outer curved boundary description is
+ * assumed to be given by another manifold (e.g. a polar manifold on a circle).
+ * The mechanism to extend the boundary information is a so-called transfinite
+ * interpolation.
+ *
+ * The formula for extending such a description in 2D is inspired by
+ * <a href="https://en.wikipedia.org/wiki/Transfinite_interpolation">
+ * Wikipedia</a>
+ * Given a point $(u,v)$ on the chart, the image of this point in real space
+ * is given by
+ * @f{align*}{
+ * \bf S(u,v) &= (1-v)\bf c_0(u)+v \bf c_1(u) + (1-u)\bf c_2(v) + u \bf c_3(v) \\
+ * &\quad - \begin[(1-u)(1-v) \bf x_0 + u(1-v) \bf x_1 + (1-u)v \bf x_2 + uv \bf x_3 \right]
+ * @f}
+ * where $\bf x_0, \bf x_1, \bf x_2, \bf x_3$ denote the four bounding vertices
+ * bounding the image space and $\bf c_0, \bf c_1, \bf c_2, \bf c_3$ are the
+ * four curves describing the lines of the cell. If a curved manifold is
+ * attached to any of these lines, the evaluation is done according to
+ * Manifold::get_new_point() with the two end points of the line and
+ * appropriate weight. In 3D, the generalization of this formula is
+ * implemented, creating a weighted sum of the vertices (positive
+ * contribution), the lines (negative), and the faces (positive contribution).
+ *
+ * This manifold is usually attached to a coarse mesh and then places new
+ * points as a combination of the descriptions on the boundaries, weighted
+ * appropriately according to the position of the point on the original chart
+ * point $(u,v). Whenever possible, this manifold should be preferred over
+ * setting only a curved manifold on the boundary of a mesh, since the latter
+ * will need to switch from a curved description to a straight description in a
+ * single layer of elements, which causes an error order on cells close to the
+ * boundary that does not exceed 3 no matter how high the degree of the
+ * polynomial mapping and the finite element space. Using this class instead,
+ * the curved nature of the manifold that is originally contained in one
+ * <i>coarse</i> mesh layer will be applied to more than one <i>fine</i> mesh
+ * layer once the mesh gets refined, restoring the optimal convergence rates of
+ * the underlying finite element and mapping as appropriate.
+ *
+ * If no curved boundaries surround a coarse cell, this class reduces to a flat
+ * manifold description.
+ *
+ * To give an example of using this class, the following code attaches a
+ * transfinite manifold to a circle:
+ *
+ * @code
+ * PolarManifold<dim> polar_manifold;
+ * TransfiniteInterpolationManifold<dim> inner_manifold;
+ *
+ * Triangulation<dim> triangulation;
+ * GridGenerator::hyper_ball (triangulation);
+ *
+ * triangulation.set_all_manifold_ids(1);
+ * triangulation.set_all_manifold_ids_on_boundary(0);
+ * triangulation.set_manifold (0, polar_manifold);
+ * inner_manifold.initialize(triangulation);
+ * triangulation.set_manifold (1, inner_manifold);
+ * triangulation.refine_global(4);
+ * @endcode
+ *
+ * In this code, we first set all manifold ids to the id of the transfinite
+ * interpolation, and then re-set the manifold ids on the boundary to identify
+ * the curved boundary described by the polar manifold. With this code, one
+ * gets a really nice mesh:
+ *
+ * <p ALIGN="center">
+ * @image html circular_mesh_transfinite_interpolation.png
+ * </p>
+ *
+ * which is obviously much nicer than the polar manifold applied to just the
+ * boundary:
+ *
+ * <p ALIGN="center">
+ * @image html circular_mesh_only_boundary_manifold.png
+ * </p>
+ *
+ * @note In the implementation of this class, the manifolds surrounding a
+ * coarse cell are queried repeatedly to compute points on their interior. For
+ * optimal mesh quality, those manifolds should be compatible with a chart
+ * notion. For example, computing a point that is 0.25 along the line between
+ * two vertices using the weights 0.25 and 0.75 for the two vertices should
+ * give the same result as first computing the mid point at 0.5 and then again
+ * compute the midpoint between the first vertex and coarse mid point. This is
+ * the case for PolarManifold but not for Spherical manifold, so be careful
+ * when using the latter. In case the quality of the manifold is not good
+ * enough, upon mesh refinement it may happen that the transformation to a
+ * chart inside the get_new_point() or add_new_point() methods produces points
+ * that are outside the unit cells, then this class throws an exception of type
+ * Manifold@<dim,spacedim@>::ExcTransformationFailed. If that happens, the mesh
+ * should be refined before applying this class, as done in the following
+ * example:
+ * @code
+ * SphericalManifold<dim> spherical_manifold;
+ * TransfiniteInterpolationManifold<dim> inner_manifold;
+ * Triangulation<dim> triangulation;
+ * GridGenerator::hyper_ball (triangulation);
+ * triangulation.set_all_manifold_ids(1);
+ * triangulation.set_all_manifold_ids_on_boundary(0);
+ * triangulation.set_manifold (0, polar_manifold);
+ * inner_manifold.initialize(triangulation);
+ * triangulation.set_manifold (1, inner_manifold);
+ * triangulation.refine_global(1);
+ * // initialize the transfinite manifold again
+ * inner_manifold.initialize(triangulation);
+ * triangulation.refine_global(4);
+ * @endcode
+ *
+ * @note For performance and accuracy reasons, it is recommended to apply the
+ * transfinite manifold to as coarse a mesh as possible. Regarding accuracy,
+ * the curved description can only be applied to new points created from a
+ * given neighborhood, and the grid quality is typically higher when extending
+ * the curved description over as large a domain as possible. Regarding
+ * performance, the identification of the correct coarse cell in the
+ * get_new_point() method needs to pass all coarse cells, so expect a quadratic
+ * complexity in the number of coarse cells. Thus, the current implementation
+ * is only really economical when there are not more than a few hundred coarse
+ * cells. To make performance better for larger numbers of cells, one could
+ * extend the current implementation by a pre-identification of relevant cells
+ * with axis-aligned bounding boxes.
+ *
+ * @author Martin Kronbichler, Luca Heltai, 2017
+ */
+template <int dim, int spacedim=dim>
+class TransfiniteInterpolationManifold : public Manifold<dim,spacedim>
+{
+public:
+  /**
+   * Constructor.
+   */
+  TransfiniteInterpolationManifold();
+
+  /**
+   * Initializes the manifold with a coarse mesh. The prerequisite for using
+   * this class is that the input triangulation is uniformly refined and the
+   * manifold is later attached to the same triangulation.
+   *
+   * @note The triangulation used to construct the manifold must not be
+   * destroyed during the usage of this object.
+   */
+  void initialize (const Triangulation<dim,spacedim> &triangulation);
+
+  /**
+   * Return the point which shall become the new vertex surrounded by the
+   * given points @p surrounding_points. @p weights contains appropriate
+   * weights for the surrounding points according to which the manifold
+   * determines the new point's position.
+   *
+   * The implementation in this class overrides the method in the base class
+   * and computes the new point by a transfinite interpolation. The first step
+   * in the implementation is to identify the coarse cell on which the
+   * surrounding points are located. Then, the coordinates are transformed to
+   * the unit coordinates on the coarse cell by a Newton iteration, where the
+   * new point is then computed according to the weights. Finally, it is
+   * pushed forward to the real space according to the transfinite
+   * interpolation.
+   */
+  virtual
+  Point<spacedim>
+  get_new_point (const std::vector<Point<spacedim> > &surrounding_points,
+                 const std::vector<double>           &weights) const;
+
+  /**
+   * Compute a new set of points that interpolate between the given points
+   * @p surrounding_points. @p weights is a table with as many columns as
+   * @p surrounding_points.size(). The number of rows in @p weights determines
+   * how many new points will be computed and appended to the last input
+   * argument @p new_points. After exit of this function, the size of
+   * @p new_points equals the size at entry plus the number of rows in
+   * @p weights.
+   *
+   * The implementation in this class overrides the method in the base class
+   * and computes the new point by a transfinite interpolation. The first step
+   * in the implementation is to identify the coarse cell on which the
+   * surrounding points are located. Then, the coordinates are transformed to
+   * the unit coordinates on the coarse cell by a Newton iteration, where the
+   * new points are then computed according to the weights. Finally, the is
+   * pushed forward to the real space according to the transfinite
+   * interpolation.
+   *
+   * The implementation does not allow for @p surrounding_points and
+   * @p new_points to point to the same vector, so make sure to pass different
+   * objects into the function.
+   */
+  virtual
+  void
+  add_new_points (const std::vector<Point<spacedim> > &surrounding_points,
+                  const Table<2,double>               &weights,
+                  std::vector<Point<spacedim> >       &new_points) const;
+
+private:
+  /**
+   * Internal function to identify the most suitable cells (=charts) where the
+   * given surrounding points are located. We use a cheap algorithm to identify
+   * the cells and rank the cells by probability before we actually do the
+   * search inside the relevant cells. The cells are sorted by the distance of
+   * a Q1 approximation of the inverse mapping to the unit cell of the
+   * surrounding points. We expect at most 10 cells (it should be less than 8
+   * candidates even in 3D, typically only two or three), so get an array with
+   * 10 entries of a the indices <tt>cell->index()</tt>.
+   */
+  std::array<unsigned int, 10>
+  get_possible_cells_around_points(const std::vector<Point<spacedim> > &surrounding_points) const;
+
+  /**
+   * Finalizes the identification of the correct chart and returns the location
+   * of the surrounding points on the chart. This method internally calls
+   * @p get_possible_cells_around_points().
+   */
+  std::pair<typename Triangulation<dim,spacedim>::cell_iterator,
+      std::vector<Point<dim> > >
+      compute_chart_points(const std::vector<Point<spacedim> > &surrounding_points) const;
+
+  /**
+   * Pull back operation into the unit coordinates on the given coarse cell.
+   *
+   * @note This internal function is currently not compatible with the
+   * ChartManifold::pull_back() function because the given class represents an
+   * atlas of charts, not a single chart. Thus, the pull_back() operation is
+   * only valid with the additional information of the chart, given by a cell
+   * on the coarse grid. An alternative implementation could shift the index
+   * depending on the coarse cell for a 1-to-1 relation between the chart space
+   * and the image space.
+   */
+  Point<dim>
+  pull_back(const typename Triangulation<dim,spacedim>::cell_iterator &cell,
+            const Point<spacedim> &p) const;
+
+  /**
+   * Push forward operation.
+   *
+   * @note This internal function is currently not compatible with the
+   * ChartManifold::pull_back() function because the given class represents an
+   * atlas of charts, not a single chart. Thus, the pull_back() operation is
+   * only valid with the additional information of the chart, given by a cell
+   * on the coarse grid. An alternative implementation could shift the index
+   * depending on the coarse cell for a 1-to-1 relation between the chart space
+   * and the image space.
+   */
+  Point<spacedim>
+  push_forward(const typename Triangulation<dim,spacedim>::cell_iterator &cell,
+               const Point<dim> &chart_point) const;
+
+  /**
+   * Gradient of the push_forward method.
+   */
+  DerivativeForm<1,dim,spacedim>
+  push_forward_gradient(const typename Triangulation<dim,spacedim>::cell_iterator &cell,
+                        const Point<dim> &chart_point) const;
+
+  /**
+   * The underlying triangulation.
+   */
+  const Triangulation<dim,spacedim> *triangulation;
+
+  /**
+   * The level of the mesh cells where the transfinite approximation is
+   * applied, usually level 0.
+   */
+  int level_coarse;
+
+  /**
+   * In case there all surrounding manifolds are the transfinite manifold or
+   * have default (invalid) manifold id, the manifold degenerates to a flat
+   * manifold and we can choose cheaper algorithms for the push_forward method.
+   */
+  std::vector<bool> coarse_cell_is_flat;
+
+  /**
+   * A flat manifold used to compute new points in the chart space where it we
+   * use a FlatManifold description.
+   */
+  FlatManifold<dim> chart_manifold;
+};
+
 DEAL_II_NAMESPACE_CLOSE
 
 #endif
index 2eb5cf2f6da4e877d64c6ce0ebb2142afdeebec2..53382861240df084c8eaa3e4eaae47328b6eebd5 100644 (file)
@@ -24,6 +24,7 @@
 #include <deal.II/grid/tria.h>
 #include <deal.II/grid/tria_iterator.h>
 #include <deal.II/grid/tria_boundary.h>
+#include <deal.II/grid/manifold_lib.h>
 #include <deal.II/dofs/dof_accessor.h>
 #include <deal.II/fe/fe_tools.h>
 #include <deal.II/fe/fe.h>
@@ -1461,6 +1462,8 @@ namespace
           }
 
         // Solve  [f'(x)]d=f(x)
+        AssertThrow(determinant(df) > 0,
+                    (typename Mapping<dim,spacedim>::ExcTransformationFailed()));
         Tensor<2,spacedim> df_inverse = invert(df);
         const Tensor<1,spacedim> delta = df_inverse * static_cast<const Tensor<1,spacedim>&>(f);
 
@@ -3792,21 +3795,23 @@ compute_mapping_support_points(const typename Triangulation<dim,spacedim>::cell_
   if (this->polynomial_degree > 1)
     {
       // check if all entities have the same manifold id which is when we can
-      // simply ask the manifold for all points
+      // simply ask the manifold for all points. the transfinite manifold can
+      // do the interpolation better than this class, so if we detect that we
+      // do not have to change anything here
       Assert(dim<=3, ExcImpossibleInDim(dim));
       bool all_manifold_ids_are_equal = (dim == spacedim);
-      if (dim > 1)
-        for (unsigned int l=0; l<GeometryInfo<dim>::lines_per_cell; ++l)
-          if (&cell->line(l)->get_manifold() != &cell->get_manifold())
-            all_manifold_ids_are_equal = false;
-      if (dim == 3)
-        if (check_identical_manifolds_of_quads(cell) == false)
-          all_manifold_ids_are_equal = false;
-      if (all_manifold_ids_are_equal)
+      if (all_manifold_ids_are_equal &&
+          dynamic_cast<const TransfiniteInterpolationManifold<dim,spacedim>*>(&cell->get_manifold()) == nullptr)
         {
-          std::vector<Point<spacedim> > vertices(a);
-          cell->get_manifold().add_new_points(vertices, support_point_weights_cell, a);
+          if (dim > 1)
+            for (unsigned int l=0; l<GeometryInfo<dim>::lines_per_cell; ++l)
+              if (&cell->line(l)->get_manifold() != &cell->get_manifold())
+                all_manifold_ids_are_equal = false;
+          if (dim == 3)
+            if (check_identical_manifolds_of_quads(cell) == false)
+              all_manifold_ids_are_equal = false;
         }
+
       else
         switch (dim)
           {
index 46c46eb371d4aaab87ad9fddc3e2538b947123ab..1de6dbb3335d5417d7d919357129274d1e730aa5 100644 (file)
 //
 // ---------------------------------------------------------------------
 
+#include <deal.II/base/tensor.h>
+#include <deal.II/base/table.h>
 #include <deal.II/grid/tria.h>
 #include <deal.II/grid/tria_iterator.h>
 #include <deal.II/grid/tria_accessor.h>
 #include <deal.II/grid/manifold_lib.h>
-#include <deal.II/base/tensor.h>
 #include <deal.II/lac/vector.h>
+#include <deal.II/fe/mapping_q1.h>
 #include <cmath>
 
 DEAL_II_NAMESPACE_OPEN
@@ -549,6 +551,543 @@ TorusManifold<dim>::push_forward_gradient(const Point<3> &chart_point) const
 
 
 
+template <int dim, int spacedim>
+TransfiniteInterpolationManifold<dim,spacedim>::TransfiniteInterpolationManifold()
+  :
+  triangulation(nullptr),
+  level_coarse (-1)
+{}
+
+
+
+template <int dim, int spacedim>
+void
+TransfiniteInterpolationManifold<dim,spacedim>
+::initialize(const Triangulation<dim,spacedim> &triangulation)
+{
+  this->triangulation = &triangulation;
+  // in case the triangulatoin is cleared, remove the pointers by a signal
+  triangulation.signals.clear.connect
+  ([&]() -> void {this->triangulation = nullptr; this->level_coarse = -1;});
+  level_coarse = triangulation.last()->level();
+  coarse_cell_is_flat.resize(triangulation.n_cells(level_coarse), false);
+  typename Triangulation<dim,spacedim>::active_cell_iterator
+  cell = triangulation.begin(level_coarse),
+  endc = triangulation.end(level_coarse);
+  for ( ; cell != endc; ++cell)
+    {
+      bool cell_is_flat = true;
+      for (unsigned int l=0; l<GeometryInfo<dim>::lines_per_cell; ++l)
+        if (cell->line(l)->manifold_id() != cell->manifold_id() &&
+            cell->line(l)->manifold_id() != numbers::invalid_manifold_id)
+          cell_is_flat = false;
+      if (dim > 2)
+        for (unsigned int q=0; q<GeometryInfo<dim>::quads_per_cell; ++q)
+          if (cell->quad(q)->manifold_id() != cell->manifold_id() &&
+              cell->quad(q)->manifold_id() != numbers::invalid_manifold_id)
+            cell_is_flat = false;
+      AssertIndexRange(static_cast<unsigned int>(cell->index()),
+                       coarse_cell_is_flat.size());
+      coarse_cell_is_flat[cell->index()] = cell_is_flat;
+    }
+}
+
+
+
+namespace
+{
+  // version for 1D
+  template <typename AccessorType>
+  Point<AccessorType::space_dimension>
+  compute_transfinite_interpolation(const AccessorType &cell,
+                                    const Point<1> &chart_point,
+                                    const bool      cell_is_flat)
+  {
+    return cell.vertex(0) * (1.-chart_point[0]) + cell.vertex(1) * chart_point[0];
+  }
+
+  // version for 2D
+  template <typename AccessorType>
+  Point<AccessorType::space_dimension>
+  compute_transfinite_interpolation(const AccessorType &cell,
+                                    const Point<2>     &chart_point,
+                                    const bool          cell_is_flat)
+  {
+    const unsigned int spacedim = AccessorType::space_dimension;
+    const types::manifold_id my_manifold_id = cell.manifold_id();
+
+    // formula see wikipedia
+    // https://en.wikipedia.org/wiki/Transfinite_interpolation
+    // S(u,v) = (1-v)c_1(u)+v c_3(u) + (1-u)c_2(v) + u c_4(v) -
+    //   [(1-u)(1-v)P_0 + u(1-v) P_1 + (1-u)v P_2 + uv P_3]
+
+    Point<spacedim> new_point;
+    if (cell_is_flat)
+      for (unsigned int v=0; v<GeometryInfo<2>::vertices_per_cell; ++v)
+        new_point += GeometryInfo<2>::d_linear_shape_function(chart_point, v) *
+                     cell.vertex(v);
+    else
+      {
+        // We subtract the contribution of the vertices (second line in formula).
+        // If a line applies the same manifold as the cell, we also subtract a
+        // weighted part of the vertex, so accumulate the final weight of the
+        // the vertices while going through the faces (this is a bit artificial
+        // in 2D but it becomes clear in 3D where we avoid looking at the faces'
+        // orientation and other complications).
+        double weights_vertices[GeometryInfo<2>::vertices_per_cell];
+        for (unsigned int v=0; v<GeometryInfo<2>::vertices_per_cell; ++v)
+          weights_vertices[v] = -GeometryInfo<2>::d_linear_shape_function(chart_point, v);
+
+        // add the contribution from the lines around the cell (first line in
+        // formula)
+        std::vector<double> weights(GeometryInfo<2>::vertices_per_face);
+        std::vector<Point<spacedim> > points(GeometryInfo<2>::vertices_per_face);
+        for (unsigned int line=0; line<GeometryInfo<2>::lines_per_cell; ++line)
+          {
+            const double my_weight = line%2 ? chart_point[line/2] : 1-chart_point[line/2];
+            const double line_point = chart_point[1-line/2];
+
+            // Same manifold or invalid id which will go back to the same class ->
+            // adds to the vertices
+            if (cell.line(line)->manifold_id() == my_manifold_id ||
+                cell.line(line)->manifold_id() == numbers::invalid_manifold_id)
+              {
+                weights_vertices[GeometryInfo<2>::line_to_cell_vertices(line,0)]
+                += my_weight * (1.-line_point);
+                weights_vertices[GeometryInfo<2>::line_to_cell_vertices(line,1)]
+                += my_weight * line_point;
+              }
+            else
+              {
+                points[0] = cell.vertex(GeometryInfo<2>::line_to_cell_vertices(line,0));
+                points[1] = cell.vertex(GeometryInfo<2>::line_to_cell_vertices(line,1));
+                weights[0] = 1. - line_point;
+                weights[1] = line_point;
+                new_point += my_weight *
+                             cell.line(line)->get_manifold().get_new_point(points, weights);
+              }
+          }
+
+        // subtract contribution from the vertices (second line in formula)
+        for (unsigned int v=0; v<GeometryInfo<2>::vertices_per_cell; ++v)
+          new_point += weights_vertices[v] * cell.vertex(v);
+      }
+
+    return new_point;
+  }
+
+  // version for 3D
+  template <typename AccessorType>
+  Point<AccessorType::space_dimension>
+  compute_transfinite_interpolation(const AccessorType &cell,
+                                    const Point<3> &chart_point,
+                                    const bool      cell_is_flat)
+  {
+    const unsigned int dim = AccessorType::dimension;
+    const unsigned int spacedim = AccessorType::space_dimension;
+    const types::manifold_id my_manifold_id = cell.manifold_id();
+
+    // Same approach as in 2D, but adding the faces, subtracting the edges, and
+    // adding the vertices
+    Point<spacedim> new_point;
+
+    if (cell_is_flat)
+      for (unsigned int v=0; v<GeometryInfo<3>::vertices_per_cell; ++v)
+        new_point += GeometryInfo<3>::d_linear_shape_function(chart_point, v) *
+                     cell.vertex(v);
+    else
+      {
+        // identify the weights for the vertices and lines to be accumulated
+        double weights_vertices[GeometryInfo<3>::vertices_per_cell];
+        for (unsigned int v=0; v<GeometryInfo<3>::vertices_per_cell; ++v)
+          weights_vertices[v] = GeometryInfo<3>::d_linear_shape_function(chart_point, v);
+        double weights_lines[GeometryInfo<3>::lines_per_cell];
+        for (unsigned int line=0; line<GeometryInfo<3>::lines_per_cell; ++line)
+          weights_lines[line] = 0;
+
+        // start with the contributions of the faces
+        std::vector<double> weights;
+        std::vector<Point<spacedim> > points;
+        for (unsigned int face=0; face<GeometryInfo<3>::faces_per_cell; ++face)
+          {
+            Point<2> quad_point(chart_point[(face/2+1)%3], chart_point[(face/2+2)%3]);
+            const double my_weight = face%2 ? chart_point[face/2] : 1-chart_point[face/2];
+
+            if (std::abs(my_weight) < 1e-13)
+              continue;
+
+            // same manifold or invalid id which will go back to the same class
+            // -> face will interpolate from the surrounding lines and vertices
+            if (cell.face(face)->manifold_id() == my_manifold_id ||
+                cell.face(face)->manifold_id() == numbers::invalid_manifold_id)
+              {
+                for (unsigned int line=0; line<GeometryInfo<2>::lines_per_cell; ++line)
+                  {
+                    const double line_weight = line%2 ? quad_point[line/2] : 1-quad_point[line/2];
+                    weights_lines[GeometryInfo<3>::face_to_cell_lines(face, line)] +=
+                      my_weight * line_weight;
+                  }
+                for (unsigned int v=0; v<GeometryInfo<2>::vertices_per_cell; ++v)
+                  weights_vertices[GeometryInfo<3>::face_to_cell_vertices(face, v)]
+                  -= GeometryInfo<2>::d_linear_shape_function(quad_point, v) * my_weight;
+              }
+            else
+              {
+                points.resize(GeometryInfo<2>::vertices_per_cell);
+                weights.resize(GeometryInfo<2>::vertices_per_cell);
+                for (unsigned int v=0; v<GeometryInfo<2>::vertices_per_cell; ++v)
+                  {
+                    points[v] = cell.vertex(GeometryInfo<3>::face_to_cell_vertices(face,v));
+                    weights[v] = GeometryInfo<2>::d_linear_shape_function(quad_point, v);
+                  }
+                new_point += my_weight *
+                             cell.face(face)->get_manifold().get_new_point(points, weights);
+              }
+          }
+
+        // next subtract the contributions of the lines
+        for (unsigned int line=0; line<GeometryInfo<3>::lines_per_cell; ++line)
+          {
+            const double line_point = (line < 8 ? chart_point[1-(line%4)/2] : chart_point[2]);
+            double my_weight = 0.;
+            if (line < 8)
+              {
+                const unsigned int subline = line%4;
+                my_weight = subline % 2 ? chart_point[subline/2] : 1-chart_point[subline/2];
+                my_weight *= line/4 ? chart_point[2] : (1-chart_point[2]);
+              }
+            else
+              {
+                Point<2> xy(chart_point[0], chart_point[1]);
+                my_weight = GeometryInfo<2>::d_linear_shape_function(xy, line-8);
+              }
+            my_weight -= weights_lines[line];
+
+            if (std::abs(my_weight) < 1e-13)
+              continue;
+
+            if (cell.line(line)->manifold_id() == my_manifold_id ||
+                cell.line(line)->manifold_id() == numbers::invalid_manifold_id)
+              {
+                weights_vertices[GeometryInfo<3>::line_to_cell_vertices(line,0)]
+                -= my_weight * (1.-line_point);
+                weights_vertices[GeometryInfo<3>::line_to_cell_vertices(line,1)]
+                -= my_weight * (line_point);
+              }
+            else
+              {
+                points.resize(2);
+                weights.resize(2);
+                points[0] = cell.vertex(GeometryInfo<3>::line_to_cell_vertices(line,0));
+                points[1] = cell.vertex(GeometryInfo<3>::line_to_cell_vertices(line,1));
+                weights[0] = 1. - line_point;
+                weights[1] = line_point;
+                new_point -= my_weight *
+                             cell.line(line)->get_manifold().get_new_point(points, weights);
+              }
+          }
+
+        // finally add the contribution of the
+        for (unsigned int v=0; v<GeometryInfo<dim>::vertices_per_cell; ++v)
+          new_point += weights_vertices[v] * cell.vertex(v);
+      }
+    return new_point;
+  }
+}
+
+
+
+template <int dim, int spacedim>
+Point<spacedim>
+TransfiniteInterpolationManifold<dim,spacedim>
+::push_forward(const typename Triangulation<dim,spacedim>::cell_iterator &cell,
+               const Point<dim> &chart_point) const
+{
+  AssertDimension(cell->level(), level_coarse);
+
+  // check that the point is in the unit cell which is the current chart
+  // Tolerance 1e-6 chosen that the method also works with
+  // SphericalManifold
+  Assert(GeometryInfo<dim>::is_inside_unit_cell(chart_point, 1e-6),
+         ExcMessage("chart_point is not in unit interval"));
+
+  return compute_transfinite_interpolation(*cell, chart_point,
+                                           coarse_cell_is_flat[cell->index()]);
+}
+
+
+
+template <int dim, int spacedim>
+DerivativeForm<1,dim,spacedim>
+TransfiniteInterpolationManifold<dim,spacedim>
+::push_forward_gradient(const typename Triangulation<dim,spacedim>::cell_iterator &cell,
+                        const Point<dim> &chart_point) const
+{
+  // compute the derivative with the help of finite differences
+  Point<spacedim> point = compute_transfinite_interpolation(*cell, chart_point,
+                                                            coarse_cell_is_flat[cell->index()]);
+  DerivativeForm<1,dim,spacedim> grad;
+  for (unsigned int d=0; d<dim; ++d)
+    {
+      Point<dim> modified = chart_point;
+      const double step = chart_point[d] > 0.5 ? -1e-8 : 1e-8;
+
+      // avoid checking outside of the unit interval
+      modified[d] += step;
+      Tensor<1,spacedim> difference =
+        compute_transfinite_interpolation(*cell, modified,
+                                          coarse_cell_is_flat[cell->index()]) - point;
+      for (unsigned int e=0; e<spacedim; ++e)
+        grad[e][d] = difference[e] / step;
+    }
+  return grad;
+}
+
+
+
+template <int dim, int spacedim>
+Point<dim>
+TransfiniteInterpolationManifold<dim,spacedim>
+::pull_back(const typename Triangulation<dim,spacedim>::cell_iterator &cell,
+            const Point<spacedim> &point) const
+{
+  Point<dim> outside;
+  for (unsigned int d=0; d<dim; ++d)
+    outside[d] = 20;
+
+  // initial guess: MappingQ1
+  Point<dim> chart_point;
+  try
+    {
+      chart_point = GeometryInfo<dim>::project_to_unit_cell(StaticMappingQ1<dim,spacedim>::mapping.transform_real_to_unit_cell(cell, point));
+    }
+  catch (typename Mapping<dim,spacedim>::ExcTransformationFailed)
+    {
+      for (unsigned int d=0; d<dim; ++d)
+        chart_point[d] = 0.5;
+    }
+
+  // run Newton iteration. As opposed to the various mapping implementations,
+  // this class does not throw exception upon failure as those are relatively
+  // expensive and failure occurs quite regularly in the implementation of the
+  // compute_chart_points method.
+  Tensor<1,spacedim> residual = point - compute_transfinite_interpolation(*cell, chart_point,
+                                coarse_cell_is_flat[cell->index()]);
+  const double tolerance = 1e-21 * cell->diameter() * cell->diameter();
+  double residual_norm_square = residual.norm_square();
+  for (unsigned int i=0; i<100; ++i)
+    {
+      if (residual_norm_square < tolerance)
+        return chart_point;
+
+      // if the determinant is zero, the mapping is not invertible and we are
+      // outside the valid chart region
+      DerivativeForm<1,dim,spacedim> grad = push_forward_gradient(cell, chart_point);
+      if (grad.determinant() <= 0.0)
+        return outside;
+      DerivativeForm<1,dim,spacedim> inv_grad = grad.covariant_form();
+      Tensor<1,dim> update;
+      for (unsigned int d=0; d<spacedim; ++d)
+        for (unsigned int e=0; e<dim; ++e)
+          update[e] += inv_grad[d][e] * residual[d];
+
+      // Line search, accept step if the residual has decreased
+      double alpha = 1.;
+      while (alpha > 1e-7)
+        {
+          Point<dim> guess = chart_point + alpha*update;
+          residual = point - compute_transfinite_interpolation(*cell, guess,
+                                                               coarse_cell_is_flat[cell->index()]);
+          const double residual_norm_new = residual.norm_square();
+          if (residual_norm_new < residual_norm_square)
+            {
+              residual_norm_square = residual_norm_new;
+              chart_point += alpha*update;
+              break;
+            }
+          else
+            alpha *= 0.5;
+        }
+      if (alpha < 1e-7)
+        return outside;
+    }
+  return outside;
+}
+
+
+
+template <int dim, int spacedim>
+std::array<unsigned int, 10>
+TransfiniteInterpolationManifold<dim,spacedim>
+::get_possible_cells_around_points(const std::vector<Point<spacedim> > &points) const
+{
+  // The methods to identify cells around points in GridTools are all written
+  // for the active cells, but we are here looking at some cells at the coarse
+  // level.
+  Assert(triangulation != nullptr, ExcNotInitialized());
+  Assert(triangulation->begin_active()->level() >= level_coarse,
+         ExcMessage("The manifold was initialized with level " +
+                    Utilities::to_string(level_coarse) + " but there are now" +
+                    "active cells on a lower level. Coarsening the mesh is " +
+                    "currently not supported"));
+
+  // This computes the distance of the surrouding points transformed to the unit
+  // cell from the unit cell.
+  typename Triangulation<dim,spacedim>::cell_iterator
+  cell = triangulation->begin(level_coarse),
+  endc = triangulation->end(level_coarse),
+  best_fit = cell;
+  const double circle_gap = 1.3;
+  std::vector<std::pair<double, unsigned int> > distances_and_cells;
+  for ( ; cell != endc; ++cell)
+    {
+      // only consider cells where the current manifold is attached
+      if (&cell->get_manifold() != this)
+        continue;
+
+      // cheap check: if any of the points is not inside a circle around the
+      // center of the loop, we can skip the expensive part below (this assumes
+      // that the manifold does not deform the grid too much)
+      Point<spacedim> center;
+      for (unsigned int v=0; v<GeometryInfo<dim>::vertices_per_cell; ++v)
+        center += cell->vertex(v);
+      center *= 1./GeometryInfo<dim>::vertices_per_cell;
+      double radius_square = 0.;
+      for (unsigned int v=0; v<GeometryInfo<dim>::vertices_per_cell; ++v)
+        radius_square = std::max(radius_square, (center-cell->vertex(v)).norm_square());
+      bool inside_circle = true;
+      for (unsigned int i=0; i<points.size(); ++i)
+        if ((center-points[i]).norm_square() > radius_square*circle_gap)
+          {
+            inside_circle = false;
+            break;
+          }
+      if (inside_circle == false)
+        continue;
+
+      // expensive search
+      try
+        {
+          double current_distance = 0;
+          for (unsigned int i=0; i<points.size(); ++i)
+            {
+              Point<dim> point = StaticMappingQ1<dim,spacedim>::mapping.transform_real_to_unit_cell(cell, points[i]);
+              current_distance += GeometryInfo<dim>::distance_to_unit_cell(point);
+            }
+          distances_and_cells.emplace_back(current_distance, cell->index());
+        }
+      catch (typename Mapping<dim,spacedim>::ExcTransformationFailed &)
+        {
+          // in case the transformation by the mapping fails, use the radius
+          distances_and_cells.emplace_back(std::sqrt(radius_square), cell->index());
+        }
+    }
+  // no coarse cell could be found -> transformation failed
+  AssertThrow(distances_and_cells.size() > 0,
+              (typename Mapping<dim,spacedim>::ExcTransformationFailed()));
+  std::sort(distances_and_cells.begin(), distances_and_cells.end());
+  std::array<unsigned int,10> cells;
+  cells.fill(numbers::invalid_unsigned_int);
+  for (unsigned int i=0; i<distances_and_cells.size() && i<10; ++i)
+    cells[i] = distances_and_cells[i].second;
+
+  return cells;
+}
+
+
+
+template <int dim, int spacedim>
+std::pair<typename Triangulation<dim,spacedim>::cell_iterator,
+    std::vector<Point<dim> > >
+    TransfiniteInterpolationManifold<dim, spacedim>
+    ::compute_chart_points (const std::vector<Point<spacedim> > &surrounding_points) const
+{
+  std::pair<typename Triangulation<dim,spacedim>::cell_iterator,
+      std::vector<Point<dim> > > chart_points;
+  chart_points.second.resize(surrounding_points.size());
+
+  std::array<unsigned int,10> nearby_cells =
+    get_possible_cells_around_points(surrounding_points);
+
+  // check whether all points are inside the unit cell of the current chart
+  for (unsigned int c=0; c<nearby_cells.size(); ++c)
+    {
+      AssertThrow(nearby_cells[c] != numbers::invalid_unsigned_int,
+                  (typename Mapping<dim,spacedim>::ExcTransformationFailed()));
+      typename Triangulation<dim,spacedim>::cell_iterator cell(triangulation,
+                                                               level_coarse,
+                                                               nearby_cells[c]);
+      bool inside_unit_cell = true;
+      for (unsigned int i=0; i<surrounding_points.size(); ++i)
+        {
+          chart_points.second[i] = pull_back(cell, surrounding_points[i]);
+
+          // Tolerance 1e-6 chosen that the method also works with
+          // SphericalManifold
+          if (GeometryInfo<dim>::is_inside_unit_cell(chart_points.second[i],
+                                                     1e-6) == false)
+            {
+              inside_unit_cell = false;
+              break;
+            }
+        }
+      if (inside_unit_cell == true)
+        {
+          chart_points.first = cell;
+          return chart_points;
+        }
+    }
+
+  // a valid inversion should have returned a point above.
+  AssertThrow(false,
+              (typename Mapping<dim,spacedim>::ExcTransformationFailed()));
+  chart_points.second.clear();
+  return chart_points;
+}
+
+
+
+template <int dim, int spacedim>
+Point<spacedim>
+TransfiniteInterpolationManifold<dim, spacedim>
+::get_new_point (const std::vector<Point<spacedim> > &surrounding_points,
+                 const std::vector<double>           &weights) const
+{
+  const std::pair<typename Triangulation<dim,spacedim>::cell_iterator,
+        std::vector<Point<dim> > > chart_points =
+          compute_chart_points(surrounding_points);
+
+  const Point<dim> p_chart = chart_manifold.get_new_point(chart_points.second,weights);
+
+  return push_forward(chart_points.first, p_chart);
+}
+
+
+
+template <int dim, int spacedim>
+void
+TransfiniteInterpolationManifold<dim,spacedim>::
+add_new_points (const std::vector<Point<spacedim> > &surrounding_points,
+                const Table<2,double>               &weights,
+                std::vector<Point<spacedim> >       &new_points) const
+{
+  Assert(weights.size(0) > 0, ExcEmptyObject());
+  AssertDimension(surrounding_points.size(), weights.size(1));
+
+  const std::pair<typename Triangulation<dim,spacedim>::cell_iterator,
+        std::vector<Point<dim> > > chart_points =
+          compute_chart_points(surrounding_points);
+
+  std::vector<Point<dim> > new_points_on_chart;
+  new_points_on_chart.reserve(weights.size(0));
+  chart_manifold.add_new_points(chart_points.second, weights, new_points_on_chart);
+
+  for (unsigned int row=0; row<weights.size(0); ++row)
+    new_points.push_back(push_forward(chart_points.first, new_points_on_chart[row]));
+}
+
+
+
 // explicit instantiations
 #include "manifold_lib.inst"
 
index 63714f9ff46ecce38cbb6312eadf2f77b6d99329..911f3bf7d97052f32128ca967927ae4af54e3413 100644 (file)
@@ -20,6 +20,9 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension :  SPACE_DIMENSIONS
     template class PolarManifold<deal_II_dimension, deal_II_space_dimension>;
     template class SphericalManifold<deal_II_dimension, deal_II_space_dimension>;
     template class CylindricalManifold<deal_II_dimension, deal_II_space_dimension>;
+#if deal_II_dimension > 1
+    template class TransfiniteInterpolationManifold<deal_II_dimension, deal_II_space_dimension>;
+#endif
 #endif
 #if deal_II_dimension == deal_II_space_dimension
     template class TorusManifold<deal_II_dimension>;
diff --git a/tests/manifold/transfinite_manifold_01.cc b/tests/manifold/transfinite_manifold_01.cc
new file mode 100644 (file)
index 0000000..c5b946c
--- /dev/null
@@ -0,0 +1,152 @@
+//-------------------------------------------------------------------
+//    Copyright (C) 2017 by the deal.II authors.
+//
+//    This file is subject to LGPL and may not be  distributed
+//    without copyright and license information. Please refer
+//    to the file deal.II/doc/license.html for the  text  and
+//    further information on this license.
+//
+//-------------------------------------------------------------------
+
+
+// Test that the transfinite manifold creates good points for a circle, a
+// cylinder, and a sphere
+
+#include "../tests.h"
+#include <fstream>
+#include <deal.II/base/logstream.h>
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/tria_accessor.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/tria_boundary_lib.h>
+#include <deal.II/grid/grid_out.h>
+#include <deal.II/grid/manifold_lib.h>
+
+
+template <int dim, int spacedim>
+void do_test(const Triangulation<dim,spacedim> &tria)
+{
+  for (typename Triangulation<dim,spacedim>::cell_iterator cell = tria.begin();
+       cell != tria.end(); ++cell)
+    {
+      deallog << "Lines on cell with center: " << cell->center() << std::endl;
+      for (unsigned int line=0; line<GeometryInfo<dim>::lines_per_cell; ++line)
+        deallog << cell->line(line)->center(/*respect_manifold=*/true) << std::endl;
+      deallog << "Faces on cell with center: " << cell->center() << std::endl;
+      for (unsigned int face=0; face<GeometryInfo<dim>::faces_per_cell; ++face)
+        deallog << cell->face(face)->center(/*respect_manifold=*/true) << std::endl;
+      deallog << "Center with manifold: " << cell->center(true) << std::endl;
+      for (unsigned int face=0; face<GeometryInfo<dim>::faces_per_cell; ++face)
+        if (cell->at_boundary(face))
+          {
+            std::vector<Point<spacedim> > points;
+            points.push_back(cell->face(face)->vertex(0));
+            points.push_back(cell->face(face)->vertex(1));
+            std::vector<double> weights(2);
+            weights[0] = 0.1;
+            weights[1] = 0.9;
+            Point<spacedim> p = cell->get_manifold().get_new_point(points, weights);
+            Point<spacedim> pref = cell->face(face)->get_manifold().get_new_point(points, weights);
+            deallog << "Distance between cell manifold and face manifold: "
+                    << (pref-p) << std::endl;
+            weights[0] = 0.55;
+            weights[1] = 0.45;
+            p = cell->get_manifold().get_new_point(points, weights);
+            pref = cell->face(face)->get_manifold().get_new_point(points, weights);
+            deallog << "Distance between cell manifold and face manifold: "
+                    << (pref-p) << std::endl;
+          }
+    }
+  deallog << std::endl;
+}
+
+template <int dim, int spacedim>
+void test_polar()
+{
+  deallog << "Testing with PolarManifold dim=" << dim
+          << ", spacedim="<< spacedim << std::endl;
+
+  PolarManifold<dim,spacedim> polar_manifold;
+  TransfiniteInterpolationManifold<dim,spacedim> manifold;
+
+  Triangulation<dim,spacedim> tria;
+  GridGenerator::hyper_ball (tria);
+
+  // set all entities to the transfinite manifold except for the boundary where
+  // we put the polar manifold
+  tria.set_all_manifold_ids(1);
+  tria.set_all_manifold_ids_on_boundary(0);
+  tria.set_manifold (0, polar_manifold);
+  manifold.initialize(tria);
+  tria.set_manifold (1, manifold);
+
+  do_test(tria);
+}
+
+template <int dim, int spacedim>
+void test_spherical()
+{
+  deallog << "Testing with SphericalManifold dim=" << dim
+          << ", spacedim="<< spacedim << std::endl;
+
+  SphericalManifold<dim,spacedim> polar_manifold;
+  TransfiniteInterpolationManifold<dim,spacedim> manifold;
+
+  Triangulation<dim,spacedim> tria;
+  GridGenerator::hyper_ball (tria);
+
+  // set all entities to the transfinite manifold except for the boundary where
+  // we put the polar manifold
+  tria.set_all_manifold_ids(1);
+  tria.set_all_manifold_ids_on_boundary(0);
+  tria.set_manifold (0, polar_manifold);
+  manifold.initialize(tria);
+  tria.set_manifold (1, manifold);
+
+  do_test(tria);
+}
+
+void test_cylinder(unsigned int ref=1)
+{
+  const unsigned int dim=3, spacedim = 3;
+  deallog << "Testing with CylindricalManifold in 3d" << std::endl;
+
+  CylindricalManifold<dim,spacedim> cylinder_manifold;
+  TransfiniteInterpolationManifold<dim,spacedim> manifold;
+  Triangulation<dim,spacedim> tria;
+  GridGenerator::cylinder (tria);
+  tria.set_all_manifold_ids(1);
+
+  for (typename Triangulation<dim,spacedim>::active_cell_iterator
+       cell=tria.begin_active(); cell!=tria.end(); ++cell)
+    {
+      for (unsigned int face=0; face<GeometryInfo<dim>::faces_per_cell; ++face)
+        if (cell->at_boundary(face))
+          {
+            bool cell_at_surfaces = true;
+            for (unsigned int i=1; i<GeometryInfo<dim>::vertices_per_face; ++i)
+              if (std::abs(cell->face(face)->vertex(i)[0] - cell->face(face)->vertex(0)[0]) > 1e-10)
+                cell_at_surfaces = false;
+            if (cell_at_surfaces == false)
+              cell->face(face)->set_all_manifold_ids(0);
+          }
+    }
+  tria.set_manifold (0, cylinder_manifold);
+  manifold.initialize(tria);
+  tria.set_manifold (1, manifold);
+
+  do_test(tria);
+}
+
+int main ()
+{
+  initlog();
+
+  test_polar<2,2>();
+  test_spherical<2,2>();
+  test_spherical<3,3>();
+  test_cylinder();
+
+  return 0;
+}
+
diff --git a/tests/manifold/transfinite_manifold_01.output b/tests/manifold/transfinite_manifold_01.output
new file mode 100644 (file)
index 0000000..602efe0
--- /dev/null
@@ -0,0 +1,540 @@
+
+DEAL::Testing with PolarManifold dim=2, spacedim=2
+DEAL::Lines on cell with center: 0.00000 -0.500000
+DEAL::-0.500000 -0.500000
+DEAL::0.500000 -0.500000
+DEAL::-1.83697e-16 -1.00000
+DEAL::0.00000 -0.292893
+DEAL::Faces on cell with center: 0.00000 -0.500000
+DEAL::-0.500000 -0.500000
+DEAL::0.500000 -0.500000
+DEAL::-1.83697e-16 -1.00000
+DEAL::0.00000 -0.292893
+DEAL::Center with manifold: -9.71445e-17 -0.646447
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000
+DEAL::Lines on cell with center: -0.500000 -1.38778e-17
+DEAL::-1.00000 1.22465e-16
+DEAL::-0.292893 0.00000
+DEAL::-0.500000 -0.500000
+DEAL::-0.500000 0.500000
+DEAL::Faces on cell with center: -0.500000 -1.38778e-17
+DEAL::-1.00000 1.22465e-16
+DEAL::-0.292893 0.00000
+DEAL::-0.500000 -0.500000
+DEAL::-0.500000 0.500000
+DEAL::Center with manifold: -0.646447 5.55112e-17
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000
+DEAL::Lines on cell with center: 0.00000 0.00000
+DEAL::-0.292893 0.00000
+DEAL::0.292893 0.00000
+DEAL::0.00000 -0.292893
+DEAL::0.00000 0.292893
+DEAL::Faces on cell with center: 0.00000 0.00000
+DEAL::-0.292893 0.00000
+DEAL::0.292893 0.00000
+DEAL::0.00000 -0.292893
+DEAL::0.00000 0.292893
+DEAL::Center with manifold: 0.00000 0.00000
+DEAL::Lines on cell with center: 0.500000 0.00000
+DEAL::0.500000 -0.500000
+DEAL::0.500000 0.500000
+DEAL::1.00000 0.00000
+DEAL::0.292893 0.00000
+DEAL::Faces on cell with center: 0.500000 0.00000
+DEAL::0.500000 -0.500000
+DEAL::0.500000 0.500000
+DEAL::1.00000 0.00000
+DEAL::0.292893 0.00000
+DEAL::Center with manifold: 0.646447 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000
+DEAL::Lines on cell with center: -1.38778e-17 0.500000
+DEAL::6.12323e-17 1.00000
+DEAL::0.00000 0.292893
+DEAL::-0.500000 0.500000
+DEAL::0.500000 0.500000
+DEAL::Faces on cell with center: -1.38778e-17 0.500000
+DEAL::6.12323e-17 1.00000
+DEAL::0.00000 0.292893
+DEAL::-0.500000 0.500000
+DEAL::0.500000 0.500000
+DEAL::Center with manifold: 2.77556e-17 0.646447
+DEAL::Distance between cell manifold and face manifold: 2.22045e-16 -1.11022e-16
+DEAL::Distance between cell manifold and face manifold: 0.00000 1.11022e-16
+DEAL::
+DEAL::Testing with SphericalManifold dim=2, spacedim=2
+DEAL::Lines on cell with center: 0.00000 -0.500000
+DEAL::-0.500000 -0.500000
+DEAL::0.500000 -0.500000
+DEAL::0.00000 -1.00000
+DEAL::0.00000 -0.292893
+DEAL::Faces on cell with center: 0.00000 -0.500000
+DEAL::-0.500000 -0.500000
+DEAL::0.500000 -0.500000
+DEAL::0.00000 -1.00000
+DEAL::0.00000 -0.292893
+DEAL::Center with manifold: 0.00000 -0.646447
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000
+DEAL::Lines on cell with center: -0.500000 -1.38778e-17
+DEAL::-1.00000 0.00000
+DEAL::-0.292893 0.00000
+DEAL::-0.500000 -0.500000
+DEAL::-0.500000 0.500000
+DEAL::Faces on cell with center: -0.500000 -1.38778e-17
+DEAL::-1.00000 0.00000
+DEAL::-0.292893 0.00000
+DEAL::-0.500000 -0.500000
+DEAL::-0.500000 0.500000
+DEAL::Center with manifold: -0.646447 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000
+DEAL::Lines on cell with center: 0.00000 0.00000
+DEAL::-0.292893 0.00000
+DEAL::0.292893 0.00000
+DEAL::0.00000 -0.292893
+DEAL::0.00000 0.292893
+DEAL::Faces on cell with center: 0.00000 0.00000
+DEAL::-0.292893 0.00000
+DEAL::0.292893 0.00000
+DEAL::0.00000 -0.292893
+DEAL::0.00000 0.292893
+DEAL::Center with manifold: 0.00000 0.00000
+DEAL::Lines on cell with center: 0.500000 0.00000
+DEAL::0.500000 -0.500000
+DEAL::0.500000 0.500000
+DEAL::1.00000 0.00000
+DEAL::0.292893 0.00000
+DEAL::Faces on cell with center: 0.500000 0.00000
+DEAL::0.500000 -0.500000
+DEAL::0.500000 0.500000
+DEAL::1.00000 0.00000
+DEAL::0.292893 0.00000
+DEAL::Center with manifold: 0.646447 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000
+DEAL::Lines on cell with center: -1.38778e-17 0.500000
+DEAL::0.00000 1.00000
+DEAL::0.00000 0.292893
+DEAL::-0.500000 0.500000
+DEAL::0.500000 0.500000
+DEAL::Faces on cell with center: -1.38778e-17 0.500000
+DEAL::0.00000 1.00000
+DEAL::0.00000 0.292893
+DEAL::-0.500000 0.500000
+DEAL::0.500000 0.500000
+DEAL::Center with manifold: 0.00000 0.646447
+DEAL::Distance between cell manifold and face manifold: 1.11022e-16 -2.22045e-16
+DEAL::Distance between cell manifold and face manifold: 6.93889e-17 1.11022e-16
+DEAL::
+DEAL::Testing with SphericalManifold dim=3, spacedim=3
+DEAL::Lines on cell with center: 0.00000 0.00000 -6.93889e-18
+DEAL::-0.211325 0.00000 -0.211325
+DEAL::0.211325 -1.38778e-17 -0.211325
+DEAL::0.00000 -0.211325 -0.211325
+DEAL::0.00000 0.211325 -0.211325
+DEAL::-0.211325 -1.38778e-17 0.211325
+DEAL::0.211325 -1.38778e-17 0.211325
+DEAL::0.00000 -0.211325 0.211325
+DEAL::0.00000 0.211325 0.211325
+DEAL::-0.211325 -0.211325 0.00000
+DEAL::0.211325 -0.211325 -1.38778e-17
+DEAL::-0.211325 0.211325 -1.38778e-17
+DEAL::0.211325 0.211325 -1.38778e-17
+DEAL::Faces on cell with center: 0.00000 0.00000 -6.93889e-18
+DEAL::-0.211325 -2.08167e-17 -4.85723e-17
+DEAL::0.211325 0.00000 -4.16334e-17
+DEAL::0.00000 -0.211325 0.00000
+DEAL::-3.46945e-17 0.211325 -1.38778e-17
+DEAL::-1.38778e-17 -6.93889e-18 -0.211325
+DEAL::-1.38778e-17 -2.77556e-17 0.211325
+DEAL::Center with manifold: -9.36751e-17 -4.85723e-17 -7.63278e-17
+DEAL::Lines on cell with center: 0.00000 0.00000 -0.394338
+DEAL::-0.707107 0.00000 -0.707107
+DEAL::0.707107 0.00000 -0.707107
+DEAL::0.00000 -0.707107 -0.707107
+DEAL::0.00000 0.707107 -0.707107
+DEAL::-0.211325 0.00000 -0.211325
+DEAL::0.211325 -1.38778e-17 -0.211325
+DEAL::0.00000 -0.211325 -0.211325
+DEAL::0.00000 0.211325 -0.211325
+DEAL::-0.394338 -0.394338 -0.394338
+DEAL::0.394338 -0.394338 -0.394338
+DEAL::-0.394338 0.394338 -0.394338
+DEAL::0.394338 0.394338 -0.394338
+DEAL::Faces on cell with center: 0.00000 0.00000 -0.394338
+DEAL::-0.459216 6.93889e-18 -0.459216
+DEAL::0.459216 0.00000 -0.459216
+DEAL::6.93889e-18 -0.459216 -0.459216
+DEAL::0.00000 0.459216 -0.459216
+DEAL::0.00000 0.00000 -1.00000
+DEAL::-1.38778e-17 -6.93889e-18 -0.211325
+DEAL::Center with manifold: -1.17961e-16 3.46945e-18 -0.605662
+DEAL::Distance between cell manifold and face manifold: 1.11022e-16 -1.11022e-16 -1.11022e-16
+DEAL::Distance between cell manifold and face manifold: -1.38778e-17 -1.11022e-16 -1.11022e-16
+DEAL::Lines on cell with center: 0.394338 0.00000 2.08167e-17
+DEAL::0.394338 -0.394338 -0.394338
+DEAL::0.394338 0.394338 -0.394338
+DEAL::0.707107 0.00000 -0.707107
+DEAL::0.211325 -1.38778e-17 -0.211325
+DEAL::0.394338 -0.394338 0.394338
+DEAL::0.394338 0.394338 0.394338
+DEAL::0.707107 0.00000 0.707107
+DEAL::0.211325 -1.38778e-17 0.211325
+DEAL::0.707107 -0.707107 0.00000
+DEAL::0.707107 0.707107 0.00000
+DEAL::0.211325 -0.211325 -1.38778e-17
+DEAL::0.211325 0.211325 -1.38778e-17
+DEAL::Faces on cell with center: 0.394338 0.00000 2.08167e-17
+DEAL::0.459216 -0.459216 0.00000
+DEAL::0.459216 0.459216 0.00000
+DEAL::1.00000 0.00000 0.00000
+DEAL::0.211325 0.00000 -4.16334e-17
+DEAL::0.459216 0.00000 -0.459216
+DEAL::0.459216 0.00000 0.459216
+DEAL::Center with manifold: 0.605662 2.42861e-17 -6.24500e-17
+DEAL::Distance between cell manifold and face manifold: 1.11022e-16 -1.11022e-16 5.55112e-17
+DEAL::Distance between cell manifold and face manifold: 1.11022e-16 -1.11022e-16 -1.38778e-17
+DEAL::Lines on cell with center: 0.00000 2.08167e-17 0.394338
+DEAL::-0.394338 -0.394338 0.394338
+DEAL::0.394338 -0.394338 0.394338
+DEAL::0.00000 -0.707107 0.707107
+DEAL::0.00000 -0.211325 0.211325
+DEAL::-0.394338 0.394338 0.394338
+DEAL::0.394338 0.394338 0.394338
+DEAL::0.00000 0.707107 0.707107
+DEAL::0.00000 0.211325 0.211325
+DEAL::-0.707107 0.00000 0.707107
+DEAL::0.707107 0.00000 0.707107
+DEAL::-0.211325 -1.38778e-17 0.211325
+DEAL::0.211325 -1.38778e-17 0.211325
+DEAL::Faces on cell with center: 0.00000 2.08167e-17 0.394338
+DEAL::-0.459216 0.00000 0.459216
+DEAL::0.459216 0.00000 0.459216
+DEAL::0.00000 0.00000 1.00000
+DEAL::-1.38778e-17 -2.77556e-17 0.211325
+DEAL::6.93889e-18 -0.459216 0.459216
+DEAL::0.00000 0.459216 0.459216
+DEAL::Center with manifold: 2.42861e-17 -6.24500e-17 0.605662
+DEAL::Distance between cell manifold and face manifold: -1.11022e-16 5.55112e-17 1.11022e-16
+DEAL::Distance between cell manifold and face manifold: -1.11022e-16 -1.94289e-16 1.11022e-16
+DEAL::Lines on cell with center: -0.394338 3.46945e-18 1.73472e-17
+DEAL::-0.707107 0.00000 -0.707107
+DEAL::-0.211325 0.00000 -0.211325
+DEAL::-0.394338 -0.394338 -0.394338
+DEAL::-0.394338 0.394338 -0.394338
+DEAL::-0.707107 0.00000 0.707107
+DEAL::-0.211325 -1.38778e-17 0.211325
+DEAL::-0.394338 -0.394338 0.394338
+DEAL::-0.394338 0.394338 0.394338
+DEAL::-0.707107 -0.707107 0.00000
+DEAL::-0.211325 -0.211325 0.00000
+DEAL::-0.707107 0.707107 0.00000
+DEAL::-0.211325 0.211325 -1.38778e-17
+DEAL::Faces on cell with center: -0.394338 3.46945e-18 1.73472e-17
+DEAL::-1.00000 0.00000 0.00000
+DEAL::-0.211325 -2.08167e-17 -4.85723e-17
+DEAL::-0.459216 -0.459216 6.93889e-18
+DEAL::-0.459216 0.459216 0.00000
+DEAL::-0.459216 6.93889e-18 -0.459216
+DEAL::-0.459216 0.00000 0.459216
+DEAL::Center with manifold: -0.605662 -3.46945e-17 -1.04083e-17
+DEAL::Distance between cell manifold and face manifold: -1.11022e-16 1.11022e-16 -1.11022e-16
+DEAL::Distance between cell manifold and face manifold: -1.11022e-16 -1.38778e-17 -1.11022e-16
+DEAL::Lines on cell with center: 0.00000 -0.394338 2.08167e-17
+DEAL::-0.394338 -0.394338 -0.394338
+DEAL::0.394338 -0.394338 -0.394338
+DEAL::0.00000 -0.707107 -0.707107
+DEAL::0.00000 -0.211325 -0.211325
+DEAL::-0.394338 -0.394338 0.394338
+DEAL::0.394338 -0.394338 0.394338
+DEAL::0.00000 -0.707107 0.707107
+DEAL::0.00000 -0.211325 0.211325
+DEAL::-0.707107 -0.707107 0.00000
+DEAL::0.707107 -0.707107 0.00000
+DEAL::-0.211325 -0.211325 0.00000
+DEAL::0.211325 -0.211325 -1.38778e-17
+DEAL::Faces on cell with center: 0.00000 -0.394338 2.08167e-17
+DEAL::-0.459216 -0.459216 6.93889e-18
+DEAL::0.459216 -0.459216 0.00000
+DEAL::0.00000 -1.00000 0.00000
+DEAL::0.00000 -0.211325 0.00000
+DEAL::6.93889e-18 -0.459216 -0.459216
+DEAL::6.93889e-18 -0.459216 0.459216
+DEAL::Center with manifold: 2.42861e-17 -0.605662 -6.24500e-17
+DEAL::Distance between cell manifold and face manifold: 1.11022e-16 -2.22045e-16 -5.55112e-17
+DEAL::Distance between cell manifold and face manifold: -2.22045e-16 -2.22045e-16 -9.71445e-17
+DEAL::Lines on cell with center: 3.46945e-18 0.394338 1.73472e-17
+DEAL::0.00000 0.707107 -0.707107
+DEAL::0.00000 0.211325 -0.211325
+DEAL::-0.394338 0.394338 -0.394338
+DEAL::0.394338 0.394338 -0.394338
+DEAL::0.00000 0.707107 0.707107
+DEAL::0.00000 0.211325 0.211325
+DEAL::-0.394338 0.394338 0.394338
+DEAL::0.394338 0.394338 0.394338
+DEAL::-0.707107 0.707107 0.00000
+DEAL::-0.211325 0.211325 -1.38778e-17
+DEAL::0.707107 0.707107 0.00000
+DEAL::0.211325 0.211325 -1.38778e-17
+DEAL::Faces on cell with center: 3.46945e-18 0.394338 1.73472e-17
+DEAL::0.00000 1.00000 0.00000
+DEAL::-3.46945e-17 0.211325 -1.38778e-17
+DEAL::-0.459216 0.459216 0.00000
+DEAL::0.459216 0.459216 0.00000
+DEAL::0.00000 0.459216 -0.459216
+DEAL::0.00000 0.459216 0.459216
+DEAL::Center with manifold: -3.46945e-17 0.605662 -1.04083e-17
+DEAL::Distance between cell manifold and face manifold: 1.11022e-16 1.11022e-16 -1.11022e-16
+DEAL::Distance between cell manifold and face manifold: -1.38778e-17 1.11022e-16 -1.11022e-16
+DEAL::
+DEAL::Testing with CylindricalManifold in 3d
+DEAL::Lines on cell with center: -0.500000 0.00000 -0.500000
+DEAL::-0.500000 0.707107 -0.707107
+DEAL::-0.500000 -0.707107 -0.707107
+DEAL::-1.00000 0.00000 -1.00000
+DEAL::0.00000 0.00000 -1.00000
+DEAL::-0.500000 0.292893 -0.292893
+DEAL::-0.500000 -0.292893 -0.292893
+DEAL::-1.00000 0.00000 -0.292893
+DEAL::0.00000 0.00000 -0.292893
+DEAL::-1.00000 0.500000 -0.500000
+DEAL::-1.00000 -0.500000 -0.500000
+DEAL::0.00000 0.500000 -0.500000
+DEAL::0.00000 -0.500000 -0.500000
+DEAL::Faces on cell with center: -0.500000 0.00000 -0.500000
+DEAL::-0.500000 0.500000 -0.500000
+DEAL::-0.500000 -0.500000 -0.500000
+DEAL::-1.00000 0.00000 -0.646447
+DEAL::0.00000 0.00000 -0.646447
+DEAL::-0.500000 0.00000 -1.00000
+DEAL::-0.500000 0.00000 -0.292893
+DEAL::Center with manifold: -0.500000 0.00000 -0.646447
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Lines on cell with center: -0.500000 0.500000 -6.93889e-18
+DEAL::-0.500000 0.707107 -0.707107
+DEAL::-0.500000 0.292893 -0.292893
+DEAL::-1.00000 0.500000 -0.500000
+DEAL::0.00000 0.500000 -0.500000
+DEAL::-0.500000 0.707107 0.707107
+DEAL::-0.500000 0.292893 0.292893
+DEAL::-1.00000 0.500000 0.500000
+DEAL::0.00000 0.500000 0.500000
+DEAL::-1.00000 1.00000 0.00000
+DEAL::-1.00000 0.292893 0.00000
+DEAL::0.00000 1.00000 0.00000
+DEAL::0.00000 0.292893 0.00000
+DEAL::Faces on cell with center: -0.500000 0.500000 -6.93889e-18
+DEAL::-0.500000 1.00000 0.00000
+DEAL::-0.500000 0.292893 0.00000
+DEAL::-1.00000 0.646447 0.00000
+DEAL::0.00000 0.646447 0.00000
+DEAL::-0.500000 0.500000 -0.500000
+DEAL::-0.500000 0.500000 0.500000
+DEAL::Center with manifold: -0.500000 0.646447 1.66533e-16
+DEAL::Distance between cell manifold and face manifold: -2.77556e-17 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Lines on cell with center: -0.500000 0.00000 1.38778e-17
+DEAL::-0.500000 0.292893 -0.292893
+DEAL::-0.500000 -0.292893 -0.292893
+DEAL::-1.00000 0.00000 -0.292893
+DEAL::0.00000 0.00000 -0.292893
+DEAL::-0.500000 0.292893 0.292893
+DEAL::-0.500000 -0.292893 0.292893
+DEAL::-1.00000 0.00000 0.292893
+DEAL::0.00000 0.00000 0.292893
+DEAL::-1.00000 0.292893 0.00000
+DEAL::-1.00000 -0.292893 0.00000
+DEAL::0.00000 0.292893 0.00000
+DEAL::0.00000 -0.292893 0.00000
+DEAL::Faces on cell with center: -0.500000 0.00000 1.38778e-17
+DEAL::-0.500000 0.292893 0.00000
+DEAL::-0.500000 -0.292893 0.00000
+DEAL::-1.00000 0.00000 0.00000
+DEAL::0.00000 0.00000 0.00000
+DEAL::-0.500000 0.00000 -0.292893
+DEAL::-0.500000 0.00000 0.292893
+DEAL::Center with manifold: -0.500000 0.00000 6.93889e-17
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Lines on cell with center: -0.500000 -0.500000 0.00000
+DEAL::-0.500000 -0.707107 -0.707107
+DEAL::-0.500000 -0.707107 0.707107
+DEAL::-1.00000 -1.00000 0.00000
+DEAL::0.00000 -1.00000 0.00000
+DEAL::-0.500000 -0.292893 -0.292893
+DEAL::-0.500000 -0.292893 0.292893
+DEAL::-1.00000 -0.292893 0.00000
+DEAL::0.00000 -0.292893 0.00000
+DEAL::-1.00000 -0.500000 -0.500000
+DEAL::-1.00000 -0.500000 0.500000
+DEAL::0.00000 -0.500000 -0.500000
+DEAL::0.00000 -0.500000 0.500000
+DEAL::Faces on cell with center: -0.500000 -0.500000 0.00000
+DEAL::-0.500000 -0.500000 -0.500000
+DEAL::-0.500000 -0.500000 0.500000
+DEAL::-1.00000 -0.646447 0.00000
+DEAL::0.00000 -0.646447 0.00000
+DEAL::-0.500000 -1.00000 0.00000
+DEAL::-0.500000 -0.292893 0.00000
+DEAL::Center with manifold: -0.500000 -0.646447 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Lines on cell with center: -0.500000 6.93889e-18 0.500000
+DEAL::-0.500000 0.707107 0.707107
+DEAL::-0.500000 0.292893 0.292893
+DEAL::-1.00000 0.500000 0.500000
+DEAL::0.00000 0.500000 0.500000
+DEAL::-0.500000 -0.707107 0.707107
+DEAL::-0.500000 -0.292893 0.292893
+DEAL::-1.00000 -0.500000 0.500000
+DEAL::0.00000 -0.500000 0.500000
+DEAL::-1.00000 0.00000 1.00000
+DEAL::-1.00000 0.00000 0.292893
+DEAL::0.00000 0.00000 1.00000
+DEAL::0.00000 0.00000 0.292893
+DEAL::Faces on cell with center: -0.500000 6.93889e-18 0.500000
+DEAL::-0.500000 0.00000 1.00000
+DEAL::-0.500000 0.00000 0.292893
+DEAL::-1.00000 0.00000 0.646447
+DEAL::0.00000 0.00000 0.646447
+DEAL::-0.500000 0.500000 0.500000
+DEAL::-0.500000 -0.500000 0.500000
+DEAL::Center with manifold: -0.500000 -1.66533e-16 0.646447
+DEAL::Distance between cell manifold and face manifold: -2.77556e-17 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Lines on cell with center: 0.500000 0.00000 -0.500000
+DEAL::0.500000 0.707107 -0.707107
+DEAL::0.500000 -0.707107 -0.707107
+DEAL::0.00000 0.00000 -1.00000
+DEAL::1.00000 0.00000 -1.00000
+DEAL::0.500000 0.292893 -0.292893
+DEAL::0.500000 -0.292893 -0.292893
+DEAL::0.00000 0.00000 -0.292893
+DEAL::1.00000 0.00000 -0.292893
+DEAL::0.00000 0.500000 -0.500000
+DEAL::0.00000 -0.500000 -0.500000
+DEAL::1.00000 0.500000 -0.500000
+DEAL::1.00000 -0.500000 -0.500000
+DEAL::Faces on cell with center: 0.500000 0.00000 -0.500000
+DEAL::0.500000 0.500000 -0.500000
+DEAL::0.500000 -0.500000 -0.500000
+DEAL::0.00000 0.00000 -0.646447
+DEAL::1.00000 0.00000 -0.646447
+DEAL::0.500000 0.00000 -1.00000
+DEAL::0.500000 0.00000 -0.292893
+DEAL::Center with manifold: 0.500000 0.00000 -0.646447
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Lines on cell with center: 0.500000 0.500000 -6.93889e-18
+DEAL::0.500000 0.707107 -0.707107
+DEAL::0.500000 0.292893 -0.292893
+DEAL::0.00000 0.500000 -0.500000
+DEAL::1.00000 0.500000 -0.500000
+DEAL::0.500000 0.707107 0.707107
+DEAL::0.500000 0.292893 0.292893
+DEAL::0.00000 0.500000 0.500000
+DEAL::1.00000 0.500000 0.500000
+DEAL::0.00000 1.00000 0.00000
+DEAL::0.00000 0.292893 0.00000
+DEAL::1.00000 1.00000 0.00000
+DEAL::1.00000 0.292893 0.00000
+DEAL::Faces on cell with center: 0.500000 0.500000 -6.93889e-18
+DEAL::0.500000 1.00000 0.00000
+DEAL::0.500000 0.292893 0.00000
+DEAL::0.00000 0.646447 0.00000
+DEAL::1.00000 0.646447 0.00000
+DEAL::0.500000 0.500000 -0.500000
+DEAL::0.500000 0.500000 0.500000
+DEAL::Center with manifold: 0.500000 0.646447 1.66533e-16
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Lines on cell with center: 0.500000 0.00000 1.38778e-17
+DEAL::0.500000 0.292893 -0.292893
+DEAL::0.500000 -0.292893 -0.292893
+DEAL::0.00000 0.00000 -0.292893
+DEAL::1.00000 0.00000 -0.292893
+DEAL::0.500000 0.292893 0.292893
+DEAL::0.500000 -0.292893 0.292893
+DEAL::0.00000 0.00000 0.292893
+DEAL::1.00000 0.00000 0.292893
+DEAL::0.00000 0.292893 0.00000
+DEAL::0.00000 -0.292893 0.00000
+DEAL::1.00000 0.292893 0.00000
+DEAL::1.00000 -0.292893 0.00000
+DEAL::Faces on cell with center: 0.500000 0.00000 1.38778e-17
+DEAL::0.500000 0.292893 0.00000
+DEAL::0.500000 -0.292893 0.00000
+DEAL::0.00000 0.00000 0.00000
+DEAL::1.00000 0.00000 0.00000
+DEAL::0.500000 0.00000 -0.292893
+DEAL::0.500000 0.00000 0.292893
+DEAL::Center with manifold: 0.500000 0.00000 6.93889e-17
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Lines on cell with center: 0.500000 -0.500000 0.00000
+DEAL::0.500000 -0.707107 -0.707107
+DEAL::0.500000 -0.707107 0.707107
+DEAL::0.00000 -1.00000 0.00000
+DEAL::1.00000 -1.00000 0.00000
+DEAL::0.500000 -0.292893 -0.292893
+DEAL::0.500000 -0.292893 0.292893
+DEAL::0.00000 -0.292893 0.00000
+DEAL::1.00000 -0.292893 0.00000
+DEAL::0.00000 -0.500000 -0.500000
+DEAL::0.00000 -0.500000 0.500000
+DEAL::1.00000 -0.500000 -0.500000
+DEAL::1.00000 -0.500000 0.500000
+DEAL::Faces on cell with center: 0.500000 -0.500000 0.00000
+DEAL::0.500000 -0.500000 -0.500000
+DEAL::0.500000 -0.500000 0.500000
+DEAL::0.00000 -0.646447 0.00000
+DEAL::1.00000 -0.646447 0.00000
+DEAL::0.500000 -1.00000 0.00000
+DEAL::0.500000 -0.292893 0.00000
+DEAL::Center with manifold: 0.500000 -0.646447 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Lines on cell with center: 0.500000 6.93889e-18 0.500000
+DEAL::0.500000 0.707107 0.707107
+DEAL::0.500000 0.292893 0.292893
+DEAL::0.00000 0.500000 0.500000
+DEAL::1.00000 0.500000 0.500000
+DEAL::0.500000 -0.707107 0.707107
+DEAL::0.500000 -0.292893 0.292893
+DEAL::0.00000 -0.500000 0.500000
+DEAL::1.00000 -0.500000 0.500000
+DEAL::0.00000 0.00000 1.00000
+DEAL::0.00000 0.00000 0.292893
+DEAL::1.00000 0.00000 1.00000
+DEAL::1.00000 0.00000 0.292893
+DEAL::Faces on cell with center: 0.500000 6.93889e-18 0.500000
+DEAL::0.500000 0.00000 1.00000
+DEAL::0.500000 0.00000 0.292893
+DEAL::0.00000 0.00000 0.646447
+DEAL::1.00000 0.00000 0.646447
+DEAL::0.500000 0.500000 0.500000
+DEAL::0.500000 -0.500000 0.500000
+DEAL::Center with manifold: 0.500000 -1.66533e-16 0.646447
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::
diff --git a/tests/manifold/transfinite_manifold_02.cc b/tests/manifold/transfinite_manifold_02.cc
new file mode 100644 (file)
index 0000000..2e9f755
--- /dev/null
@@ -0,0 +1,161 @@
+//-------------------------------------------------------------------
+//    Copyright (C) 2017 by the deal.II authors.
+//
+//    This file is subject to LGPL and may not be  distributed
+//    without copyright and license information. Please refer
+//    to the file deal.II/doc/license.html for the  text  and
+//    further information on this license.
+//
+//-------------------------------------------------------------------
+
+
+// Similar to transfinite_manifold_01 but now only setting the transfinite
+// manifold on the cells at the boundary, skipping the interior ones
+
+#include "../tests.h"
+#include <fstream>
+#include <deal.II/base/logstream.h>
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/tria_accessor.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/tria_boundary_lib.h>
+#include <deal.II/grid/grid_out.h>
+#include <deal.II/grid/manifold_lib.h>
+
+
+template <int dim, int spacedim>
+void do_test(const Triangulation<dim,spacedim> &tria)
+{
+  for (typename Triangulation<dim,spacedim>::cell_iterator cell = tria.begin();
+       cell != tria.end(); ++cell)
+    {
+      deallog << "Lines on cell with center: " << cell->center() << std::endl;
+      for (unsigned int line=0; line<GeometryInfo<dim>::lines_per_cell; ++line)
+        deallog << cell->line(line)->center(/*respect_manifold=*/true) << std::endl;
+      deallog << "Faces on cell with center: " << cell->center() << std::endl;
+      for (unsigned int face=0; face<GeometryInfo<dim>::faces_per_cell; ++face)
+        deallog << cell->face(face)->center(/*respect_manifold=*/true) << std::endl;
+      deallog << "Center with manifold: " << cell->center(true) << std::endl;
+      for (unsigned int face=0; face<GeometryInfo<dim>::faces_per_cell; ++face)
+        if (cell->at_boundary(face))
+          {
+            std::vector<Point<spacedim> > points;
+            points.push_back(cell->face(face)->vertex(0));
+            points.push_back(cell->face(face)->vertex(1));
+            std::vector<double> weights(2);
+            weights[0] = 0.1;
+            weights[1] = 0.9;
+            Point<spacedim> p = cell->get_manifold().get_new_point(points, weights);
+            Point<spacedim> pref = cell->face(face)->get_manifold().get_new_point(points, weights);
+            deallog << "Distance between cell manifold and face manifold: "
+                    << (pref-p) << std::endl;
+            weights[0] = 0.55;
+            weights[1] = 0.45;
+            p = cell->get_manifold().get_new_point(points, weights);
+            pref = cell->face(face)->get_manifold().get_new_point(points, weights);
+            deallog << "Distance between cell manifold and face manifold: "
+                    << (pref-p) << std::endl;
+          }
+    }
+  deallog << std::endl;
+}
+
+template <int dim, int spacedim>
+void test_polar()
+{
+  deallog << "Testing with PolarManifold dim=" << dim
+          << ", spacedim="<< spacedim << std::endl;
+
+  PolarManifold<dim,spacedim> polar_manifold;
+  TransfiniteInterpolationManifold<dim,spacedim> manifold;
+
+  Triangulation<dim,spacedim> tria;
+  GridGenerator::hyper_ball (tria);
+
+  // set all cells at the boundary to the transfinite manifold except for the
+  // boundaries where we put the polar manifold
+  for (typename Triangulation<dim>::active_cell_iterator cell = tria.begin();
+       cell != tria.end(); ++cell)
+    if (cell->at_boundary())
+      cell->set_manifold_id(1);
+
+  tria.set_all_manifold_ids_on_boundary(0);
+  tria.set_manifold (0, polar_manifold);
+  manifold.initialize(tria);
+  tria.set_manifold (1, manifold);
+
+  do_test(tria);
+}
+
+template <int dim, int spacedim>
+void test_spherical()
+{
+  deallog << "Testing with SphericalManifold dim=" << dim
+          << ", spacedim="<< spacedim << std::endl;
+
+  SphericalManifold<dim,spacedim> polar_manifold;
+  TransfiniteInterpolationManifold<dim,spacedim> manifold;
+
+  Triangulation<dim,spacedim> tria;
+  GridGenerator::hyper_ball (tria);
+
+  // set all cells at the boundary to the transfinite manifold except for the
+  // boundaries where we put the polar manifold
+  for (typename Triangulation<dim>::active_cell_iterator cell = tria.begin();
+       cell != tria.end(); ++cell)
+    if (cell->at_boundary())
+      cell->set_manifold_id(1);
+
+  tria.set_all_manifold_ids_on_boundary(0);
+  tria.set_manifold (0, polar_manifold);
+  manifold.initialize(tria);
+  tria.set_manifold (1, manifold);
+
+  do_test(tria);
+}
+
+void test_cylinder(unsigned int ref=1)
+{
+  const unsigned int dim=3, spacedim = 3;
+  deallog << "Testing with CylindricalManifold in 3d" << std::endl;
+
+  CylindricalManifold<dim,spacedim> cylinder_manifold;
+  TransfiniteInterpolationManifold<dim,spacedim> manifold;
+  Triangulation<dim,spacedim> tria;
+  GridGenerator::cylinder (tria);
+
+  for (typename Triangulation<dim,spacedim>::active_cell_iterator
+       cell=tria.begin_active(); cell!=tria.end(); ++cell)
+    {
+      if (cell->at_boundary())
+        cell->set_manifold_id(1);
+      for (unsigned int face=0; face<GeometryInfo<dim>::faces_per_cell; ++face)
+        if (cell->at_boundary(face))
+          {
+            bool cell_at_surfaces = true;
+            for (unsigned int i=1; i<GeometryInfo<dim>::vertices_per_face; ++i)
+              if (std::abs(cell->face(face)->vertex(i)[0] - cell->face(face)->vertex(0)[0]) > 1e-10)
+                cell_at_surfaces = false;
+            if (cell_at_surfaces == false)
+              cell->face(face)->set_all_manifold_ids(0);
+          }
+    }
+  tria.set_manifold (0, cylinder_manifold);
+  manifold.initialize(tria);
+  tria.set_manifold (1, manifold);
+
+  do_test(tria);
+}
+
+int main ()
+{
+  initlog();
+
+  test_polar<2,2>();
+  test_spherical<2,2>();
+  test_spherical<3,3>();
+  test_cylinder();
+
+  return 0;
+}
+
diff --git a/tests/manifold/transfinite_manifold_02.output b/tests/manifold/transfinite_manifold_02.output
new file mode 100644 (file)
index 0000000..75b400a
--- /dev/null
@@ -0,0 +1,540 @@
+
+DEAL::Testing with PolarManifold dim=2, spacedim=2
+DEAL::Lines on cell with center: 0.00000 -0.500000
+DEAL::-0.500000 -0.500000
+DEAL::0.500000 -0.500000
+DEAL::-1.83697e-16 -1.00000
+DEAL::0.00000 -0.292893
+DEAL::Faces on cell with center: 0.00000 -0.500000
+DEAL::-0.500000 -0.500000
+DEAL::0.500000 -0.500000
+DEAL::-1.83697e-16 -1.00000
+DEAL::0.00000 -0.292893
+DEAL::Center with manifold: -9.71445e-17 -0.646447
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000
+DEAL::Lines on cell with center: -0.500000 -1.38778e-17
+DEAL::-1.00000 1.22465e-16
+DEAL::-0.292893 0.00000
+DEAL::-0.500000 -0.500000
+DEAL::-0.500000 0.500000
+DEAL::Faces on cell with center: -0.500000 -1.38778e-17
+DEAL::-1.00000 1.22465e-16
+DEAL::-0.292893 0.00000
+DEAL::-0.500000 -0.500000
+DEAL::-0.500000 0.500000
+DEAL::Center with manifold: -0.646447 5.55112e-17
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000
+DEAL::Lines on cell with center: 0.00000 0.00000
+DEAL::-0.292893 0.00000
+DEAL::0.292893 0.00000
+DEAL::0.00000 -0.292893
+DEAL::0.00000 0.292893
+DEAL::Faces on cell with center: 0.00000 0.00000
+DEAL::-0.292893 0.00000
+DEAL::0.292893 0.00000
+DEAL::0.00000 -0.292893
+DEAL::0.00000 0.292893
+DEAL::Center with manifold: 0.00000 0.00000
+DEAL::Lines on cell with center: 0.500000 0.00000
+DEAL::0.500000 -0.500000
+DEAL::0.500000 0.500000
+DEAL::1.00000 0.00000
+DEAL::0.292893 0.00000
+DEAL::Faces on cell with center: 0.500000 0.00000
+DEAL::0.500000 -0.500000
+DEAL::0.500000 0.500000
+DEAL::1.00000 0.00000
+DEAL::0.292893 0.00000
+DEAL::Center with manifold: 0.646447 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000
+DEAL::Lines on cell with center: -1.38778e-17 0.500000
+DEAL::6.12323e-17 1.00000
+DEAL::0.00000 0.292893
+DEAL::-0.500000 0.500000
+DEAL::0.500000 0.500000
+DEAL::Faces on cell with center: -1.38778e-17 0.500000
+DEAL::6.12323e-17 1.00000
+DEAL::0.00000 0.292893
+DEAL::-0.500000 0.500000
+DEAL::0.500000 0.500000
+DEAL::Center with manifold: 2.77556e-17 0.646447
+DEAL::Distance between cell manifold and face manifold: 2.22045e-16 -1.11022e-16
+DEAL::Distance between cell manifold and face manifold: 0.00000 1.11022e-16
+DEAL::
+DEAL::Testing with SphericalManifold dim=2, spacedim=2
+DEAL::Lines on cell with center: 0.00000 -0.500000
+DEAL::-0.500000 -0.500000
+DEAL::0.500000 -0.500000
+DEAL::0.00000 -1.00000
+DEAL::0.00000 -0.292893
+DEAL::Faces on cell with center: 0.00000 -0.500000
+DEAL::-0.500000 -0.500000
+DEAL::0.500000 -0.500000
+DEAL::0.00000 -1.00000
+DEAL::0.00000 -0.292893
+DEAL::Center with manifold: 0.00000 -0.646447
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000
+DEAL::Lines on cell with center: -0.500000 -1.38778e-17
+DEAL::-1.00000 0.00000
+DEAL::-0.292893 0.00000
+DEAL::-0.500000 -0.500000
+DEAL::-0.500000 0.500000
+DEAL::Faces on cell with center: -0.500000 -1.38778e-17
+DEAL::-1.00000 0.00000
+DEAL::-0.292893 0.00000
+DEAL::-0.500000 -0.500000
+DEAL::-0.500000 0.500000
+DEAL::Center with manifold: -0.646447 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000
+DEAL::Lines on cell with center: 0.00000 0.00000
+DEAL::-0.292893 0.00000
+DEAL::0.292893 0.00000
+DEAL::0.00000 -0.292893
+DEAL::0.00000 0.292893
+DEAL::Faces on cell with center: 0.00000 0.00000
+DEAL::-0.292893 0.00000
+DEAL::0.292893 0.00000
+DEAL::0.00000 -0.292893
+DEAL::0.00000 0.292893
+DEAL::Center with manifold: 0.00000 0.00000
+DEAL::Lines on cell with center: 0.500000 0.00000
+DEAL::0.500000 -0.500000
+DEAL::0.500000 0.500000
+DEAL::1.00000 0.00000
+DEAL::0.292893 0.00000
+DEAL::Faces on cell with center: 0.500000 0.00000
+DEAL::0.500000 -0.500000
+DEAL::0.500000 0.500000
+DEAL::1.00000 0.00000
+DEAL::0.292893 0.00000
+DEAL::Center with manifold: 0.646447 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000
+DEAL::Lines on cell with center: -1.38778e-17 0.500000
+DEAL::0.00000 1.00000
+DEAL::0.00000 0.292893
+DEAL::-0.500000 0.500000
+DEAL::0.500000 0.500000
+DEAL::Faces on cell with center: -1.38778e-17 0.500000
+DEAL::0.00000 1.00000
+DEAL::0.00000 0.292893
+DEAL::-0.500000 0.500000
+DEAL::0.500000 0.500000
+DEAL::Center with manifold: 0.00000 0.646447
+DEAL::Distance between cell manifold and face manifold: 1.11022e-16 -2.22045e-16
+DEAL::Distance between cell manifold and face manifold: 6.93889e-17 1.11022e-16
+DEAL::
+DEAL::Testing with SphericalManifold dim=3, spacedim=3
+DEAL::Lines on cell with center: 0.00000 0.00000 -6.93889e-18
+DEAL::-0.211325 0.00000 -0.211325
+DEAL::0.211325 0.00000 -0.211325
+DEAL::0.00000 -0.211325 -0.211325
+DEAL::0.00000 0.211325 -0.211325
+DEAL::-0.211325 0.00000 0.211325
+DEAL::0.211325 0.00000 0.211325
+DEAL::0.00000 -0.211325 0.211325
+DEAL::0.00000 0.211325 0.211325
+DEAL::-0.211325 -0.211325 0.00000
+DEAL::0.211325 -0.211325 0.00000
+DEAL::-0.211325 0.211325 0.00000
+DEAL::0.211325 0.211325 0.00000
+DEAL::Faces on cell with center: 0.00000 0.00000 -6.93889e-18
+DEAL::-0.211325 0.00000 0.00000
+DEAL::0.211325 0.00000 0.00000
+DEAL::0.00000 -0.211325 0.00000
+DEAL::0.00000 0.211325 0.00000
+DEAL::0.00000 0.00000 -0.211325
+DEAL::0.00000 0.00000 0.211325
+DEAL::Center with manifold: 0.00000 0.00000 0.00000
+DEAL::Lines on cell with center: 0.00000 0.00000 -0.394338
+DEAL::-0.707107 0.00000 -0.707107
+DEAL::0.707107 0.00000 -0.707107
+DEAL::0.00000 -0.707107 -0.707107
+DEAL::0.00000 0.707107 -0.707107
+DEAL::-0.211325 0.00000 -0.211325
+DEAL::0.211325 0.00000 -0.211325
+DEAL::0.00000 -0.211325 -0.211325
+DEAL::0.00000 0.211325 -0.211325
+DEAL::-0.394338 -0.394338 -0.394338
+DEAL::0.394338 -0.394338 -0.394338
+DEAL::-0.394338 0.394338 -0.394338
+DEAL::0.394338 0.394338 -0.394338
+DEAL::Faces on cell with center: 0.00000 0.00000 -0.394338
+DEAL::-0.394338 0.00000 -0.394338
+DEAL::0.394338 0.00000 -0.394338
+DEAL::6.93889e-18 -0.394338 -0.394338
+DEAL::6.93889e-18 0.394338 -0.394338
+DEAL::0.00000 0.00000 -1.00000
+DEAL::0.00000 0.00000 -0.211325
+DEAL::Center with manifold: -1.11022e-16 -1.11022e-16 -0.571260
+DEAL::Distance between cell manifold and face manifold: 1.11022e-16 -1.11022e-16 -1.11022e-16
+DEAL::Distance between cell manifold and face manifold: -1.38778e-17 -1.11022e-16 -1.11022e-16
+DEAL::Lines on cell with center: 0.394338 0.00000 2.08167e-17
+DEAL::0.394338 -0.394338 -0.394338
+DEAL::0.394338 0.394338 -0.394338
+DEAL::0.707107 0.00000 -0.707107
+DEAL::0.211325 0.00000 -0.211325
+DEAL::0.394338 -0.394338 0.394338
+DEAL::0.394338 0.394338 0.394338
+DEAL::0.707107 0.00000 0.707107
+DEAL::0.211325 0.00000 0.211325
+DEAL::0.707107 -0.707107 0.00000
+DEAL::0.707107 0.707107 0.00000
+DEAL::0.211325 -0.211325 0.00000
+DEAL::0.211325 0.211325 0.00000
+DEAL::Faces on cell with center: 0.394338 0.00000 2.08167e-17
+DEAL::0.394338 -0.394338 6.93889e-18
+DEAL::0.394338 0.394338 6.93889e-18
+DEAL::1.00000 0.00000 0.00000
+DEAL::0.211325 0.00000 0.00000
+DEAL::0.394338 0.00000 -0.394338
+DEAL::0.394338 0.00000 0.394338
+DEAL::Center with manifold: 0.571260 2.08167e-17 -5.55112e-17
+DEAL::Distance between cell manifold and face manifold: 1.11022e-16 -1.11022e-16 5.55112e-17
+DEAL::Distance between cell manifold and face manifold: 1.11022e-16 -1.11022e-16 -1.38778e-17
+DEAL::Lines on cell with center: 0.00000 2.08167e-17 0.394338
+DEAL::-0.394338 -0.394338 0.394338
+DEAL::0.394338 -0.394338 0.394338
+DEAL::0.00000 -0.707107 0.707107
+DEAL::0.00000 -0.211325 0.211325
+DEAL::-0.394338 0.394338 0.394338
+DEAL::0.394338 0.394338 0.394338
+DEAL::0.00000 0.707107 0.707107
+DEAL::0.00000 0.211325 0.211325
+DEAL::-0.707107 0.00000 0.707107
+DEAL::0.707107 0.00000 0.707107
+DEAL::-0.211325 0.00000 0.211325
+DEAL::0.211325 0.00000 0.211325
+DEAL::Faces on cell with center: 0.00000 2.08167e-17 0.394338
+DEAL::-0.394338 6.93889e-18 0.394338
+DEAL::0.394338 0.00000 0.394338
+DEAL::0.00000 0.00000 1.00000
+DEAL::0.00000 0.00000 0.211325
+DEAL::0.00000 -0.394338 0.394338
+DEAL::0.00000 0.394338 0.394338
+DEAL::Center with manifold: 2.08167e-17 -5.55112e-17 0.571260
+DEAL::Distance between cell manifold and face manifold: -1.11022e-16 5.55112e-17 1.11022e-16
+DEAL::Distance between cell manifold and face manifold: -1.11022e-16 -1.94289e-16 1.11022e-16
+DEAL::Lines on cell with center: -0.394338 3.46945e-18 1.73472e-17
+DEAL::-0.707107 0.00000 -0.707107
+DEAL::-0.211325 0.00000 -0.211325
+DEAL::-0.394338 -0.394338 -0.394338
+DEAL::-0.394338 0.394338 -0.394338
+DEAL::-0.707107 0.00000 0.707107
+DEAL::-0.211325 0.00000 0.211325
+DEAL::-0.394338 -0.394338 0.394338
+DEAL::-0.394338 0.394338 0.394338
+DEAL::-0.707107 -0.707107 0.00000
+DEAL::-0.211325 -0.211325 0.00000
+DEAL::-0.707107 0.707107 0.00000
+DEAL::-0.211325 0.211325 0.00000
+DEAL::Faces on cell with center: -0.394338 3.46945e-18 1.73472e-17
+DEAL::-1.00000 0.00000 0.00000
+DEAL::-0.211325 0.00000 0.00000
+DEAL::-0.394338 -0.394338 0.00000
+DEAL::-0.394338 0.394338 0.00000
+DEAL::-0.394338 0.00000 -0.394338
+DEAL::-0.394338 6.93889e-18 0.394338
+DEAL::Center with manifold: -0.571260 -3.81639e-17 -1.73472e-17
+DEAL::Distance between cell manifold and face manifold: -1.11022e-16 1.11022e-16 -1.11022e-16
+DEAL::Distance between cell manifold and face manifold: -1.11022e-16 -1.38778e-17 -1.11022e-16
+DEAL::Lines on cell with center: 0.00000 -0.394338 2.08167e-17
+DEAL::-0.394338 -0.394338 -0.394338
+DEAL::0.394338 -0.394338 -0.394338
+DEAL::0.00000 -0.707107 -0.707107
+DEAL::0.00000 -0.211325 -0.211325
+DEAL::-0.394338 -0.394338 0.394338
+DEAL::0.394338 -0.394338 0.394338
+DEAL::0.00000 -0.707107 0.707107
+DEAL::0.00000 -0.211325 0.211325
+DEAL::-0.707107 -0.707107 0.00000
+DEAL::0.707107 -0.707107 0.00000
+DEAL::-0.211325 -0.211325 0.00000
+DEAL::0.211325 -0.211325 0.00000
+DEAL::Faces on cell with center: 0.00000 -0.394338 2.08167e-17
+DEAL::-0.394338 -0.394338 0.00000
+DEAL::0.394338 -0.394338 6.93889e-18
+DEAL::0.00000 -1.00000 0.00000
+DEAL::0.00000 -0.211325 0.00000
+DEAL::6.93889e-18 -0.394338 -0.394338
+DEAL::0.00000 -0.394338 0.394338
+DEAL::Center with manifold: 2.08167e-17 -0.571260 -5.55112e-17
+DEAL::Distance between cell manifold and face manifold: 1.11022e-16 -2.22045e-16 -5.55112e-17
+DEAL::Distance between cell manifold and face manifold: -2.22045e-16 -2.22045e-16 -9.71445e-17
+DEAL::Lines on cell with center: 3.46945e-18 0.394338 1.73472e-17
+DEAL::0.00000 0.707107 -0.707107
+DEAL::0.00000 0.211325 -0.211325
+DEAL::-0.394338 0.394338 -0.394338
+DEAL::0.394338 0.394338 -0.394338
+DEAL::0.00000 0.707107 0.707107
+DEAL::0.00000 0.211325 0.211325
+DEAL::-0.394338 0.394338 0.394338
+DEAL::0.394338 0.394338 0.394338
+DEAL::-0.707107 0.707107 0.00000
+DEAL::-0.211325 0.211325 0.00000
+DEAL::0.707107 0.707107 0.00000
+DEAL::0.211325 0.211325 0.00000
+DEAL::Faces on cell with center: 3.46945e-18 0.394338 1.73472e-17
+DEAL::0.00000 1.00000 0.00000
+DEAL::0.00000 0.211325 0.00000
+DEAL::-0.394338 0.394338 0.00000
+DEAL::0.394338 0.394338 6.93889e-18
+DEAL::6.93889e-18 0.394338 -0.394338
+DEAL::0.00000 0.394338 0.394338
+DEAL::Center with manifold: -3.81639e-17 0.571260 -1.73472e-17
+DEAL::Distance between cell manifold and face manifold: 1.11022e-16 1.11022e-16 -1.11022e-16
+DEAL::Distance between cell manifold and face manifold: -1.38778e-17 1.11022e-16 -1.11022e-16
+DEAL::
+DEAL::Testing with CylindricalManifold in 3d
+DEAL::Lines on cell with center: -0.500000 0.00000 -0.500000
+DEAL::-0.500000 0.707107 -0.707107
+DEAL::-0.500000 -0.707107 -0.707107
+DEAL::-1.00000 0.00000 -1.00000
+DEAL::0.00000 0.00000 -1.00000
+DEAL::-0.500000 0.292893 -0.292893
+DEAL::-0.500000 -0.292893 -0.292893
+DEAL::-1.00000 0.00000 -0.292893
+DEAL::0.00000 0.00000 -0.292893
+DEAL::-1.00000 0.500000 -0.500000
+DEAL::-1.00000 -0.500000 -0.500000
+DEAL::0.00000 0.500000 -0.500000
+DEAL::0.00000 -0.500000 -0.500000
+DEAL::Faces on cell with center: -0.500000 0.00000 -0.500000
+DEAL::-0.500000 0.500000 -0.500000
+DEAL::-0.500000 -0.500000 -0.500000
+DEAL::-1.00000 0.00000 -0.646447
+DEAL::0.00000 0.00000 -0.500000
+DEAL::-0.500000 0.00000 -1.00000
+DEAL::-0.500000 0.00000 -0.292893
+DEAL::Center with manifold: -0.500000 0.00000 -0.634243
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Lines on cell with center: -0.500000 0.500000 -6.93889e-18
+DEAL::-0.500000 0.707107 -0.707107
+DEAL::-0.500000 0.292893 -0.292893
+DEAL::-1.00000 0.500000 -0.500000
+DEAL::0.00000 0.500000 -0.500000
+DEAL::-0.500000 0.707107 0.707107
+DEAL::-0.500000 0.292893 0.292893
+DEAL::-1.00000 0.500000 0.500000
+DEAL::0.00000 0.500000 0.500000
+DEAL::-1.00000 1.00000 0.00000
+DEAL::-1.00000 0.292893 0.00000
+DEAL::0.00000 1.00000 0.00000
+DEAL::0.00000 0.292893 0.00000
+DEAL::Faces on cell with center: -0.500000 0.500000 -6.93889e-18
+DEAL::-0.500000 1.00000 0.00000
+DEAL::-0.500000 0.292893 0.00000
+DEAL::-1.00000 0.646447 0.00000
+DEAL::0.00000 0.500000 0.00000
+DEAL::-0.500000 0.500000 -0.500000
+DEAL::-0.500000 0.500000 0.500000
+DEAL::Center with manifold: -0.500000 0.634243 1.66533e-16
+DEAL::Distance between cell manifold and face manifold: -2.77556e-17 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Lines on cell with center: -0.500000 0.00000 1.38778e-17
+DEAL::-0.500000 0.292893 -0.292893
+DEAL::-0.500000 -0.292893 -0.292893
+DEAL::-1.00000 0.00000 -0.292893
+DEAL::0.00000 0.00000 -0.292893
+DEAL::-0.500000 0.292893 0.292893
+DEAL::-0.500000 -0.292893 0.292893
+DEAL::-1.00000 0.00000 0.292893
+DEAL::0.00000 0.00000 0.292893
+DEAL::-1.00000 0.292893 0.00000
+DEAL::-1.00000 -0.292893 0.00000
+DEAL::0.00000 0.292893 0.00000
+DEAL::0.00000 -0.292893 0.00000
+DEAL::Faces on cell with center: -0.500000 0.00000 1.38778e-17
+DEAL::-0.500000 0.292893 0.00000
+DEAL::-0.500000 -0.292893 0.00000
+DEAL::-1.00000 0.00000 0.00000
+DEAL::0.00000 0.00000 0.00000
+DEAL::-0.500000 0.00000 -0.292893
+DEAL::-0.500000 0.00000 0.292893
+DEAL::Center with manifold: -0.500000 0.00000 6.93889e-17
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Lines on cell with center: -0.500000 -0.500000 0.00000
+DEAL::-0.500000 -0.707107 -0.707107
+DEAL::-0.500000 -0.707107 0.707107
+DEAL::-1.00000 -1.00000 0.00000
+DEAL::0.00000 -1.00000 0.00000
+DEAL::-0.500000 -0.292893 -0.292893
+DEAL::-0.500000 -0.292893 0.292893
+DEAL::-1.00000 -0.292893 0.00000
+DEAL::0.00000 -0.292893 0.00000
+DEAL::-1.00000 -0.500000 -0.500000
+DEAL::-1.00000 -0.500000 0.500000
+DEAL::0.00000 -0.500000 -0.500000
+DEAL::0.00000 -0.500000 0.500000
+DEAL::Faces on cell with center: -0.500000 -0.500000 0.00000
+DEAL::-0.500000 -0.500000 -0.500000
+DEAL::-0.500000 -0.500000 0.500000
+DEAL::-1.00000 -0.646447 0.00000
+DEAL::0.00000 -0.500000 0.00000
+DEAL::-0.500000 -1.00000 0.00000
+DEAL::-0.500000 -0.292893 0.00000
+DEAL::Center with manifold: -0.500000 -0.634243 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Lines on cell with center: -0.500000 6.93889e-18 0.500000
+DEAL::-0.500000 0.707107 0.707107
+DEAL::-0.500000 0.292893 0.292893
+DEAL::-1.00000 0.500000 0.500000
+DEAL::0.00000 0.500000 0.500000
+DEAL::-0.500000 -0.707107 0.707107
+DEAL::-0.500000 -0.292893 0.292893
+DEAL::-1.00000 -0.500000 0.500000
+DEAL::0.00000 -0.500000 0.500000
+DEAL::-1.00000 0.00000 1.00000
+DEAL::-1.00000 0.00000 0.292893
+DEAL::0.00000 0.00000 1.00000
+DEAL::0.00000 0.00000 0.292893
+DEAL::Faces on cell with center: -0.500000 6.93889e-18 0.500000
+DEAL::-0.500000 0.00000 1.00000
+DEAL::-0.500000 0.00000 0.292893
+DEAL::-1.00000 0.00000 0.646447
+DEAL::0.00000 0.00000 0.500000
+DEAL::-0.500000 0.500000 0.500000
+DEAL::-0.500000 -0.500000 0.500000
+DEAL::Center with manifold: -0.500000 -1.66533e-16 0.634243
+DEAL::Distance between cell manifold and face manifold: -2.77556e-17 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Lines on cell with center: 0.500000 0.00000 -0.500000
+DEAL::0.500000 0.707107 -0.707107
+DEAL::0.500000 -0.707107 -0.707107
+DEAL::0.00000 0.00000 -1.00000
+DEAL::1.00000 0.00000 -1.00000
+DEAL::0.500000 0.292893 -0.292893
+DEAL::0.500000 -0.292893 -0.292893
+DEAL::0.00000 0.00000 -0.292893
+DEAL::1.00000 0.00000 -0.292893
+DEAL::0.00000 0.500000 -0.500000
+DEAL::0.00000 -0.500000 -0.500000
+DEAL::1.00000 0.500000 -0.500000
+DEAL::1.00000 -0.500000 -0.500000
+DEAL::Faces on cell with center: 0.500000 0.00000 -0.500000
+DEAL::0.500000 0.500000 -0.500000
+DEAL::0.500000 -0.500000 -0.500000
+DEAL::0.00000 0.00000 -0.500000
+DEAL::1.00000 0.00000 -0.500000
+DEAL::0.500000 0.00000 -1.00000
+DEAL::0.500000 0.00000 -0.292893
+DEAL::Center with manifold: 0.500000 0.00000 -0.622039
+DEAL::Distance between cell manifold and face manifold: 0.00000 5.55112e-17 -5.55112e-17
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Lines on cell with center: 0.500000 0.500000 -6.93889e-18
+DEAL::0.500000 0.707107 -0.707107
+DEAL::0.500000 0.292893 -0.292893
+DEAL::0.00000 0.500000 -0.500000
+DEAL::1.00000 0.500000 -0.500000
+DEAL::0.500000 0.707107 0.707107
+DEAL::0.500000 0.292893 0.292893
+DEAL::0.00000 0.500000 0.500000
+DEAL::1.00000 0.500000 0.500000
+DEAL::0.00000 1.00000 0.00000
+DEAL::0.00000 0.292893 0.00000
+DEAL::1.00000 1.00000 0.00000
+DEAL::1.00000 0.292893 0.00000
+DEAL::Faces on cell with center: 0.500000 0.500000 -6.93889e-18
+DEAL::0.500000 1.00000 0.00000
+DEAL::0.500000 0.292893 0.00000
+DEAL::0.00000 0.500000 0.00000
+DEAL::1.00000 0.500000 0.00000
+DEAL::0.500000 0.500000 -0.500000
+DEAL::0.500000 0.500000 0.500000
+DEAL::Center with manifold: 0.500000 0.622039 1.52656e-16
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 -0.0737620 -0.0590096
+DEAL::Distance between cell manifold and face manifold: 0.00000 -0.287930 0.0287930
+DEAL::Lines on cell with center: 0.500000 0.00000 1.38778e-17
+DEAL::0.500000 0.292893 -0.292893
+DEAL::0.500000 -0.292893 -0.292893
+DEAL::0.00000 0.00000 -0.292893
+DEAL::1.00000 0.00000 -0.292893
+DEAL::0.500000 0.292893 0.292893
+DEAL::0.500000 -0.292893 0.292893
+DEAL::0.00000 0.00000 0.292893
+DEAL::1.00000 0.00000 0.292893
+DEAL::0.00000 0.292893 0.00000
+DEAL::0.00000 -0.292893 0.00000
+DEAL::1.00000 0.292893 0.00000
+DEAL::1.00000 -0.292893 0.00000
+DEAL::Faces on cell with center: 0.500000 0.00000 1.38778e-17
+DEAL::0.500000 0.292893 0.00000
+DEAL::0.500000 -0.292893 0.00000
+DEAL::0.00000 0.00000 0.00000
+DEAL::1.00000 0.00000 0.00000
+DEAL::0.500000 0.00000 -0.292893
+DEAL::0.500000 0.00000 0.292893
+DEAL::Center with manifold: 0.500000 0.00000 6.93889e-17
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Lines on cell with center: 0.500000 -0.500000 0.00000
+DEAL::0.500000 -0.707107 -0.707107
+DEAL::0.500000 -0.707107 0.707107
+DEAL::0.00000 -1.00000 0.00000
+DEAL::1.00000 -1.00000 0.00000
+DEAL::0.500000 -0.292893 -0.292893
+DEAL::0.500000 -0.292893 0.292893
+DEAL::0.00000 -0.292893 0.00000
+DEAL::1.00000 -0.292893 0.00000
+DEAL::0.00000 -0.500000 -0.500000
+DEAL::0.00000 -0.500000 0.500000
+DEAL::1.00000 -0.500000 -0.500000
+DEAL::1.00000 -0.500000 0.500000
+DEAL::Faces on cell with center: 0.500000 -0.500000 0.00000
+DEAL::0.500000 -0.500000 -0.500000
+DEAL::0.500000 -0.500000 0.500000
+DEAL::0.00000 -0.500000 0.00000
+DEAL::1.00000 -0.500000 0.00000
+DEAL::0.500000 -1.00000 0.00000
+DEAL::0.500000 -0.292893 0.00000
+DEAL::Center with manifold: 0.500000 -0.622039 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 -5.55112e-17 -5.55112e-17
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Lines on cell with center: 0.500000 6.93889e-18 0.500000
+DEAL::0.500000 0.707107 0.707107
+DEAL::0.500000 0.292893 0.292893
+DEAL::0.00000 0.500000 0.500000
+DEAL::1.00000 0.500000 0.500000
+DEAL::0.500000 -0.707107 0.707107
+DEAL::0.500000 -0.292893 0.292893
+DEAL::0.00000 -0.500000 0.500000
+DEAL::1.00000 -0.500000 0.500000
+DEAL::0.00000 0.00000 1.00000
+DEAL::0.00000 0.00000 0.292893
+DEAL::1.00000 0.00000 1.00000
+DEAL::1.00000 0.00000 0.292893
+DEAL::Faces on cell with center: 0.500000 6.93889e-18 0.500000
+DEAL::0.500000 0.00000 1.00000
+DEAL::0.500000 0.00000 0.292893
+DEAL::0.00000 0.00000 0.500000
+DEAL::1.00000 0.00000 0.500000
+DEAL::0.500000 0.500000 0.500000
+DEAL::0.500000 -0.500000 0.500000
+DEAL::Center with manifold: 0.500000 -1.52656e-16 0.622039
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.0590096 -0.0737620
+DEAL::Distance between cell manifold and face manifold: 0.00000 -0.0287930 -0.287930
+DEAL::
diff --git a/tests/manifold/transfinite_manifold_03.cc b/tests/manifold/transfinite_manifold_03.cc
new file mode 100644 (file)
index 0000000..a101c16
--- /dev/null
@@ -0,0 +1,143 @@
+//-------------------------------------------------------------------
+//    Copyright (C) 2017 by the deal.II authors.
+//
+//    This file is subject to LGPL and may not be  distributed
+//    without copyright and license information. Please refer
+//    to the file deal.II/doc/license.html for the  text  and
+//    further information on this license.
+//
+//-------------------------------------------------------------------
+
+
+// Similar to transfinite_manifold_01 but now applying a MappingQGeneric and
+// computing some areas
+
+#include "../tests.h"
+#include <fstream>
+#include <deal.II/base/logstream.h>
+#include <deal.II/base/quadrature_lib.h>
+#include <deal.II/fe/mapping_q_generic.h>
+#include <deal.II/fe/fe_nothing.h>
+#include <deal.II/fe/fe_values.h>
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/tria_accessor.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/tria_boundary_lib.h>
+#include <deal.II/grid/grid_out.h>
+#include <deal.II/grid/manifold_lib.h>
+
+
+template <int dim, int spacedim>
+void do_test(const Triangulation<dim,spacedim> &tria)
+{
+  for (unsigned int degree = 1; degree < 5; ++degree)
+    {
+      MappingQGeneric<dim,spacedim> mapping(degree);
+      FE_Nothing<dim,spacedim> fe;
+      QGauss<dim> gauss(degree+1);
+      FEValues<dim,spacedim> fe_values(mapping, fe, gauss, update_JxW_values);
+      double volume = 0;
+      for (typename Triangulation<dim,spacedim>::cell_iterator cell = tria.begin_active();
+           cell != tria.end(); ++cell)
+        {
+          fe_values.reinit(cell);
+          for (unsigned int q=0; q<gauss.size(); ++q)
+            volume += fe_values.JxW(q);
+        }
+
+      deallog << std::setw(4) << tria.n_active_cells() << " cells, degree=" << degree << ": "
+              << std::setprecision(12) << std::setw(14) << volume << std::endl;
+    }
+  deallog << std::endl;
+}
+
+template <int dim, int spacedim>
+void test_polar()
+{
+  deallog << "Testing with PolarManifold dim=" << dim
+          << ", spacedim="<< spacedim << std::endl;
+
+  PolarManifold<dim,spacedim> polar_manifold;
+  TransfiniteInterpolationManifold<dim,spacedim> manifold;
+
+  Triangulation<dim,spacedim> tria;
+  GridGenerator::hyper_ball (tria);
+
+  tria.set_all_manifold_ids(1);
+  tria.set_all_manifold_ids_on_boundary(0);
+  tria.set_manifold (0, polar_manifold);
+  manifold.initialize(tria);
+  tria.set_manifold (1, manifold);
+
+  tria.refine_global(2);
+
+  do_test(tria);
+}
+
+template <int dim, int spacedim>
+void test_spherical()
+{
+  deallog << "Testing with SphericalManifold dim=" << dim
+          << ", spacedim="<< spacedim << std::endl;
+
+  SphericalManifold<dim,spacedim> spherical_manifold;
+  TransfiniteInterpolationManifold<dim,spacedim> manifold;
+
+  Triangulation<dim,spacedim> tria;
+  GridGenerator::hyper_ball (tria);
+
+  tria.set_all_manifold_ids(1);
+  tria.set_all_manifold_ids_on_boundary(0);
+  tria.set_manifold (0, spherical_manifold);
+  manifold.initialize(tria);
+  tria.set_manifold (1, manifold);
+  tria.refine_global(4-dim);
+
+  do_test(tria);
+}
+
+void test_cylinder(unsigned int ref=1)
+{
+  const unsigned int dim=3, spacedim = 3;
+  deallog << "Testing with CylindricalManifold in 3d" << std::endl;
+
+  CylindricalManifold<dim,spacedim> cylinder_manifold;
+  TransfiniteInterpolationManifold<dim,spacedim> manifold;
+  Triangulation<dim,spacedim> tria;
+  GridGenerator::cylinder (tria);
+
+  tria.set_all_manifold_ids(1);
+  for (typename Triangulation<dim,spacedim>::active_cell_iterator
+       cell=tria.begin_active(); cell!=tria.end(); ++cell)
+    {
+      for (unsigned int face=0; face<GeometryInfo<dim>::faces_per_cell; ++face)
+        if (cell->at_boundary(face))
+          {
+            bool cell_at_surfaces = true;
+            for (unsigned int i=1; i<GeometryInfo<dim>::vertices_per_face; ++i)
+              if (std::abs(cell->face(face)->vertex(i)[0] - cell->face(face)->vertex(0)[0]) > 1e-10)
+                cell_at_surfaces = false;
+            if (cell_at_surfaces == false)
+              cell->face(face)->set_all_manifold_ids(0);
+          }
+    }
+  tria.set_manifold (0, cylinder_manifold);
+  manifold.initialize(tria);
+  tria.set_manifold (1, manifold);
+  tria.refine_global(1);
+
+  do_test(tria);
+}
+
+int main ()
+{
+  initlog();
+
+  test_polar<2,2>();
+  test_spherical<2,2>();
+  test_spherical<3,3>();
+  test_cylinder();
+
+  return 0;
+}
+
diff --git a/tests/manifold/transfinite_manifold_03.output b/tests/manifold/transfinite_manifold_03.output
new file mode 100644 (file)
index 0000000..afd60ee
--- /dev/null
@@ -0,0 +1,25 @@
+
+DEAL::Testing with PolarManifold dim=2, spacedim=2
+DEAL::80   cells, degree=1: 3.06146745892 
+DEAL::80   cells, degree=2: 3.14143771670 
+DEAL::80   cells, degree=3: 3.14159250174 
+DEAL::80   cells, degree=4: 3.14159265350 
+DEAL::
+DEAL::Testing with SphericalManifold dim=2, spacedim=2
+DEAL::80   cells, degree=1: 3.06146745892 
+DEAL::80   cells, degree=2: 3.14143771670 
+DEAL::80   cells, degree=3: 3.14158706712 
+DEAL::80   cells, degree=4: 3.14159263998 
+DEAL::
+DEAL::Testing with SphericalManifold dim=3, spacedim=3
+DEAL::56   cells, degree=1: 3.21034343105 
+DEAL::56   cells, degree=2: 4.18080434784 
+DEAL::56   cells, degree=3: 4.18759484733 
+DEAL::56   cells, degree=4: 4.18877458327 
+DEAL::
+DEAL::Testing with CylindricalManifold in 3d
+DEAL::80   cells, degree=1: 5.65685424949 
+DEAL::80   cells, degree=2: 6.27829514062 
+DEAL::80   cells, degree=3: 6.28283288022 
+DEAL::80   cells, degree=4: 6.28318111226 
+DEAL::
diff --git a/tests/manifold/transfinite_manifold_04.cc b/tests/manifold/transfinite_manifold_04.cc
new file mode 100644 (file)
index 0000000..c58ce66
--- /dev/null
@@ -0,0 +1,124 @@
+//-------------------------------------------------------------------
+//    Copyright (C) 2017 by the deal.II authors.
+//
+//    This file is subject to LGPL and may not be  distributed
+//    without copyright and license information. Please refer
+//    to the file deal.II/doc/license.html for the  text  and
+//    further information on this license.
+//
+//-------------------------------------------------------------------
+
+
+// Test that clearing a triangulation is behaving correctly for the transfinite
+// interpolation on a test that is similar to transfinite_manifold_01.cc
+
+#include "../tests.h"
+#include <fstream>
+#include <deal.II/base/logstream.h>
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/tria_accessor.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/tria_boundary_lib.h>
+#include <deal.II/grid/grid_out.h>
+#include <deal.II/grid/manifold_lib.h>
+
+
+template <int dim, int spacedim>
+void do_test(const Triangulation<dim,spacedim> &tria)
+{
+  for (typename Triangulation<dim,spacedim>::cell_iterator cell = tria.begin();
+       cell != tria.end(); ++cell)
+    {
+      deallog << "Lines on cell with center: " << cell->center() << std::endl;
+      for (unsigned int line=0; line<GeometryInfo<dim>::lines_per_cell; ++line)
+        deallog << cell->line(line)->center(/*respect_manifold=*/true) << std::endl;
+      deallog << "Faces on cell with center: " << cell->center() << std::endl;
+      for (unsigned int face=0; face<GeometryInfo<dim>::faces_per_cell; ++face)
+        deallog << cell->face(face)->center(/*respect_manifold=*/true) << std::endl;
+      deallog << "Center with manifold: " << cell->center(true) << std::endl;
+      for (unsigned int face=0; face<GeometryInfo<dim>::faces_per_cell; ++face)
+        if (cell->at_boundary(face))
+          {
+            std::vector<Point<spacedim> > points;
+            points.push_back(cell->face(face)->vertex(0));
+            points.push_back(cell->face(face)->vertex(1));
+            std::vector<double> weights(2);
+            weights[0] = 0.1;
+            weights[1] = 0.9;
+            Point<spacedim> p = cell->get_manifold().get_new_point(points, weights);
+            Point<spacedim> pref = cell->face(face)->get_manifold().get_new_point(points, weights);
+            deallog << "Distance between cell manifold and face manifold: "
+                    << (pref-p) << std::endl;
+            weights[0] = 0.55;
+            weights[1] = 0.45;
+            p = cell->get_manifold().get_new_point(points, weights);
+            pref = cell->face(face)->get_manifold().get_new_point(points, weights);
+            deallog << "Distance between cell manifold and face manifold: "
+                    << (pref-p) << std::endl;
+          }
+    }
+  deallog << std::endl;
+}
+
+template <int dim, int spacedim>
+void test_polar()
+{
+  deallog << "Testing with PolarManifold dim=" << dim
+          << ", spacedim="<< spacedim << std::endl;
+
+  PolarManifold<dim,spacedim> polar_manifold;
+  TransfiniteInterpolationManifold<dim,spacedim> manifold;
+
+  {
+    Triangulation<dim,spacedim> tria;
+    GridGenerator::hyper_ball (tria);
+
+    // set all entities to the transfinite manifold except for the boundary where
+    // we put the polar manifold
+    tria.set_all_manifold_ids(1);
+    tria.set_all_manifold_ids_on_boundary(0);
+    tria.set_manifold (0, polar_manifold);
+    manifold.initialize(tria);
+    tria.set_manifold (1, manifold);
+
+    do_test(tria);
+
+    // clear the triangulation, set up a similar problem again and check that it
+    // is sane.
+    tria.clear();
+    GridGenerator::hyper_ball (tria);
+    tria.set_all_manifold_ids(1);
+    tria.set_all_manifold_ids_on_boundary(0);
+    tria.set_manifold (0, polar_manifold);
+    manifold.initialize(tria);
+    tria.set_manifold (1, manifold);
+
+    do_test(tria);
+  }
+  {
+    Triangulation<dim,spacedim> tria;
+    GridGenerator::hyper_ball (tria);
+
+    // set all entities to the transfinite manifold except for the boundary where
+    // we put the polar manifold
+    tria.set_all_manifold_ids(1);
+    tria.set_all_manifold_ids_on_boundary(0);
+    tria.set_manifold (0, polar_manifold);
+    manifold.initialize(tria);
+    tria.set_manifold (1, manifold);
+
+    do_test(tria);
+  }
+}
+
+
+int main ()
+{
+  initlog();
+
+  test_polar<2,2>();
+
+  return 0;
+}
+
+
diff --git a/tests/manifold/transfinite_manifold_04.output b/tests/manifold/transfinite_manifold_04.output
new file mode 100644 (file)
index 0000000..9679a2b
--- /dev/null
@@ -0,0 +1,194 @@
+
+DEAL::Testing with PolarManifold dim=2, spacedim=2
+DEAL::Lines on cell with center: 0.00000 -0.500000
+DEAL::-0.500000 -0.500000
+DEAL::0.500000 -0.500000
+DEAL::-1.83697e-16 -1.00000
+DEAL::0.00000 -0.292893
+DEAL::Faces on cell with center: 0.00000 -0.500000
+DEAL::-0.500000 -0.500000
+DEAL::0.500000 -0.500000
+DEAL::-1.83697e-16 -1.00000
+DEAL::0.00000 -0.292893
+DEAL::Center with manifold: -9.71445e-17 -0.646447
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000
+DEAL::Lines on cell with center: -0.500000 -1.38778e-17
+DEAL::-1.00000 1.22465e-16
+DEAL::-0.292893 0.00000
+DEAL::-0.500000 -0.500000
+DEAL::-0.500000 0.500000
+DEAL::Faces on cell with center: -0.500000 -1.38778e-17
+DEAL::-1.00000 1.22465e-16
+DEAL::-0.292893 0.00000
+DEAL::-0.500000 -0.500000
+DEAL::-0.500000 0.500000
+DEAL::Center with manifold: -0.646447 5.55112e-17
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000
+DEAL::Lines on cell with center: 0.00000 0.00000
+DEAL::-0.292893 0.00000
+DEAL::0.292893 0.00000
+DEAL::0.00000 -0.292893
+DEAL::0.00000 0.292893
+DEAL::Faces on cell with center: 0.00000 0.00000
+DEAL::-0.292893 0.00000
+DEAL::0.292893 0.00000
+DEAL::0.00000 -0.292893
+DEAL::0.00000 0.292893
+DEAL::Center with manifold: 0.00000 0.00000
+DEAL::Lines on cell with center: 0.500000 0.00000
+DEAL::0.500000 -0.500000
+DEAL::0.500000 0.500000
+DEAL::1.00000 0.00000
+DEAL::0.292893 0.00000
+DEAL::Faces on cell with center: 0.500000 0.00000
+DEAL::0.500000 -0.500000
+DEAL::0.500000 0.500000
+DEAL::1.00000 0.00000
+DEAL::0.292893 0.00000
+DEAL::Center with manifold: 0.646447 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000
+DEAL::Lines on cell with center: -1.38778e-17 0.500000
+DEAL::6.12323e-17 1.00000
+DEAL::0.00000 0.292893
+DEAL::-0.500000 0.500000
+DEAL::0.500000 0.500000
+DEAL::Faces on cell with center: -1.38778e-17 0.500000
+DEAL::6.12323e-17 1.00000
+DEAL::0.00000 0.292893
+DEAL::-0.500000 0.500000
+DEAL::0.500000 0.500000
+DEAL::Center with manifold: 2.77556e-17 0.646447
+DEAL::Distance between cell manifold and face manifold: 2.22045e-16 -1.11022e-16
+DEAL::Distance between cell manifold and face manifold: 0.00000 1.11022e-16
+DEAL::
+DEAL::Lines on cell with center: 0.00000 -0.500000
+DEAL::-0.500000 -0.500000
+DEAL::0.500000 -0.500000
+DEAL::-1.83697e-16 -1.00000
+DEAL::0.00000 -0.292893
+DEAL::Faces on cell with center: 0.00000 -0.500000
+DEAL::-0.500000 -0.500000
+DEAL::0.500000 -0.500000
+DEAL::-1.83697e-16 -1.00000
+DEAL::0.00000 -0.292893
+DEAL::Center with manifold: -9.71445e-17 -0.646447
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000
+DEAL::Lines on cell with center: -0.500000 -1.38778e-17
+DEAL::-1.00000 1.22465e-16
+DEAL::-0.292893 0.00000
+DEAL::-0.500000 -0.500000
+DEAL::-0.500000 0.500000
+DEAL::Faces on cell with center: -0.500000 -1.38778e-17
+DEAL::-1.00000 1.22465e-16
+DEAL::-0.292893 0.00000
+DEAL::-0.500000 -0.500000
+DEAL::-0.500000 0.500000
+DEAL::Center with manifold: -0.646447 5.55112e-17
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000
+DEAL::Lines on cell with center: 0.00000 0.00000
+DEAL::-0.292893 0.00000
+DEAL::0.292893 0.00000
+DEAL::0.00000 -0.292893
+DEAL::0.00000 0.292893
+DEAL::Faces on cell with center: 0.00000 0.00000
+DEAL::-0.292893 0.00000
+DEAL::0.292893 0.00000
+DEAL::0.00000 -0.292893
+DEAL::0.00000 0.292893
+DEAL::Center with manifold: 0.00000 0.00000
+DEAL::Lines on cell with center: 0.500000 0.00000
+DEAL::0.500000 -0.500000
+DEAL::0.500000 0.500000
+DEAL::1.00000 0.00000
+DEAL::0.292893 0.00000
+DEAL::Faces on cell with center: 0.500000 0.00000
+DEAL::0.500000 -0.500000
+DEAL::0.500000 0.500000
+DEAL::1.00000 0.00000
+DEAL::0.292893 0.00000
+DEAL::Center with manifold: 0.646447 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000
+DEAL::Lines on cell with center: -1.38778e-17 0.500000
+DEAL::6.12323e-17 1.00000
+DEAL::0.00000 0.292893
+DEAL::-0.500000 0.500000
+DEAL::0.500000 0.500000
+DEAL::Faces on cell with center: -1.38778e-17 0.500000
+DEAL::6.12323e-17 1.00000
+DEAL::0.00000 0.292893
+DEAL::-0.500000 0.500000
+DEAL::0.500000 0.500000
+DEAL::Center with manifold: 2.77556e-17 0.646447
+DEAL::Distance between cell manifold and face manifold: 2.22045e-16 -1.11022e-16
+DEAL::Distance between cell manifold and face manifold: 0.00000 1.11022e-16
+DEAL::
+DEAL::Lines on cell with center: 0.00000 -0.500000
+DEAL::-0.500000 -0.500000
+DEAL::0.500000 -0.500000
+DEAL::-1.83697e-16 -1.00000
+DEAL::0.00000 -0.292893
+DEAL::Faces on cell with center: 0.00000 -0.500000
+DEAL::-0.500000 -0.500000
+DEAL::0.500000 -0.500000
+DEAL::-1.83697e-16 -1.00000
+DEAL::0.00000 -0.292893
+DEAL::Center with manifold: -9.71445e-17 -0.646447
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000
+DEAL::Lines on cell with center: -0.500000 -1.38778e-17
+DEAL::-1.00000 1.22465e-16
+DEAL::-0.292893 0.00000
+DEAL::-0.500000 -0.500000
+DEAL::-0.500000 0.500000
+DEAL::Faces on cell with center: -0.500000 -1.38778e-17
+DEAL::-1.00000 1.22465e-16
+DEAL::-0.292893 0.00000
+DEAL::-0.500000 -0.500000
+DEAL::-0.500000 0.500000
+DEAL::Center with manifold: -0.646447 5.55112e-17
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000
+DEAL::Lines on cell with center: 0.00000 0.00000
+DEAL::-0.292893 0.00000
+DEAL::0.292893 0.00000
+DEAL::0.00000 -0.292893
+DEAL::0.00000 0.292893
+DEAL::Faces on cell with center: 0.00000 0.00000
+DEAL::-0.292893 0.00000
+DEAL::0.292893 0.00000
+DEAL::0.00000 -0.292893
+DEAL::0.00000 0.292893
+DEAL::Center with manifold: 0.00000 0.00000
+DEAL::Lines on cell with center: 0.500000 0.00000
+DEAL::0.500000 -0.500000
+DEAL::0.500000 0.500000
+DEAL::1.00000 0.00000
+DEAL::0.292893 0.00000
+DEAL::Faces on cell with center: 0.500000 0.00000
+DEAL::0.500000 -0.500000
+DEAL::0.500000 0.500000
+DEAL::1.00000 0.00000
+DEAL::0.292893 0.00000
+DEAL::Center with manifold: 0.646447 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000
+DEAL::Distance between cell manifold and face manifold: 0.00000 0.00000
+DEAL::Lines on cell with center: -1.38778e-17 0.500000
+DEAL::6.12323e-17 1.00000
+DEAL::0.00000 0.292893
+DEAL::-0.500000 0.500000
+DEAL::0.500000 0.500000
+DEAL::Faces on cell with center: -1.38778e-17 0.500000
+DEAL::6.12323e-17 1.00000
+DEAL::0.00000 0.292893
+DEAL::-0.500000 0.500000
+DEAL::0.500000 0.500000
+DEAL::Center with manifold: 2.77556e-17 0.646447
+DEAL::Distance between cell manifold and face manifold: 2.22045e-16 -1.11022e-16
+DEAL::Distance between cell manifold and face manifold: 0.00000 1.11022e-16
+DEAL::

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.