]> https://gitweb.dealii.org/ - dealii.git/commitdiff
New GridOut::write_vtk function.
authorLuca Heltai <luca.heltai@sissa.it>
Wed, 7 Nov 2018 20:59:43 +0000 (21:59 +0100)
committerLuca Heltai <luca.heltai@sissa.it>
Thu, 15 Nov 2018 19:37:57 +0000 (20:37 +0100)
This version generates a grid that can be read back into deal.II

include/deal.II/grid/grid_out.h
source/grid/grid_out.cc
tests/grid/grid_out_vtk_01.output
tests/grid/grid_out_vtk_02.with_p4est=true.mpirun=3.output

index b221fd6a5e9501158c22deb095d63ecc76bda9bb..c8ef27c336c966745b409a24799d7940ec2ac500 100644 (file)
@@ -1146,21 +1146,42 @@ public:
                std::ostream &                      out) const;
 
   /**
-   * Write triangulation in VTK format.
+   * Write triangulation in VTK format. This function writes a
+   * UNSTRUCTURED_GRID file, that contains the following VTK cell types:
+   * VTK_HEXAHEDRON, VTK_QUAD, and VTK_LINE, depending on the template
+   * dimension.
    *
-   * Due to the way this function writes data to the output stream,
-   * the resulting output files correspond to a faithful representation
-   * of the mesh in that all cells are visible for visualization. However,
-   * the data is not in a format that allows reading this file in again
-   * through the GridIn class. This is because every vertex of the mesh is
-   * duplicated as many times as there are adjacent cells. In other words,
-   * every cell has its own, separate set of vertices that are at the
-   * same location as the vertices of other cells, but are separately
-   * numbered. If such a file is read in through the GridIn class, then
-   * that will result in a mesh that has the correct cells and vertex
-   * locations, but because the vertices are logically separate (though at
-   * the same locations) all cells are unconnected and have no neighbors
-   * across faces.
+   * In three dimensions, this function writes a file that contains
+   *
+   * - VTK_HEXAHEDRON cell types, containing the cell information of the
+   *   Triangulation
+   * - VTK_QUAD cell types, containing all boundary faces with non-zero
+   *   boundary ids, and all faces with non-flat manifold ids
+   * - VTK_LINE cell types, containing all boundary edges with non-zero
+   *   boundary ids, and all edges with non-flat manifold ids
+   *
+   * In two dimensions:
+   *
+   * - VTK_QUAD cell types, containing the cell information of the
+   *   Triangulation
+   * - VTK_LINE cell types, containing all boundary faces with non-zero
+   *   boundary ids, and all faces with non-flat manifold ids
+   *
+   * In one dimension
+   *
+   * - VTK_LINE cell types, containing the cell information of the
+   *   Triangulation
+   *
+   * The output file will contain two CELL_DATA sections, `MaterialID` and
+   * `ManifoldID`, recording for each VTK cell type the material or boundary id,
+   * and the manifold. See the
+   * [VTK file format](http://www.vtk.org/VTK/img/file-formats.pdf)
+   * documentation for an explanation of the generated output.
+   *
+   * The companion GridIn::read_vtk function can be used to read VTK files
+   * generated with this method.
+   *
+   * @author Luca Heltai, 2018
    */
   template <int dim, int spacedim>
   void
index b536dbd01adfbf8e93a34f559e1770c8a5eb9e27..43cbcac71d3b70d35a5e832b015795373b795402 100644 (file)
@@ -2921,6 +2921,8 @@ namespace
       }
   }
 
+
+
   std::vector<std::string>
   triangulation_patch_data_names()
   {
@@ -2932,6 +2934,74 @@ namespace
     v[4] = "level_subdomain";
     return v;
   }
+
+
+  /**
+   * Return all lines of a face in three dimension that have a non-standard
+   * boundary indicator (!=0), or a non-flat manifold indicator.
+   */
+  std::vector<typename Triangulation<3, 3>::active_line_iterator>
+  relevant_co_faces(const Triangulation<3, 3> &tria)
+  {
+    std::vector<typename Triangulation<3, 3>::active_line_iterator> res;
+
+    std::vector<bool> flags;
+    tria.save_user_flags_line(flags);
+    const_cast<Triangulation<3, 3> &>(tria).clear_user_flags_line();
+
+    for (auto face = tria.begin_active_face(); face != tria.end_face(); ++face)
+      for (unsigned int l = 0; l < GeometryInfo<3>::lines_per_face; ++l)
+        {
+          auto line = face->line(l);
+          if (line->user_flag_set())
+            continue;
+          else
+            line->set_user_flag();
+          if (line->manifold_id() != numbers::flat_manifold_id ||
+              (line->boundary_id() != 0 &&
+               line->boundary_id() != numbers::invalid_boundary_id))
+            res.push_back(line);
+        }
+    const_cast<Triangulation<3, 3> &>(tria).load_user_flags_line(flags);
+    return res;
+  }
+
+
+  /**
+   * Same as above, for 1 and 2 dimensional grids. Does nothing.
+   */
+  template <int dim, int spacedim>
+  std::vector<typename Triangulation<dim, spacedim>::active_line_iterator>
+  relevant_co_faces(const Triangulation<dim, spacedim> &)
+  {
+    std::vector<typename Triangulation<dim, spacedim>::active_line_iterator>
+      res;
+    return res;
+  }
+
+
+
+  /**
+   * Return all faces of a triangulation that have a non-standard
+   * boundary indicator (!=0), or a non-flat manifold indicator.
+   */
+  template <int dim, int spacedim>
+  std::vector<typename Triangulation<dim, spacedim>::active_face_iterator>
+  relevant_faces(const Triangulation<dim, spacedim> &tria)
+  {
+    std::vector<typename Triangulation<dim, spacedim>::active_face_iterator>
+      res;
+    if (dim == 1)
+      return res;
+    for (auto face = tria.begin_active_face(); face != tria.end_face(); ++face)
+      {
+        if (face->manifold_id() != numbers::flat_manifold_id ||
+            (face->boundary_id() != 0 &&
+             face->boundary_id() != numbers::invalid_boundary_id))
+          res.push_back(face);
+      }
+    return res;
+  }
 } // namespace
 
 
@@ -2943,23 +3013,142 @@ GridOut::write_vtk(const Triangulation<dim, spacedim> &tria,
 {
   AssertThrow(out, ExcIO());
 
-  // convert the cells of the triangulation into a set of patches
-  // and then have them output. since there is no data attached to
-  // the geometry, we also do not have to provide any names, identifying
-  // information, etc.
-  std::vector<DataOutBase::Patch<dim, spacedim>> patches;
-  patches.reserve(tria.n_active_cells());
-  generate_triangulation_patches(patches, tria.begin_active(), tria.end());
-  DataOutBase::write_vtk(
-    patches,
-    triangulation_patch_data_names(),
-    std::vector<
-      std::tuple<unsigned int,
-                 unsigned int,
-                 std::string,
-                 DataComponentInterpretation::DataComponentInterpretation>>(),
-    vtk_flags,
-    out);
+  // get the positions of the vertices
+  const std::vector<Point<spacedim>> &vertices = tria.get_vertices();
+
+  const auto n_vertices = vertices.size();
+
+  out << "# vtk DataFile Version 3.0" << std::endl
+      << "Triangulation generated with deal.II" << std::endl
+      << "ASCII" << std::endl
+      << "DATASET UNSTRUCTURED_GRID" << std::endl
+      << "POINTS " << n_vertices << " double" << std::endl;
+
+  // actually write the vertices.
+  for (auto v : vertices)
+    {
+      out << v;
+      for (unsigned int d = spacedim + 1; d <= 3; ++d)
+        out << " 0"; // fill with zeroes
+      out << '\n';
+    }
+
+  const auto faces    = relevant_faces(tria);
+  const auto co_faces = relevant_co_faces(tria);
+
+  // Write cells preamble
+  const int n_cells = tria.n_active_cells() + faces.size() + co_faces.size();
+
+  const int cells_size =
+    tria.n_active_cells() * (GeometryInfo<dim>::vertices_per_cell + 1) +
+    faces.size() * (GeometryInfo<dim>::vertices_per_face + 1) +
+    co_faces.size() * (3); // only in 3d, otherwise it is always zero.
+
+
+  out << std::endl << "CELLS " << n_cells << " " << cells_size << std::endl;
+  /*
+   * VTK cells:
+   *
+   * 1 VTK_VERTEX
+   * 3 VTK_LINE
+   * 9 VTK_QUAD
+   * 12 VTK_HEXAHEDRON
+   * ...
+   */
+  const int cell_type    = (dim == 1 ? 3 : dim == 2 ? 9 : 12);
+  const int face_type    = (dim == 1 ? 1 : dim == 2 ? 3 : 9);
+  const int co_face_type = (dim == 1 ? -1 : dim == 2 ? -1 : 3);
+
+  // write cells.
+  for (auto cell : tria.active_cell_iterators())
+    {
+      out << GeometryInfo<dim>::vertices_per_cell;
+      for (unsigned int i = 0; i < GeometryInfo<dim>::vertices_per_cell; ++i)
+        {
+          out << " " << cell->vertex_index(GeometryInfo<dim>::ucd_to_deal[i]);
+        }
+      out << std::endl;
+    }
+  for (auto face : faces)
+    {
+      out << GeometryInfo<dim>::vertices_per_face;
+      for (unsigned int i = 0; i < GeometryInfo<dim>::vertices_per_face; ++i)
+        {
+          out << " "
+              << face->vertex_index(
+                   GeometryInfo < (dim > 1) ? dim - 1 : dim > ::ucd_to_deal[i]);
+        }
+      out << std::endl;
+    }
+  for (auto co_face : co_faces)
+    {
+      out << 2;
+      for (unsigned int i = 0; i < 2; ++i)
+        out << " " << co_face->vertex_index(i);
+      out << std::endl;
+    }
+
+  // write cell types
+  out << std::endl << "CELL_TYPES " << n_cells << std::endl;
+  for (unsigned int i = 0; i < tria.n_active_cells(); ++i)
+    {
+      out << cell_type << " ";
+    }
+  out << std::endl;
+  for (unsigned int i = 0; i < faces.size(); ++i)
+    {
+      out << face_type << " ";
+    }
+  out << std::endl;
+  for (unsigned int i = 0; i < co_faces.size(); ++i)
+    {
+      out << co_face_type << " ";
+    }
+  out << std::endl
+      << std::endl
+      << "CELL_DATA " << n_cells << std::endl
+      << "SCALARS MaterialID int 1" << std::endl
+      << "LOOKUP_TABLE default" << std::endl;
+
+  // Now material id and boundary id
+  for (auto cell : tria.active_cell_iterators())
+    {
+      out << static_cast<int>(cell->material_id()) << " ";
+    }
+  out << std::endl;
+  for (auto face : faces)
+    {
+      out << static_cast<int>(face->boundary_id()) << " ";
+    }
+  out << std::endl;
+  for (auto co_face : co_faces)
+    {
+      out << static_cast<int>(co_face->boundary_id()) << " ";
+    }
+
+  out << std::endl
+      << std::endl
+      << "SCALARS ManifoldID int 1" << std::endl
+      << "LOOKUP_TABLE default" << std::endl;
+
+  // Now material id and boundary id
+  for (auto cell : tria.active_cell_iterators())
+    {
+      out << static_cast<int>(cell->manifold_id()) << " ";
+    }
+  out << std::endl;
+  for (auto face : faces)
+    {
+      out << static_cast<int>(face->manifold_id()) << " ";
+    }
+  out << std::endl;
+  for (auto co_face : co_faces)
+    {
+      out << static_cast<int>(co_face->manifold_id()) << " ";
+    }
+  out << std::endl;
+
+  out.flush();
 
   AssertThrow(out, ExcIO());
 }
index 2fa9835ccdb74a7917becf44355f0f87d63a8824..d1f049c4357d5a8c879e451ec302517804e1bfc3 100644 (file)
 
 # vtk DataFile Version 3.0
-#This file was generated 
+Triangulation generated with deal.II
 ASCII
 DATASET UNSTRUCTURED_GRID
-
-POINTS 6 double
+POINTS 4 double
+0.00000 0 0
 0.500000 0 0
 1.00000 0 0
-0.00000 0 0
-0.250000 0 0
 0.250000 0 0
-0.500000 0 0
 
 CELLS 3 9
-2      0       1
-2      2       3
-2      4       5
+2 1 2
+2 0 3
+2 3 1
 
 CELL_TYPES 3
- 3 3 3
-POINT_DATA 6
-SCALARS level double 1
-LOOKUP_TABLE default
-0.00000 0.00000 1.00000 1.00000 1.00000 1.00000 
-SCALARS manifold double 1
-LOOKUP_TABLE default
--1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 
-SCALARS material double 1
-LOOKUP_TABLE default
-0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 
-SCALARS subdomain double 1
+3 3 3 
+
+
+
+CELL_DATA 3
+SCALARS MaterialID int 1
 LOOKUP_TABLE default
-0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 
-SCALARS level_subdomain double 1
+0 0 0 
+
+
+
+SCALARS ManifoldID int 1
 LOOKUP_TABLE default
-0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 
+-1 -1 -1 
+
+
 # vtk DataFile Version 3.0
-#This file was generated 
+Triangulation generated with deal.II
 ASCII
 DATASET UNSTRUCTURED_GRID
-
-POINTS 6 double
+POINTS 4 double
+0.00000 0.00000 0
 0.500000 0.00000 0
 1.00000 0.00000 0
-0.00000 0.00000 0
 0.250000 0.00000 0
-0.250000 0.00000 0
-0.500000 0.00000 0
 
 CELLS 3 9
-2      0       1
-2      2       3
-2      4       5
+2 1 2
+2 0 3
+2 3 1
 
 CELL_TYPES 3
- 3 3 3
-POINT_DATA 6
-SCALARS level double 1
-LOOKUP_TABLE default
-0.00000 0.00000 1.00000 1.00000 1.00000 1.00000 
-SCALARS manifold double 1
-LOOKUP_TABLE default
--1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 
-SCALARS material double 1
-LOOKUP_TABLE default
-0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 
-SCALARS subdomain double 1
+3 3 3 
+
+
+
+CELL_DATA 3
+SCALARS MaterialID int 1
 LOOKUP_TABLE default
-0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 
-SCALARS level_subdomain double 1
+0 0 0 
+
+
+
+SCALARS ManifoldID int 1
 LOOKUP_TABLE default
-0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 
+-1 -1 -1 
+
+
 # vtk DataFile Version 3.0
-#This file was generated 
+Triangulation generated with deal.II
 ASCII
 DATASET UNSTRUCTURED_GRID
-
-POINTS 32 double
--1.50000 -0.500000 0
+POINTS 17 double
 -0.500000 -0.500000 0
--1.50000 0.500000 0
--0.500000 0.500000 0
 0.500000 -0.500000 0
-1.50000 -0.500000 0
+-0.500000 0.500000 0
 0.500000 0.500000 0
+-1.50000 -0.500000 0
+-1.50000 0.500000 0
+1.50000 -0.500000 0
 1.50000 0.500000 0
 -0.500000 -1.50000 0
 0.500000 -1.50000 0
--0.500000 -0.500000 0
-0.500000 -0.500000 0
--0.500000 0.500000 0
-0.500000 0.500000 0
 -0.500000 1.50000 0
 0.500000 1.50000 0
--0.500000 -0.500000 0
 0.00000 -0.500000 0
 -0.500000 0.00000 0
-0.00000 0.00000 0
-0.00000 -0.500000 0
-0.500000 -0.500000 0
-0.00000 0.00000 0
 0.500000 0.00000 0
--0.500000 0.00000 0
-0.00000 0.00000 0
--0.500000 0.500000 0
 0.00000 0.500000 0
 0.00000 0.00000 0
-0.500000 0.00000 0
-0.00000 0.500000 0
-0.500000 0.500000 0
 
 CELLS 8 40
-4      0       1       3       2
-4      4       5       7       6
-4      8       9       11      10
-4      12      13      15      14
-4      16      17      19      18
-4      20      21      23      22
-4      24      25      27      26
-4      28      29      31      30
+4 4 0 2 5
+4 1 6 7 3
+4 8 9 1 0
+4 2 3 11 10
+4 0 12 16 13
+4 12 1 14 16
+4 13 16 15 2
+4 16 14 3 15
 
 CELL_TYPES 8
- 9 9 9 9 9 9 9 9
-POINT_DATA 32
-SCALARS level double 1
-LOOKUP_TABLE default
-0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 
-SCALARS manifold double 1
-LOOKUP_TABLE default
--1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 
-SCALARS material double 1
-LOOKUP_TABLE default
-1.00000 1.00000 1.00000 1.00000 2.00000 2.00000 2.00000 2.00000 3.00000 3.00000 3.00000 3.00000 4.00000 4.00000 4.00000 4.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 
-SCALARS subdomain double 1
+9 9 9 9 9 9 9 9 
+
+
+
+CELL_DATA 8
+SCALARS MaterialID int 1
 LOOKUP_TABLE default
-0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 
-SCALARS level_subdomain double 1
+1 2 3 4 0 0 0 0 
+
+
+
+SCALARS ManifoldID int 1
 LOOKUP_TABLE default
-0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 
+-1 -1 -1 -1 -1 -1 -1 -1 
+
+
 # vtk DataFile Version 3.0
-#This file was generated 
+Triangulation generated with deal.II
 ASCII
 DATASET UNSTRUCTURED_GRID
-
-POINTS 32 double
--1.50000 -0.500000 0.00000
+POINTS 17 double
 -0.500000 -0.500000 0.00000
--1.50000 0.500000 0.00000
--0.500000 0.500000 0.00000
 0.500000 -0.500000 0.00000
-1.50000 -0.500000 0.00000
+-0.500000 0.500000 0.00000
 0.500000 0.500000 0.00000
+-1.50000 -0.500000 0.00000
+-1.50000 0.500000 0.00000
+1.50000 -0.500000 0.00000
 1.50000 0.500000 0.00000
 -0.500000 -1.50000 0.00000
 0.500000 -1.50000 0.00000
--0.500000 -0.500000 0.00000
-0.500000 -0.500000 0.00000
--0.500000 0.500000 0.00000
-0.500000 0.500000 0.00000
 -0.500000 1.50000 0.00000
 0.500000 1.50000 0.00000
--0.500000 -0.500000 0.00000
 0.00000 -0.500000 0.00000
 -0.500000 0.00000 0.00000
-0.00000 0.00000 0.00000
-0.00000 -0.500000 0.00000
-0.500000 -0.500000 0.00000
-0.00000 0.00000 0.00000
 0.500000 0.00000 0.00000
--0.500000 0.00000 0.00000
-0.00000 0.00000 0.00000
--0.500000 0.500000 0.00000
 0.00000 0.500000 0.00000
 0.00000 0.00000 0.00000
-0.500000 0.00000 0.00000
-0.00000 0.500000 0.00000
-0.500000 0.500000 0.00000
 
 CELLS 8 40
-4      0       1       3       2
-4      4       5       7       6
-4      8       9       11      10
-4      12      13      15      14
-4      16      17      19      18
-4      20      21      23      22
-4      24      25      27      26
-4      28      29      31      30
+4 4 0 2 5
+4 1 6 7 3
+4 8 9 1 0
+4 2 3 11 10
+4 0 12 16 13
+4 12 1 14 16
+4 13 16 15 2
+4 16 14 3 15
 
 CELL_TYPES 8
- 9 9 9 9 9 9 9 9
-POINT_DATA 32
-SCALARS level double 1
-LOOKUP_TABLE default
-0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 
-SCALARS manifold double 1
-LOOKUP_TABLE default
--1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 
-SCALARS material double 1
-LOOKUP_TABLE default
-1.00000 1.00000 1.00000 1.00000 2.00000 2.00000 2.00000 2.00000 3.00000 3.00000 3.00000 3.00000 4.00000 4.00000 4.00000 4.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 
-SCALARS subdomain double 1
+9 9 9 9 9 9 9 9 
+
+
+
+CELL_DATA 8
+SCALARS MaterialID int 1
 LOOKUP_TABLE default
-0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 
-SCALARS level_subdomain double 1
+1 2 3 4 0 0 0 0 
+
+
+
+SCALARS ManifoldID int 1
 LOOKUP_TABLE default
-0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 
+-1 -1 -1 -1 -1 -1 -1 -1 
+
+
 # vtk DataFile Version 3.0
-#This file was generated 
+Triangulation generated with deal.II
 ASCII
 DATASET UNSTRUCTURED_GRID
-
-POINTS 112 double
--1.50000 -0.500000 -0.500000
+POINTS 51 double
 -0.500000 -0.500000 -0.500000
--1.50000 0.500000 -0.500000
+0.500000 -0.500000 -0.500000
 -0.500000 0.500000 -0.500000
--1.50000 -0.500000 0.500000
+0.500000 0.500000 -0.500000
 -0.500000 -0.500000 0.500000
--1.50000 0.500000 0.500000
+0.500000 -0.500000 0.500000
 -0.500000 0.500000 0.500000
-0.500000 -0.500000 -0.500000
+0.500000 0.500000 0.500000
+-1.50000 -0.500000 -0.500000
+-1.50000 0.500000 -0.500000
+-1.50000 -0.500000 0.500000
+-1.50000 0.500000 0.500000
 1.50000 -0.500000 -0.500000
-0.500000 0.500000 -0.500000
 1.50000 0.500000 -0.500000
-0.500000 -0.500000 0.500000
 1.50000 -0.500000 0.500000
-0.500000 0.500000 0.500000
 1.50000 0.500000 0.500000
 -0.500000 -1.50000 -0.500000
-0.500000 -1.50000 -0.500000
--0.500000 -0.500000 -0.500000
-0.500000 -0.500000 -0.500000
 -0.500000 -1.50000 0.500000
+0.500000 -1.50000 -0.500000
 0.500000 -1.50000 0.500000
--0.500000 -0.500000 0.500000
-0.500000 -0.500000 0.500000
--0.500000 0.500000 -0.500000
-0.500000 0.500000 -0.500000
 -0.500000 1.50000 -0.500000
-0.500000 1.50000 -0.500000
--0.500000 0.500000 0.500000
-0.500000 0.500000 0.500000
 -0.500000 1.50000 0.500000
+0.500000 1.50000 -0.500000
 0.500000 1.50000 0.500000
 -0.500000 -0.500000 -1.50000
 0.500000 -0.500000 -1.50000
 -0.500000 0.500000 -1.50000
 0.500000 0.500000 -1.50000
--0.500000 -0.500000 -0.500000
-0.500000 -0.500000 -0.500000
--0.500000 0.500000 -0.500000
-0.500000 0.500000 -0.500000
--0.500000 -0.500000 0.500000
-0.500000 -0.500000 0.500000
--0.500000 0.500000 0.500000
-0.500000 0.500000 0.500000
 -0.500000 -0.500000 1.50000
 0.500000 -0.500000 1.50000
 -0.500000 0.500000 1.50000
 0.500000 0.500000 1.50000
--0.500000 -0.500000 -0.500000
 0.00000 -0.500000 -0.500000
 -0.500000 0.00000 -0.500000
-0.00000 0.00000 -0.500000
 -0.500000 -0.500000 0.00000
-0.00000 -0.500000 0.00000
--0.500000 0.00000 0.00000
-0.00000 0.00000 6.93889e-18
-0.00000 -0.500000 -0.500000
-0.500000 -0.500000 -0.500000
-0.00000 0.00000 -0.500000
 0.500000 0.00000 -0.500000
-0.00000 -0.500000 0.00000
 0.500000 -0.500000 0.00000
-0.00000 0.00000 6.93889e-18
-0.500000 0.00000 0.00000
--0.500000 0.00000 -0.500000
-0.00000 0.00000 -0.500000
--0.500000 0.500000 -0.500000
 0.00000 0.500000 -0.500000
--0.500000 0.00000 0.00000
-0.00000 0.00000 6.93889e-18
 -0.500000 0.500000 0.00000
-0.00000 0.500000 0.00000
-0.00000 0.00000 -0.500000
-0.500000 0.00000 -0.500000
-0.00000 0.500000 -0.500000
-0.500000 0.500000 -0.500000
-0.00000 0.00000 6.93889e-18
-0.500000 0.00000 0.00000
-0.00000 0.500000 0.00000
 0.500000 0.500000 0.00000
--0.500000 -0.500000 0.00000
-0.00000 -0.500000 0.00000
--0.500000 0.00000 0.00000
-0.00000 0.00000 6.93889e-18
--0.500000 -0.500000 0.500000
 0.00000 -0.500000 0.500000
 -0.500000 0.00000 0.500000
-0.00000 0.00000 0.500000
-0.00000 -0.500000 0.00000
-0.500000 -0.500000 0.00000
-0.00000 0.00000 6.93889e-18
-0.500000 0.00000 0.00000
-0.00000 -0.500000 0.500000
-0.500000 -0.500000 0.500000
-0.00000 0.00000 0.500000
 0.500000 0.00000 0.500000
--0.500000 0.00000 0.00000
-0.00000 0.00000 6.93889e-18
--0.500000 0.500000 0.00000
-0.00000 0.500000 0.00000
--0.500000 0.00000 0.500000
-0.00000 0.00000 0.500000
--0.500000 0.500000 0.500000
 0.00000 0.500000 0.500000
-0.00000 0.00000 6.93889e-18
+0.00000 -0.500000 0.00000
+0.00000 0.00000 -0.500000
+-0.500000 0.00000 0.00000
 0.500000 0.00000 0.00000
 0.00000 0.500000 0.00000
-0.500000 0.500000 0.00000
 0.00000 0.00000 0.500000
-0.500000 0.00000 0.500000
-0.00000 0.500000 0.500000
-0.500000 0.500000 0.500000
+0.00000 0.00000 0.00000
 
 CELLS 14 126
-8      0       1       3       2       4       5       7       6
-8      8       9       11      10      12      13      15      14
-8      16      17      19      18      20      21      23      22
-8      24      25      27      26      28      29      31      30
-8      32      33      35      34      36      37      39      38
-8      40      41      43      42      44      45      47      46
-8      48      49      51      50      52      53      55      54
-8      56      57      59      58      60      61      63      62
-8      64      65      67      66      68      69      71      70
-8      72      73      75      74      76      77      79      78
-8      80      81      83      82      84      85      87      86
-8      88      89      91      90      92      93      95      94
-8      96      97      99      98      100     101     103     102
-8      104     105     107     106     108     109     111     110
+8 8 0 4 10 9 2 6 11
+8 1 12 14 5 3 13 15 7
+8 16 18 19 17 0 1 5 4
+8 2 3 7 6 20 22 23 21
+8 24 25 1 0 26 27 3 2
+8 4 5 29 28 6 7 31 30
+8 0 32 44 34 33 45 50 46
+8 32 1 36 44 45 35 47 50
+8 33 45 50 46 2 37 48 38
+8 45 35 47 50 37 3 39 48
+8 34 44 40 4 46 50 49 41
+8 44 36 5 40 50 47 42 49
+8 46 50 49 41 38 48 43 6
+8 50 47 42 49 48 39 7 43
 
 CELL_TYPES 14
- 12 12 12 12 12 12 12 12 12 12 12 12 12 12
-POINT_DATA 112
-SCALARS level double 1
-LOOKUP_TABLE default
-0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 
-SCALARS manifold double 1
-LOOKUP_TABLE default
--1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 
-SCALARS material double 1
-LOOKUP_TABLE default
-1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 3.00000 3.00000 3.00000 3.00000 3.00000 3.00000 3.00000 3.00000 4.00000 4.00000 4.00000 4.00000 4.00000 4.00000 4.00000 4.00000 5.00000 5.00000 5.00000 5.00000 5.00000 5.00000 5.00000 5.00000 6.00000 6.00000 6.00000 6.00000 6.00000 6.00000 6.00000 6.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 
-SCALARS subdomain double 1
+12 12 12 12 12 12 12 12 12 12 12 12 12 12 
+
+
+
+CELL_DATA 14
+SCALARS MaterialID int 1
 LOOKUP_TABLE default
-0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 
-SCALARS level_subdomain double 1
+1 2 3 4 5 6 0 0 0 0 0 0 0 0 
+
+
+
+SCALARS ManifoldID int 1
 LOOKUP_TABLE default
-0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 
+-1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 
+
+
index 023f39c87b26930a457d5f33388db939473c07fd..0c40d5d644a757993a7da196a671586c60485496 100644 (file)
 
 DEAL:0:2d::hyper_cube
 # vtk DataFile Version 3.0
-#This file was generated 
+Triangulation generated with deal.II
 ASCII
 DATASET UNSTRUCTURED_GRID
-
-POINTS 64 double
+POINTS 25 double
 0.00000 0.00000 0
-0.250000 0.00000 0
-0.00000 0.250000 0
-0.250000 0.250000 0
-0.250000 0.00000 0
+1.00000 0.00000 0
+0.00000 1.00000 0
+1.00000 1.00000 0
 0.500000 0.00000 0
-0.250000 0.250000 0
-0.500000 0.250000 0
-0.00000 0.250000 0
-0.250000 0.250000 0
 0.00000 0.500000 0
-0.250000 0.500000 0
-0.250000 0.250000 0
-0.500000 0.250000 0
-0.250000 0.500000 0
+1.00000 0.500000 0
+0.500000 1.00000 0
 0.500000 0.500000 0
-0.500000 0.00000 0
-0.750000 0.00000 0
-0.500000 0.250000 0
-0.750000 0.250000 0
+0.250000 0.00000 0
 0.750000 0.00000 0
-1.00000 0.00000 0
-0.750000 0.250000 0
+0.00000 0.250000 0
+0.00000 0.750000 0
 1.00000 0.250000 0
+1.00000 0.750000 0
+0.250000 1.00000 0
+0.750000 1.00000 0
 0.500000 0.250000 0
-0.750000 0.250000 0
-0.500000 0.500000 0
+0.500000 0.750000 0
+0.250000 0.500000 0
 0.750000 0.500000 0
+0.250000 0.250000 0
 0.750000 0.250000 0
-1.00000 0.250000 0
-0.750000 0.500000 0
-1.00000 0.500000 0
-0.00000 0.500000 0
-0.250000 0.500000 0
-0.00000 0.750000 0
-0.250000 0.750000 0
-0.250000 0.500000 0
-0.500000 0.500000 0
-0.250000 0.750000 0
-0.500000 0.750000 0
-0.00000 0.750000 0
-0.250000 0.750000 0
-0.00000 1.00000 0
-0.250000 1.00000 0
 0.250000 0.750000 0
-0.500000 0.750000 0
-0.250000 1.00000 0
-0.500000 1.00000 0
-0.500000 0.500000 0
-0.750000 0.500000 0
-0.500000 0.750000 0
 0.750000 0.750000 0
-0.750000 0.500000 0
-1.00000 0.500000 0
-0.750000 0.750000 0
-1.00000 0.750000 0
-0.500000 0.750000 0
-0.750000 0.750000 0
-0.500000 1.00000 0
-0.750000 1.00000 0
-0.750000 0.750000 0
-1.00000 0.750000 0
-0.750000 1.00000 0
-1.00000 1.00000 0
 
 CELLS 16 80
-4      0       1       3       2
-4      4       5       7       6
-4      8       9       11      10
-4      12      13      15      14
-4      16      17      19      18
-4      20      21      23      22
-4      24      25      27      26
-4      28      29      31      30
-4      32      33      35      34
-4      36      37      39      38
-4      40      41      43      42
-4      44      45      47      46
-4      48      49      51      50
-4      52      53      55      54
-4      56      57      59      58
-4      60      61      63      62
+4 0 9 21 11
+4 9 4 17 21
+4 11 21 19 5
+4 21 17 8 19
+4 4 10 22 17
+4 10 1 13 22
+4 17 22 20 8
+4 22 13 6 20
+4 5 19 23 12
+4 19 8 18 23
+4 12 23 15 2
+4 23 18 7 15
+4 8 20 24 18
+4 20 6 14 24
+4 18 24 16 7
+4 24 14 3 16
 
 CELL_TYPES 16
- 9 9 9 9 9 9 9 9 9 9 9 9 9 9 9 9
-POINT_DATA 64
-SCALARS level double 1
-LOOKUP_TABLE default
-2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 
-SCALARS manifold double 1
-LOOKUP_TABLE default
--1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 
-SCALARS material double 1
-LOOKUP_TABLE default
-0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 
-SCALARS subdomain double 1
+9 9 9 9 9 9 9 9 9 9 9 9 9 9 9 9 
+
+
+
+CELL_DATA 16
+SCALARS MaterialID int 1
 LOOKUP_TABLE default
-0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 1.00000 1.00000 1.00000 1.00000 -2.00000 -2.00000 -2.00000 -2.00000 1.00000 1.00000 1.00000 1.00000 -2.00000 -2.00000 -2.00000 -2.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 -2.00000 -2.00000 -2.00000 -2.00000 -2.00000 -2.00000 -2.00000 -2.00000 2.00000 2.00000 2.00000 2.00000 -2.00000 -2.00000 -2.00000 -2.00000 -2.00000 -2.00000 -2.00000 -2.00000 -2.00000 -2.00000 -2.00000 -2.00000 
-SCALARS level_subdomain double 1
+0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 
+
+
+
+SCALARS ManifoldID int 1
 LOOKUP_TABLE default
-0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 1.00000 1.00000 1.00000 1.00000 -2.00000 -2.00000 -2.00000 -2.00000 1.00000 1.00000 1.00000 1.00000 -2.00000 -2.00000 -2.00000 -2.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 -2.00000 -2.00000 -2.00000 -2.00000 -2.00000 -2.00000 -2.00000 -2.00000 2.00000 2.00000 2.00000 2.00000 -2.00000 -2.00000 -2.00000 -2.00000 -2.00000 -2.00000 -2.00000 -2.00000 -2.00000 -2.00000 -2.00000 -2.00000 
+-1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 
+
+
 DEAL:0:2d::ID = 0
 DEAL:0:2d::0
 DEAL:0:2d::0112
@@ -115,112 +74,71 @@ DEAL:0:2d::
 DEAL:0:2d::my levels: 3<= global levels:3
 
 # vtk DataFile Version 3.0
-#This file was generated 
+Triangulation generated with deal.II
 ASCII
 DATASET UNSTRUCTURED_GRID
-
-POINTS 64 double
+POINTS 25 double
 0.00000 0.00000 0
-0.250000 0.00000 0
-0.00000 0.250000 0
-0.250000 0.250000 0
-0.250000 0.00000 0
+1.00000 0.00000 0
+0.00000 1.00000 0
+1.00000 1.00000 0
 0.500000 0.00000 0
-0.250000 0.250000 0
-0.500000 0.250000 0
-0.00000 0.250000 0
-0.250000 0.250000 0
 0.00000 0.500000 0
-0.250000 0.500000 0
-0.250000 0.250000 0
-0.500000 0.250000 0
-0.250000 0.500000 0
+1.00000 0.500000 0
+0.500000 1.00000 0
 0.500000 0.500000 0
-0.500000 0.00000 0
-0.750000 0.00000 0
-0.500000 0.250000 0
-0.750000 0.250000 0
+0.250000 0.00000 0
 0.750000 0.00000 0
-1.00000 0.00000 0
-0.750000 0.250000 0
+0.00000 0.250000 0
+0.00000 0.750000 0
 1.00000 0.250000 0
+1.00000 0.750000 0
+0.250000 1.00000 0
+0.750000 1.00000 0
 0.500000 0.250000 0
-0.750000 0.250000 0
-0.500000 0.500000 0
+0.500000 0.750000 0
+0.250000 0.500000 0
 0.750000 0.500000 0
+0.250000 0.250000 0
 0.750000 0.250000 0
-1.00000 0.250000 0
-0.750000 0.500000 0
-1.00000 0.500000 0
-0.00000 0.500000 0
-0.250000 0.500000 0
-0.00000 0.750000 0
-0.250000 0.750000 0
-0.250000 0.500000 0
-0.500000 0.500000 0
 0.250000 0.750000 0
-0.500000 0.750000 0
-0.00000 0.750000 0
-0.250000 0.750000 0
-0.00000 1.00000 0
-0.250000 1.00000 0
-0.250000 0.750000 0
-0.500000 0.750000 0
-0.250000 1.00000 0
-0.500000 1.00000 0
-0.500000 0.500000 0
-0.750000 0.500000 0
-0.500000 0.750000 0
-0.750000 0.750000 0
-0.750000 0.500000 0
-1.00000 0.500000 0
 0.750000 0.750000 0
-1.00000 0.750000 0
-0.500000 0.750000 0
-0.750000 0.750000 0
-0.500000 1.00000 0
-0.750000 1.00000 0
-0.750000 0.750000 0
-1.00000 0.750000 0
-0.750000 1.00000 0
-1.00000 1.00000 0
 
 CELLS 16 80
-4      0       1       3       2
-4      4       5       7       6
-4      8       9       11      10
-4      12      13      15      14
-4      16      17      19      18
-4      20      21      23      22
-4      24      25      27      26
-4      28      29      31      30
-4      32      33      35      34
-4      36      37      39      38
-4      40      41      43      42
-4      44      45      47      46
-4      48      49      51      50
-4      52      53      55      54
-4      56      57      59      58
-4      60      61      63      62
+4 0 9 21 11
+4 9 4 17 21
+4 11 21 19 5
+4 21 17 8 19
+4 4 10 22 17
+4 10 1 13 22
+4 17 22 20 8
+4 22 13 6 20
+4 5 19 23 12
+4 19 8 18 23
+4 12 23 15 2
+4 23 18 7 15
+4 8 20 24 18
+4 20 6 14 24
+4 18 24 16 7
+4 24 14 3 16
 
 CELL_TYPES 16
- 9 9 9 9 9 9 9 9 9 9 9 9 9 9 9 9
-POINT_DATA 64
-SCALARS level double 1
-LOOKUP_TABLE default
-2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 
-SCALARS manifold double 1
-LOOKUP_TABLE default
--1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 
-SCALARS material double 1
-LOOKUP_TABLE default
-0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 
-SCALARS subdomain double 1
+9 9 9 9 9 9 9 9 9 9 9 9 9 9 9 9 
+
+
+
+CELL_DATA 16
+SCALARS MaterialID int 1
 LOOKUP_TABLE default
--2.00000 -2.00000 -2.00000 -2.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 -2.00000 -2.00000 -2.00000 -2.00000 
-SCALARS level_subdomain double 1
+0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 
+
+
+
+SCALARS ManifoldID int 1
 LOOKUP_TABLE default
--2.00000 -2.00000 -2.00000 -2.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 -2.00000 -2.00000 -2.00000 -2.00000 
+-1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 
+
+
 DEAL:1:2d::ID = 1
 DEAL:1:2d::0
 DEAL:1:2d::0112
@@ -229,112 +147,71 @@ DEAL:1:2d::
 
 
 # vtk DataFile Version 3.0
-#This file was generated 
+Triangulation generated with deal.II
 ASCII
 DATASET UNSTRUCTURED_GRID
-
-POINTS 64 double
+POINTS 25 double
 0.00000 0.00000 0
-0.250000 0.00000 0
-0.00000 0.250000 0
-0.250000 0.250000 0
-0.250000 0.00000 0
+1.00000 0.00000 0
+0.00000 1.00000 0
+1.00000 1.00000 0
 0.500000 0.00000 0
-0.250000 0.250000 0
-0.500000 0.250000 0
-0.00000 0.250000 0
-0.250000 0.250000 0
 0.00000 0.500000 0
-0.250000 0.500000 0
-0.250000 0.250000 0
-0.500000 0.250000 0
-0.250000 0.500000 0
+1.00000 0.500000 0
+0.500000 1.00000 0
 0.500000 0.500000 0
-0.500000 0.00000 0
-0.750000 0.00000 0
-0.500000 0.250000 0
-0.750000 0.250000 0
+0.250000 0.00000 0
 0.750000 0.00000 0
-1.00000 0.00000 0
-0.750000 0.250000 0
+0.00000 0.250000 0
+0.00000 0.750000 0
 1.00000 0.250000 0
+1.00000 0.750000 0
+0.250000 1.00000 0
+0.750000 1.00000 0
 0.500000 0.250000 0
-0.750000 0.250000 0
-0.500000 0.500000 0
+0.500000 0.750000 0
+0.250000 0.500000 0
 0.750000 0.500000 0
+0.250000 0.250000 0
 0.750000 0.250000 0
-1.00000 0.250000 0
-0.750000 0.500000 0
-1.00000 0.500000 0
-0.00000 0.500000 0
-0.250000 0.500000 0
-0.00000 0.750000 0
-0.250000 0.750000 0
-0.250000 0.500000 0
-0.500000 0.500000 0
-0.250000 0.750000 0
-0.500000 0.750000 0
-0.00000 0.750000 0
-0.250000 0.750000 0
-0.00000 1.00000 0
-0.250000 1.00000 0
 0.250000 0.750000 0
-0.500000 0.750000 0
-0.250000 1.00000 0
-0.500000 1.00000 0
-0.500000 0.500000 0
-0.750000 0.500000 0
-0.500000 0.750000 0
-0.750000 0.750000 0
-0.750000 0.500000 0
-1.00000 0.500000 0
-0.750000 0.750000 0
-1.00000 0.750000 0
-0.500000 0.750000 0
-0.750000 0.750000 0
-0.500000 1.00000 0
-0.750000 1.00000 0
 0.750000 0.750000 0
-1.00000 0.750000 0
-0.750000 1.00000 0
-1.00000 1.00000 0
 
 CELLS 16 80
-4      0       1       3       2
-4      4       5       7       6
-4      8       9       11      10
-4      12      13      15      14
-4      16      17      19      18
-4      20      21      23      22
-4      24      25      27      26
-4      28      29      31      30
-4      32      33      35      34
-4      36      37      39      38
-4      40      41      43      42
-4      44      45      47      46
-4      48      49      51      50
-4      52      53      55      54
-4      56      57      59      58
-4      60      61      63      62
+4 0 9 21 11
+4 9 4 17 21
+4 11 21 19 5
+4 21 17 8 19
+4 4 10 22 17
+4 10 1 13 22
+4 17 22 20 8
+4 22 13 6 20
+4 5 19 23 12
+4 19 8 18 23
+4 12 23 15 2
+4 23 18 7 15
+4 8 20 24 18
+4 20 6 14 24
+4 18 24 16 7
+4 24 14 3 16
 
 CELL_TYPES 16
- 9 9 9 9 9 9 9 9 9 9 9 9 9 9 9 9
-POINT_DATA 64
-SCALARS level double 1
-LOOKUP_TABLE default
-2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 
-SCALARS manifold double 1
-LOOKUP_TABLE default
--1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 -1.00000 
-SCALARS material double 1
-LOOKUP_TABLE default
-0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 
-SCALARS subdomain double 1
+9 9 9 9 9 9 9 9 9 9 9 9 9 9 9 9 
+
+
+
+CELL_DATA 16
+SCALARS MaterialID int 1
 LOOKUP_TABLE default
--2.00000 -2.00000 -2.00000 -2.00000 -2.00000 -2.00000 -2.00000 -2.00000 -2.00000 -2.00000 -2.00000 -2.00000 0.00000 0.00000 0.00000 0.00000 -2.00000 -2.00000 -2.00000 -2.00000 -2.00000 -2.00000 -2.00000 -2.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 -2.00000 -2.00000 -2.00000 -2.00000 1.00000 1.00000 1.00000 1.00000 -2.00000 -2.00000 -2.00000 -2.00000 1.00000 1.00000 1.00000 1.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 
-SCALARS level_subdomain double 1
+0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 
+
+
+
+SCALARS ManifoldID int 1
 LOOKUP_TABLE default
--2.00000 -2.00000 -2.00000 -2.00000 -2.00000 -2.00000 -2.00000 -2.00000 -2.00000 -2.00000 -2.00000 -2.00000 0.00000 0.00000 0.00000 0.00000 -2.00000 -2.00000 -2.00000 -2.00000 -2.00000 -2.00000 -2.00000 -2.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 1.00000 -2.00000 -2.00000 -2.00000 -2.00000 1.00000 1.00000 1.00000 1.00000 -2.00000 -2.00000 -2.00000 -2.00000 1.00000 1.00000 1.00000 1.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 2.00000 
+-1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 
+
+
 DEAL:2:2d::ID = 2
 DEAL:2:2d::0
 DEAL:2:2d::0112

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.