]> https://gitweb.dealii.org/ - dealii.git/commitdiff
add GridOut::write_mesh_per_processor_as_vtu
authortcclevenger <tcleven@clemson.edu>
Tue, 5 Apr 2016 18:15:09 +0000 (14:15 -0400)
committertcclevenger <tcleven@clemson.edu>
Thu, 28 Apr 2016 15:16:25 +0000 (11:16 -0400)
Writes vtu with values "level", "subdomain", "lvl_subdomain", and "proc_writing". Also combines these into a .pvtu file. The motivation for this function is to visually debug geometric multigrid for parallel computations.

Added entry in changes.h

doc/doxygen/images/write_mesh_vtu_active.png [new file with mode: 0644]
doc/doxygen/images/write_mesh_vtu_levels.png [new file with mode: 0644]
doc/news/changes.h
include/deal.II/grid/grid_out.h
source/grid/grid_out.cc
source/grid/grid_out.inst.in
tests/grid/grid_out_per_processor_vtu_01.cc [new file with mode: 0644]
tests/grid/grid_out_per_processor_vtu_01.mpirun=3.output [new file with mode: 0644]

diff --git a/doc/doxygen/images/write_mesh_vtu_active.png b/doc/doxygen/images/write_mesh_vtu_active.png
new file mode 100644 (file)
index 0000000..13528f7
Binary files /dev/null and b/doc/doxygen/images/write_mesh_vtu_active.png differ
diff --git a/doc/doxygen/images/write_mesh_vtu_levels.png b/doc/doxygen/images/write_mesh_vtu_levels.png
new file mode 100644 (file)
index 0000000..d5d0c59
Binary files /dev/null and b/doc/doxygen/images/write_mesh_vtu_levels.png differ
index 0965897dc091df92b3318d0975efc8b1a2d28291..68767f7749083a8b237adae24ca167d9a3a3968f 100644 (file)
@@ -253,6 +253,14 @@ inconvenience this causes.
  (Daniel Shapero, 2016/04/19)
  </li>
 
+ <li> New: added function GridOut::write_mesh_per_processor_as_vtu. This allows 
+ the visualization of a parallel finite element mesh that can be seperated into each 
+ processor's owned and ghost cells. It also allows for the visualization of each level
+ of a multilevel mesh. 
+ <br>
+ (Conrad Clevenger, 2016/04/28)
+ </li>
+
  <li> Improved: The parallel loops in the deal.II Vector class for
  vector-vector operations have been revised for performance. This includes
  adjusting the minimum parallel grain size to 4096 vector entries and using an
index 93b91ddf9134ac74e7c8ae47aa3b5089c960003f..bb9c89dc32b1f25b5e7f0256988d80d57fc15eab 100644 (file)
@@ -1079,6 +1079,40 @@ public:
   void write_vtu (const Triangulation<dim,spacedim> &tria,
                   std::ostream                      &out) const;
 
+  /**
+   * Write triangulation in VTU format for each processor, and add a .pvtu file for
+   * visualization in Visit or Paraview that describes the collection of VTU files
+   * as all part of the same simulation. The output is in the form
+   * <tt>filename_without_extension.proc000*.vtu</tt> where * is 0,1,...,n_proc-1 and
+   * <tt>filename_without_extension.pvtu</tt>. The input <tt>view_levels</tt> can be
+   * set as true to view each level of a multilevel method. The input
+   * <tt>include_artificial</tt> can be set as true to view the artificial cells for
+   * each processor. Each .vtu and .pvtu file will have the attributes subdomain,
+   * level_subdomain, level, and proc_writing. The level value can be used to seperate the
+   * image into the view of the grid on each level of a multilevel method and the
+   * proc_writing value can be used to seperate the image into each processor's owned and
+   * ghost cells.
+   * This is accomplished by applying the "warp by scalar" filter in paraview
+   * to each of the values. After opening the .pvtu file of a mesh where the input
+   * <tt>view_levels</tt> is set to true, select the "warp by scalar"
+   * filter. For the "Scalars" input select <tt>proc_writing</tt> and for the "Normal" input
+   * enter in 1 0 0, then click Apply. Next select the "warp by scalar" filter again. For the
+   * "Scalars" input select <tt>level</tt> and for the "Normal" input enter in 0 1 0,
+   * then click Apply. This will give you the following image.
+   * @image html write_mesh_vtu_levels.png
+   * If the <tt>view_levels</tt> remains at false, thereby only giving the mesh for the active
+   * level, it is enough to seperate the image into the views written by different processors.
+   * This is shown in the following image where the <tt>include_artificial</tt> input is set as true.
+   * @image html write_mesh_vtu_active.png
+   * Note: Depending on the size of the mesh you may need to increase the "Scale Factor" input
+   * so that each piece does not overlap.
+   */
+  template <int dim, int spacedim>
+  void write_mesh_per_processor_as_vtu (const Triangulation<dim,spacedim> &tria,
+                                        const std::string                 &filename_without_extension,
+                                        const bool                        view_levels=false,
+                                        const bool                        include_artificial=false) const;
+
   /**
    * Write grid to @p out according to the given data format. This function
    * simply calls the appropriate <tt>write_*</tt> function.
index 324b6ea29c918181984c3f7e9facfe7b21bab9b0..f1d6f9b7890a105c20a89c0deaef34b8afe938c2 100644 (file)
 #include <deal.II/base/point.h>
 #include <deal.II/base/quadrature.h>
 #include <deal.II/base/qprojector.h>
+#include <deal.II/base/geometry_info.h>
 #include <deal.II/grid/tria.h>
 #include <deal.II/grid/tria_accessor.h>
 #include <deal.II/grid/tria_iterator.h>
 #include <deal.II/fe/mapping.h>
+#include <deal.II/numerics/data_out.h>
 
+#include <fstream>
 #include <cstring>
 #include <iomanip>
 #include <algorithm>
@@ -2526,6 +2529,105 @@ void GridOut::write_vtu (const Triangulation<dim,spacedim> &tria,
 
 
 
+template <int dim, int spacedim>
+void GridOut::write_mesh_per_processor_as_vtu (const Triangulation<dim,spacedim> &tria,
+                                               const std::string                 &filename_without_extension,
+                                               const bool                        view_levels,
+                                               const bool                        include_artificial) const
+{
+  std::vector<DataOutBase::Patch<dim,spacedim> > patches;
+  const unsigned int n_datasets=4;
+  std::vector<std::string> data_names;
+  data_names.push_back("level");
+  data_names.push_back("subdomain");
+  data_names.push_back("level_subdomain");
+  data_names.push_back("proc_writing");
+
+  const unsigned int n_q_points = GeometryInfo<dim>::vertices_per_cell;
+
+  typename Triangulation<dim, spacedim>::cell_iterator cell, endc;
+  for (cell=tria.begin(), endc=tria.end();
+       cell != endc; ++cell)
+    {
+      if (!include_artificial && cell->level_subdomain_id() ==
+          numbers::artificial_subdomain_id)
+        continue;
+      if (!view_levels && cell->has_children())
+        continue;
+
+      DataOutBase::Patch<dim,spacedim> patch;
+      patch.data.reinit(n_datasets, n_q_points);
+      patch.points_are_available = false;
+
+      for (unsigned int vertex=0; vertex<n_q_points; ++vertex)
+        {
+          patch.vertices[vertex] = cell->vertex(vertex);
+          patch.data(0,vertex) = cell->level();
+          if (!cell->has_children())
+            patch.data(1,vertex) =
+              (double)static_cast<int>(cell->subdomain_id());
+          else
+            patch.data(1,vertex) = -1.0;
+          patch.data(2,vertex) =
+            (double)static_cast<int>(cell->level_subdomain_id());
+          patch.data(3,vertex) = tria.locally_owned_subdomain();
+        }
+
+      for (unsigned int f=0; f<GeometryInfo<dim>::faces_per_cell; ++f)
+        patch.neighbors[f] = numbers::invalid_unsigned_int;
+      patches.push_back(patch);
+    }
+
+  const std::string new_file = (filename_without_extension + ".proc" +
+                                Utilities::int_to_string (tria.locally_owned_subdomain(), 4) +
+                                ".vtu");
+  std::ofstream out(new_file.c_str());
+  std::vector<std_cxx1x::tuple<unsigned int, unsigned int, std::string> > vector_data_ranges;
+  DataOutBase::VtkFlags flags;
+  DataOutBase::write_vtu (patches,
+                          data_names,
+                          vector_data_ranges,
+                          flags,
+                          out);
+  //create .pvtu record
+  if (tria.locally_owned_subdomain() == 0)
+    {
+      std::vector<std::string> filenames;
+
+      //.pvtu needs to reference the files without a relative path because it will be written
+      //in the same directory. For this, remove any paths from filename.
+      std::size_t pos = filename_without_extension.find_last_of('/');
+      if (pos == std::string::npos)
+        pos = 0;
+      else
+        pos += 1;
+      for (unsigned int i=0; i<Utilities::MPI::n_mpi_processes(MPI_COMM_WORLD); ++i)
+        filenames.push_back (filename_without_extension.substr(pos) +
+                             ".proc" + Utilities::int_to_string(i, 4) +
+                             ".vtu");
+
+      const std::string pvtu_master_filename = (filename_without_extension + ".pvtu");
+      std::ofstream pvtu_master (pvtu_master_filename.c_str());
+
+      DataOut<dim, DoFHandler<dim,spacedim> > data_out;
+      data_out.attach_triangulation (tria);
+
+      //We need a dummy vector with the names of the data values in the .vtu files
+      //in order that the .pvtu contains reference these values
+      Vector<float> dummy_vector (tria.n_active_cells());
+      data_out.add_data_vector (dummy_vector, "level");
+      data_out.add_data_vector (dummy_vector, "subdomain");
+      data_out.add_data_vector (dummy_vector, "level_subdomain");
+      data_out.add_data_vector (dummy_vector, "proc_writing");
+
+      data_out.build_patches ();
+
+      data_out.write_pvtu_record (pvtu_master, filenames);
+    }
+}
+
+
+
 unsigned int GridOut::n_boundary_faces (const Triangulation<1> &) const
 {
   return 0;
index 420ba1d8704f88da5cc168d86db96e88306c39bd..c36139e9d2aa09b1feff03d26f5b0bff50cabdcf 100644 (file)
@@ -53,8 +53,13 @@ for (deal_II_dimension : DIMENSIONS)
        std::ostream&) const;       
     template void GridOut::write_vtu
       (const Triangulation<deal_II_dimension>&,
-       std::ostream&) const;       
-       
+       std::ostream&) const;
+    template void GridOut::write_mesh_per_processor_as_vtu
+      (const Triangulation<deal_II_dimension>&,
+       const std::string&,
+       const bool,
+       const bool) const;
+
     template void GridOut::write<deal_II_dimension>
       (const Triangulation<deal_II_dimension> &,
        std::ostream &, const OutputFormat,
@@ -83,6 +88,12 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension :  SPACE_DIMENSIONS
    template void GridOut::write_vtu
       (const Triangulation<deal_II_dimension,deal_II_space_dimension>&,
        std::ostream&) const;
+   template void GridOut::write_mesh_per_processor_as_vtu
+      (const Triangulation<deal_II_dimension,deal_II_space_dimension>&,
+       const std::string&,
+       const bool,
+       const bool) const;
+
 
     template void GridOut::write<deal_II_dimension,deal_II_space_dimension>
       (const Triangulation<deal_II_dimension,deal_II_space_dimension> &,
diff --git a/tests/grid/grid_out_per_processor_vtu_01.cc b/tests/grid/grid_out_per_processor_vtu_01.cc
new file mode 100644 (file)
index 0000000..5a505b0
--- /dev/null
@@ -0,0 +1,92 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2009 - 2016 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE at
+// the top level of the deal.II distribution.
+//
+// ---------------------------------------------------------------------
+
+
+
+// check GriOut::write_mesh_per_processor_as_vtu()
+
+#include "../tests.h"
+#include <deal.II/base/logstream.h>
+#include <deal.II/base/tensor.h>
+#include <deal.II/base/geometry_info.h>
+#include <deal.II/grid/tria.h>
+#include <deal.II/distributed/tria.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/grid_out.h>
+#include <deal.II/dofs/dof_handler.h>
+#include <deal.II/dofs/dof_tools.h>
+#include <deal.II/grid/tria_accessor.h>
+#include <deal.II/grid/tria_iterator.h>
+#include <deal.II/dofs/dof_accessor.h>
+#include <deal.II/grid/grid_out.h>
+#include <deal.II/numerics/data_out.h>
+#include <deal.II/fe/fe_q.h>
+
+#include <fstream>
+
+template<int dim>
+void output(const parallel::distributed::Triangulation<dim> &tr,
+            const std::string                               &filename,
+            const bool                                      view_levels,
+            const bool                                      include_artificial)
+{
+  GridOut out;
+  out.write_mesh_per_processor_as_vtu(tr, filename, view_levels, include_artificial);
+
+  // copy the .pvtu and .vtu files
+  // into the logstream
+  int myid = Utilities::MPI::this_mpi_process (MPI_COMM_WORLD);
+  if (myid==0)
+    {
+      cat_file((std::string(filename) + ".pvtu").c_str());
+      cat_file((std::string(filename) + ".proc0000.vtu").c_str());
+    }
+  else if (myid==1)
+    cat_file((std::string(filename) + ".proc0001.vtu").c_str());
+  else if (myid==2)
+    cat_file((std::string(filename) + ".proc0002.vtu").c_str());
+  else
+    AssertThrow(false, ExcNotImplemented());
+}
+
+template<int dim>
+void test()
+{
+  unsigned int myid = Utilities::MPI::this_mpi_process (MPI_COMM_WORLD);
+  if (myid == 0)
+    deallog << "hyper_cube" << std::endl;
+
+  parallel::distributed::Triangulation<dim> tr(MPI_COMM_WORLD,
+                                               Triangulation<dim>::limit_level_difference_at_vertices,
+                                               parallel::distributed::Triangulation<dim>::construct_multigrid_hierarchy);
+  GridGenerator::hyper_cube(tr);
+  tr.refine_global(3);
+  DoFHandler<dim> dofh(tr);
+
+  output(tr, "file1", false, true);
+  output(tr, "file2", true, false);
+}
+
+
+int main(int argc, char *argv[])
+{
+  Utilities::MPI::MPI_InitFinalize mpi_initialization (argc, argv, 1);
+
+  MPILogInitAll init;
+
+  deallog.push("2d");
+  test<2>();
+  deallog.pop();
+}
diff --git a/tests/grid/grid_out_per_processor_vtu_01.mpirun=3.output b/tests/grid/grid_out_per_processor_vtu_01.mpirun=3.output
new file mode 100644 (file)
index 0000000..30fb43b
--- /dev/null
@@ -0,0 +1,278 @@
+
+DEAL:0:2d::hyper_cube
+<?xml version="1.0"?>
+<!--
+#This file was generated 
+-->
+<VTKFile type="PUnstructuredGrid" version="0.1" byte_order="LittleEndian">
+  <PUnstructuredGrid GhostLevel="0">
+    <PPointData Scalars="scalars">
+    <PDataArray type="Float64" Name="subdomain" format="ascii"/>
+    <PDataArray type="Float64" Name="proc_writing" format="ascii"/>
+    </PPointData>
+    <PPoints>
+      <PDataArray type="Float64" NumberOfComponents="3"/>
+    </PPoints>
+    <Piece Source="file1.proc0000.vtu"/>
+    <Piece Source="file1.proc0001.vtu"/>
+    <Piece Source="file1.proc0002.vtu"/>
+  </PUnstructuredGrid>
+</VTKFile>
+
+<?xml version="1.0" ?> 
+<!-- 
+# vtk DataFile Version 3.0
+#This file was generated 
+-->
+<VTKFile type="UnstructuredGrid" version="0.1" compressor="vtkZLibDataCompressor" byte_order="LittleEndian">
+<UnstructuredGrid>
+<Piece NumberOfPoints="196" NumberOfCells="49" >
+  <Points>
+    <DataArray type="Float64" NumberOfComponents="3" format="binary">
+AQAAAGASAABgEgAAbAEAAA==eNqlWDGuQzEI6824eo7AyMiQgSFDjtAp0v9RHYNge65rLEJ5pJ/P3wj59/hR+Y2f2IC/kzou77wb8G/8fN/ltz8H/AD8rM8gfORzN/0gnQ1wFENqn59nhN+h8tYd4NyH9Pyg/Ien0vODdJTgdxjxfec1gt/hJO8AfZutgxG+Ar5Kzw/ScXnXO9sPqK6Hb81+QDoq7/PM+kHnP+SNZ/0gHSV4th9QXiN4th+QjhfrYITvgO9NP0jHCX7HJH1z/+4mwe8IkneA91p2DkzCV8BX6flBOkHwOxap2513EZztZ0xnF+uwCF8BX6XnB+ns4hyYhG+Ab81+QDpRnAOT8B3wvelnkr3am/2A8i6CZ/sB6exiHRbhO+B7088q3neq+wCa64c/m/sA0lF5v0+yftDcHfLGs36G1O7L1X0A5TWCZ/cBpOPFOhjhB+BH048V7+nVeYjyToJn5yHSiWIdJuEH4EfTz0z/H/IFub2Npg==
+    </DataArray>
+  </Points>
+
+  <Cells>
+    <DataArray type="Int32" Name="connectivity" format="binary">
+AQAAABADAAAQAwAAPgEAAA==eNoNwwdXiAEAAMBPJVRm9iikZFdGJSuzorI32aPIXmVTMkIkI0KECNmJrP6Yu/cuCIKgg6GGGGZHOxluZ7sYaYRRdrW73exhT6PtZW/72M++9neAgxzoYIcYa4xDHWacwx1hvCNNMNFRjnG0Yx3nBMebZLITTXGSk011immmm+FUpzndmc5wlpnOcbZznecC55tltgvNcZG55pvnYpe4zKUud4WrXOlq17jOta53gwVudJOb3eoWt7ndne5wl7ststA97nWfxe73gIc86GGPeMyjHveEpZZ40lOe8bRnPecFz3vRMi9ZboWXveoVr1npDa970ypve8tq73jXGu9531of+NBHPrbOJz71mfU+t8GXvvCVjb7xtU2+9b3vbPaDn/zoZ7/4za+2+N0ftvrTNn/7yz/+td1//gef8Eqn
+    </DataArray>
+    <DataArray type="Int32" Name="offsets" format="binary">
+AQAAAMQAAADEAAAAcwAAAA==eNoNwxEXwmAAAMAPgiAIgiAYBEEQBEEQBINgEAyCYBAEQTAYBMEgCIIgCAZBEARBEOyHdffeNUIITVu27di1Z2TfgUNHjp04debc2IWJS1NXrs3cuHXn3tzCg0dLT569ePXm3cqHT1++/fj1Z+0fFp8TJQ==
+    </DataArray>
+    <DataArray type="UInt8" Name="types" format="binary">
+AQAAADEAAAAxAAAADAAAAA==eNrj5CQRAAArQgG6
+    </DataArray>
+  </Cells>
+  <PointData Scalars="scalars">
+    <DataArray type="Float64" Name="level" format="binary">
+AQAAACAGAAAgBgAAGwAAAA==eNpjYAADB4ZBTXOM0qP0KD1Kj9J0pgGhvDaB    </DataArray>
+    <DataArray type="Float64" Name="subdomain" format="binary">
+AQAAACAGAAAgBgAAMQAAAA==eNpjYACDAwyDmx4Fo2AAwQd7/DSh9Etr/ZTStHY/pTSt/TfY7WdwIEDTLPwBSdZmWQ==    </DataArray>
+    <DataArray type="Float64" Name="level_subdomain" format="binary">
+AQAAACAGAAAgBgAAMQAAAA==eNpjYACDAwyDmx4Fo2AAwQd7/DSh9Etr/ZTStHY/pTSt/TfY7WdwIEDTLPwBSdZmWQ==    </DataArray>
+    <DataArray type="Float64" Name="proc_writing" format="binary">
+AQAAACAGAAAgBgAAFQAAAA==eNpjYBgFo2AUjIJRMAowAQAGIAAB    </DataArray>
+  </PointData>
+ </Piece>
+ </UnstructuredGrid>
+</VTKFile>
+
+<?xml version="1.0"?>
+<!--
+#This file was generated 
+-->
+<VTKFile type="PUnstructuredGrid" version="0.1" byte_order="LittleEndian">
+  <PUnstructuredGrid GhostLevel="0">
+    <PPointData Scalars="scalars">
+    <PDataArray type="Float64" Name="subdomain" format="ascii"/>
+    <PDataArray type="Float64" Name="proc_writing" format="ascii"/>
+    </PPointData>
+    <PPoints>
+      <PDataArray type="Float64" NumberOfComponents="3"/>
+    </PPoints>
+    <Piece Source="file2.proc0000.vtu"/>
+    <Piece Source="file2.proc0001.vtu"/>
+    <Piece Source="file2.proc0002.vtu"/>
+  </PUnstructuredGrid>
+</VTKFile>
+
+<?xml version="1.0" ?> 
+<!-- 
+# vtk DataFile Version 3.0
+#This file was generated 
+-->
+<VTKFile type="UnstructuredGrid" version="0.1" compressor="vtkZLibDataCompressor" byte_order="LittleEndian">
+<UnstructuredGrid>
+<Piece NumberOfPoints="188" NumberOfCells="47" >
+  <Points>
+    <DataArray type="Float64" NumberOfComponents="3" format="binary">
+AQAAAKARAACgEQAASwEAAA==eNqtVjsSBSEI25vl6h6BkpLCwmKLPcKrbHwTPoN2xmzIALI+j7c+PKXzvY++28tQO997hmf9M50Pvl7WD4tv8PGsH6ZTzb+gdr73DM/Wl+kYfL2sHxZf4ONZP4I7/bkC32fcFeDZ/mc6XzEPK+Ab4Rt6ftal+yuo1WfzV7MfmI4FeDafTGfB75dzDdTO957h2fnDdAS+XtYPiz/g41k/A3fmpwa+z7ga4Nn5zHSsmAcN+EL4gp4fvfR/GajVZ/O12Q9MR+DXM+uH1X/Ax7N+Bu78fzXImxK+NvuB6VgxDxrwjfCt6UcvvU9m0DfnvZsBnn3/MJ1VnAMz4AvhC3p+5qX32xvk7Yz7Fn1GfCF8aeaT9e0M8Gw+mc4q3qMZ8I3wrZlP5vO9NM/Zvdz82ZznTEcCPDsPmY4GeHYeMh0L8Gw/MJ35h/8AciNj+g==
+    </DataArray>
+  </Points>
+
+  <Cells>
+    <DataArray type="Int32" Name="connectivity" format="binary">
+AQAAAPACAADwAgAAMgEAAA==eNoNwwdXDgAAAMAvSipbA9krUdlKSZnZMlIqpBAiyiqjRMosM4RoWaGyEn/O3XsXCAQCQQ51iMGGGOowhxtmhOGOcKSjHeUYxzrecUYaZYzRTnCisU5yslOc5lSnO8NZznS2c4xzrvOMd4HzTTDRhSa5yMUudYnLXG6yK0xxpWmmusp0M1xtpmtc51rXu8EsN7rJzW51i9vc7k53mO0u97jbveaY6z7z3G+B+RZ6wEMetMjDlljsEY9a6jGPe8IyT3rK056x3LNWeM5Kz3vBS160ymqveNmrXrPWGq9b501vWO8tG23wtne8513v+8Bmm3zoI5/42Kc+87ktvvClr2z1tW98a5vvbLfTDrvs9oPv/egne/zsF7/a6zf77PeH3/3pLwf87R8H/edf/wOe8USr
+    </DataArray>
+    <DataArray type="Int32" Name="offsets" format="binary">
+AQAAALwAAAC8AAAAbwAAAA==eNoNwxEXwmAAAMAPgiAIgiAIBkEQBIMgCIJgMBgEwSAIBkEwGARBEATBIAiCYBAE/bjdvXedEELXnn0HDh05NnLi1JlzYxcuXbl2Y2Jq5taduXsPFh49WVp59uLVm3cf1j59+fZj49eff1uAeBGh
+    </DataArray>
+    <DataArray type="UInt8" Name="types" format="binary">
+AQAAAC8AAAAvAAAADAAAAA==eNrj5CQJAAAn1wGo
+    </DataArray>
+  </Cells>
+  <PointData Scalars="scalars">
+    <DataArray type="Float64" Name="level" format="binary">
+AQAAAOAFAADgBQAAHgAAAA==eNpjYCAGfLAfWJrBYZSmJc0xSo/So/QwpAFNQ0DR    </DataArray>
+    <DataArray type="Float64" Name="subdomain" format="binary">
+AQAAAOAFAADgBQAAIgAAAA==eNpjYACBD/sZRukRTI+CUTCQ4IP9KE0JzeCAiwYAjdmcGQ==    </DataArray>
+    <DataArray type="Float64" Name="level_subdomain" format="binary">
+AQAAAOAFAADgBQAAKAAAAA==eNpjYKAH+GCPn2ZwIEAPckDIfwNND/XwHQXDGwz2/DN08zcABbROwQ==    </DataArray>
+    <DataArray type="Float64" Name="proc_writing" format="binary">
+AQAAAOAFAADgBQAAFAAAAA==eNpjYBgFo2AUjIJRMBwBAAXgAAE=    </DataArray>
+  </PointData>
+ </Piece>
+ </UnstructuredGrid>
+</VTKFile>
+
+
+<?xml version="1.0" ?> 
+<!-- 
+# vtk DataFile Version 3.0
+#This file was generated 
+-->
+<VTKFile type="UnstructuredGrid" version="0.1" compressor="vtkZLibDataCompressor" byte_order="LittleEndian">
+<UnstructuredGrid>
+<Piece NumberOfPoints="232" NumberOfCells="58" >
+  <Points>
+    <DataArray type="Float64" NumberOfComponents="3" format="binary">
+AQAAAMAVAADAFQAAqgEAAA==eNqlmDGSxSAMQ3MzX50jUFK6oKCg4AhbUSy7D5mxmz/R18gaxzFOvu8W1b6n//f1iQ/7/btjAb6vF/CXxXw2u/su9j+f8DNc5D11HPBXP1TnJvCoH9Jxu9f7xIu93Z/Nb0EdqgPpVLvfz6gfuv+b70k/pFMFHu0HytsEHu0H0vHHOjTBd+B70g/puMDP6KJvzueuC/yMIfIWmKvROdAFvwK/Ws5PF+cI4WdMUbcz7xT4GUvkLXB+ReswBb8Cv1rOD+msxznQBb8BvyX7gXTG4xzogu/A96Qf0hkCj/YD5Z0Cj/YD6azHOkzBd+B70g/pLLvP++g+QHN983tyHyCdavfzJOqH5m6xOx71QzpV4NF9gPI2gUf3AdLxxzo0wR/AH0k/pON2P2+j/UDn6ubPZD+QTrX7eR71Q+d/sTse9VPEe+5K9gPlbQKP9gPp+GMdmuAv4K+kH9JxgUf3AXruusCj+wDpjMc50AV/AH8k/XTxXWgk9wHKOwUe3QdIZz3WYQr+AP5I+pmP3+te+4Geuy7waD+QznicA13wF/BX0g/p/P0u+gPINeP8
+    </DataArray>
+  </Points>
+
+  <Cells>
+    <DataArray type="Int32" Name="connectivity" format="binary">
+AQAAAKADAACgAwAAbgEAAA==eNoNwwN3FgAAAMAvm1vLtm1sGYuLW9aytYxlrsWFLbfslm1by/wb3b13gUAgkMZ0pjW9GcxkRjObxWxmNbs5zGVOc5vHIPMabD7zG2IBC1rYQhaxqMUtZglLWtpSlrGs5S1nBSta2UpWsarVrWYNa1rbWtaxrvWtZwMb2thGNrGpoTYzzOa2tIWtbG1b29jO9obbwY52soud7Wo3uxthD3va2172MdK+RtnP/g50gIMc7FCHOMzhjjDakY5yjKMd6zgnON6JTnKKk53qNGc43RhnOttZznGu853nAhca6yIXu8RlLnW5K1zlSle7xnWuNc71bjDejW5yi5vdaoLb3eYOd5pkorvc7V73uM/9HvSAyR7yiIc96jFPeNyTnvKMpz3rOVM87wUvetlLXvGq173mDW9621ve8a73vecDH/rYRz7xqc995gtf+tpXvvGt733nBz+a6ic/+8VvfvW7P/zlT3/7x3/+9T8SCWit
+    </DataArray>
+    <DataArray type="Int32" Name="offsets" format="binary">
+AQAAAOgAAADoAAAAhAAAAA==eNoNw6ESAVEAAMAXBEEQBEEQBEEQBEEQBEEQBEEQBEEQBOGCIAiCIAiCIAiCIAiCIAiCIAiCIPgMuzMbCyHETZg0ZdqMWXPmLVi0ZNmKVWvWbdi0ZduOXXv2HTh05NjIiVNnzl24dOXajVt37j149OTZi1dv3n349OXbj19//gFcIBq9
+    </DataArray>
+    <DataArray type="UInt8" Name="types" format="binary">
+AQAAADoAAAA6AAAADAAAAA==eNrj5CQXAAA8YQIL
+    </DataArray>
+  </Cells>
+  <PointData Scalars="scalars">
+    <DataArray type="Float64" Name="level" format="binary">
+AQAAAEAHAABABwAAHQAAAA==eNpjYAADBwaKaI5RepQepUfpUXqUHlE0ABj2QQE=    </DataArray>
+    <DataArray type="Float64" Name="subdomain" format="binary">
+AQAAAEAHAABABwAANQAAAA==eNpjYACDAwyDmyYEaG0+pfbTGgyw/R/sR+lReuTSDA40pgnl76FuP6XlN63Nx0kDAIRbqqE=    </DataArray>
+    <DataArray type="Float64" Name="level_subdomain" format="binary">
+AQAAAEAHAABABwAANQAAAA==eNpjYACDAwyDmyYEaG0+pfbTGgyw/R/sR+lReuTSDA40pgnl76FuP6XlN63Nx0kDAIRbqqE=    </DataArray>
+    <DataArray type="Float64" Name="proc_writing" format="binary">
+AQAAAEAHAABABwAAHQAAAA==eNpjYACBD/YMo/QoPUqP0qP0KD1KE00DAEonEqg=    </DataArray>
+  </PointData>
+ </Piece>
+ </UnstructuredGrid>
+</VTKFile>
+
+<?xml version="1.0" ?> 
+<!-- 
+# vtk DataFile Version 3.0
+#This file was generated 
+-->
+<VTKFile type="UnstructuredGrid" version="0.1" compressor="vtkZLibDataCompressor" byte_order="LittleEndian">
+<UnstructuredGrid>
+<Piece NumberOfPoints="244" NumberOfCells="61" >
+  <Points>
+    <DataArray type="Float64" NumberOfComponents="3" format="binary">
+AQAAAOAWAADgFgAAqQEAAA==eNqtlzG27iAIhP+dsXWXQElpQUFh4RJu82zyzjhwMF3Il3FCIJjf73Zs+ZWun3N23zmm1K6fcxTP+kc6W+56WT9o/Sn3eNYP0kH5V6nl//AK1lVyH9NB66vc41k/Km/qJ4jv77pB4tn6RDq7mIcg/AT8lJ6feNRfKrX3c/ho1gPSmXJ/n1k/6P0r6d+sH5U3358geQvAR7MekM4u5iEIvwH/jdu/8wF0BuCr30/7xI+uJftlyD1uoE4sqYP6qPq8BtadRZ9GeORzNv1Y8XvixfnioK6iWIdOeAW8Ss+PP5qPi+Ttu+4i8ez8RTq7mIdFeAW8Ss/PerR/cMIb4K1ZD0gnin3qhEd9Opt+/NH+apG8GeCtWQ9IZxfzsAg/AT+bftaj/ecgvAPem/MU6ajc50nWD/ruDrnHs36GvNmfG8mbA96b8xrpzGIejPAB+Gj6sUf/L0Nq8/nwq1kPSEflPs+zftD8H3KPZ/0MefN/ZyRvC/CrWQ9IZxbzYIRH/1m7Oa9RXziJZ+c10olinzrhA/DRnKfI5yLx7DxFOvvR86I69P/ifwI48yA=
+    </DataArray>
+  </Points>
+
+  <Cells>
+    <DataArray type="Int32" Name="connectivity" format="binary">
+AQAAANADAADQAwAAewEAAA==eNoNwwN3FgAAAMAv23Zr2bZrLdu2bWMLa9m2bdu2bbv+QHfvXSAQCEQxmlGNbgxjGdPYxjGecY1vAhOZ0MQmMZlJTW4KU5nS1KYxnWlNbwYzmdHMBhlsFrOazRxmN6e5zGNu85rPAua3oIUsYmGLWswSFrekpSxjactazgqWt6KVrGJlqxpiqNWsbg1rWdPa1rGeda1vAxvZ0MY2sZlNbW4LW9nS1raxnW1tbwc72dHOdrGbXe1uD3vZ0972sZ997e8ABznQwQ5xmEMd7ghHOdLRjnGcYx3vBMMNc6KTnOJkI5zqNCOd7gxnOdPZznGec53vAhe50MUucZlLXe4KV7nS1a5xnWtd7wY3udHNbnGbW93uDne5093ucZ973e8BD3nQwx7xmEc97glPedLTnvGcZz3vBS950cte8ZpXve4Nb3nT297xnne97wMf+dDHPvGZT33uC1/50te+8Z1vfe8HP/nRz37xm1/97g9/+dPf/vGff/0P0ulzzw==
+    </DataArray>
+    <DataArray type="Int32" Name="offsets" format="binary">
+AQAAAPQAAAD0AAAAiQAAAA==eNoNw6ESQQEAALAXBEEQBEEQBEEQBEEQBEEQBEEQBEEQBEEQBEEQBEEQBEEQBEEQBEEQBEEQBEEQfIDtbqEgCMJGjBozbsKkKdNmzJozb8GiJctWrFqzbsOmLdt27Nqz78ChI8dOnDpz7sKlK9du3Lpz78GjJ89evHrz7sOnL99+/PrzD7NrHY0=
+    </DataArray>
+    <DataArray type="UInt8" Name="types" format="binary">
+AQAAAD0AAAA9AAAADAAAAA==eNrj5KQAAABCuAIm
+    </DataArray>
+  </Cells>
+  <PointData Scalars="scalars">
+    <DataArray type="Float64" Name="level" format="binary">
+AQAAAKAHAACgBwAAIgAAAA==eNpjYCAGfLAfWJrBYZQeyjTHKD1Kj9Kj9ChNIg0AU5ZQMQ==    </DataArray>
+    <DataArray type="Float64" Name="subdomain" format="binary">
+AQAAAKAHAACgBwAAJQAAAA==eNpjYACBD/sZRulRmmb0KMAPPtiP0qP0yKUZHEZp7DQAAc/6lQ==    </DataArray>
+    <DataArray type="Float64" Name="level_subdomain" format="binary">
+AQAAAKAHAACgBwAALQAAAA==eNpjYKAH+GCPn2ZwIEAPMCDk/sFOEwxfWtOjYFinr1F6lB7S5dOgpQEhaaDF    </DataArray>
+    <DataArray type="Float64" Name="proc_writing" format="binary">
+AQAAAKAHAACgBwAAHQAAAA==eNpjYACBD/YMo/QoPUqP0qP0KD1KDxoaANQSINw=    </DataArray>
+  </PointData>
+ </Piece>
+ </UnstructuredGrid>
+</VTKFile>
+
+
+
+<?xml version="1.0" ?> 
+<!-- 
+# vtk DataFile Version 3.0
+#This file was generated 
+-->
+<VTKFile type="UnstructuredGrid" version="0.1" compressor="vtkZLibDataCompressor" byte_order="LittleEndian">
+<UnstructuredGrid>
+<Piece NumberOfPoints="196" NumberOfCells="49" >
+  <Points>
+    <DataArray type="Float64" NumberOfComponents="3" format="binary">
+AQAAAGASAABgEgAAXAEAAA==eNqlljGywzAIRH0zXf0fwSWlChUUKjhCKjf687QwqMv6ZdkhmOh5bucdT+n595n088xx9z2/N8fdL5uH6n/cHLmcLnzPui7084Soe/rEqPXBBG/AW9JniroGvDXz0O9oQs/mMTEnM9mHJfpmwFvSx0VdA96aeagPS+jZPOTjQj/7sEXfDHhL+oSoa8BbMw/1YQs9m4d8Ytz35an/CX4Bv5I+tAfI5wW9mof2+sd7Mw/5vELP7kOqa0LP7kPymcU+mOAdeG/mIZ857v/n2Xk437OT3815IJ8X9Goeur98fDTzkM8r9Ow8UF0TenYeyGcW+2CCD+CjmYd8ptCz9wF675bQs/cB8vHiHliCd+C9mYd8XOjZ+wDV3ULP3gfIJ4p92IJ34L2Zh3yiuAeW4DfwuzkP5OPFPbAEH8BHMw/5uNCz80B1t9Cz80A+UezDFnwAH8085BP/9B+M/6Hc
+    </DataArray>
+  </Points>
+
+  <Cells>
+    <DataArray type="Int32" Name="connectivity" format="binary">
+AQAAABADAAAQAwAAPgEAAA==eNoNwwdXiAEAAMBPJVRm9iikZFdGJSuzorI32aPIXmVTMkIkI0KECNmJrP6Yu/cuCIKgg6GGGGZHOxluZ7sYaYRRdrW73exhT6PtZW/72M++9neAgxzoYIcYa4xDHWacwx1hvCNNMNFRjnG0Yx3nBMebZLITTXGSk011immmm+FUpzndmc5wlpnOcbZznecC55tltgvNcZG55pvnYpe4zKUud4WrXOlq17jOta53gwVudJOb3eoWt7ndne5wl7ststA97nWfxe73gIc86GGPeMyjHveEpZZ40lOe8bRnPecFz3vRMi9ZboWXveoVr1npDa970ypve8tq73jXGu9531of+NBHPrbOJz71mfU+t8GXvvCVjb7xtU2+9b3vbPaDn/zoZ7/4za+2+N0ftvrTNn/7yz/+td1//gef8Eqn
+    </DataArray>
+    <DataArray type="Int32" Name="offsets" format="binary">
+AQAAAMQAAADEAAAAcwAAAA==eNoNwxEXwmAAAMAPgiAIgiAYBEEQBEEQBINgEAyCYBAEQTAYBMEgCIIgCAZBEARBEOyHdffeNUIITVu27di1Z2TfgUNHjp04debc2IWJS1NXrs3cuHXn3tzCg0dLT569ePXm3cqHT1++/fj1Z+0fFp8TJQ==
+    </DataArray>
+    <DataArray type="UInt8" Name="types" format="binary">
+AQAAADEAAAAxAAAADAAAAA==eNrj5CQRAAArQgG6
+    </DataArray>
+  </Cells>
+  <PointData Scalars="scalars">
+    <DataArray type="Float64" Name="level" format="binary">
+AQAAACAGAAAgBgAAGwAAAA==eNpjYAADB4ZBTXOM0qP0KD1Kj9J0pgGhvDaB    </DataArray>
+    <DataArray type="Float64" Name="subdomain" format="binary">
+AQAAACAGAAAgBgAALAAAAA==eNpjYACDAwwjmyYEKDT/gz1l9FC3n1L30dp/tLafYv0Oo/QoPRA0AEm6eVk=    </DataArray>
+    <DataArray type="Float64" Name="level_subdomain" format="binary">
+AQAAACAGAAAgBgAALAAAAA==eNpjYACDAwwjmyYEKDT/gz1l9FC3n1L30dp/tLafYv0Oo/QoPRA0AEm6eVk=    </DataArray>
+    <DataArray type="Float64" Name="proc_writing" format="binary">
+AQAAACAGAAAgBgAAGQAAAA==eNpjYAADB4ZRepQepUfpUXqURqIBi9sxAQ==    </DataArray>
+  </PointData>
+ </Piece>
+ </UnstructuredGrid>
+</VTKFile>
+
+<?xml version="1.0" ?> 
+<!-- 
+# vtk DataFile Version 3.0
+#This file was generated 
+-->
+<VTKFile type="UnstructuredGrid" version="0.1" compressor="vtkZLibDataCompressor" byte_order="LittleEndian">
+<UnstructuredGrid>
+<Piece NumberOfPoints="188" NumberOfCells="47" >
+  <Points>
+    <DataArray type="Float64" NumberOfComponents="3" format="binary">
+AQAAAKARAACgEQAARwEAAA==eNqtljESxDAIA/Mzf/2e4NIlBQUFBU+4ykU8owgGp8nE2ROKjhA/z9cR4ynd39fsd/uQUbu/r9F61j/SifGtl/WD6sv4Xs/6QToo/zne51NnAr6a/6ljpK4AXpI6Qep2/99JeAO8NfNEOgLWq35Qf07SP1k/c9zpfyO5GeCt2Q9IJ4o5GOED8NH0Y8X5sMb7fPpcgK/Oh1NHSV0BvCR1jNStzp9Tx0ldAbwkdYLUzfbzb7yf+3xPFfDZeVLthwXqStHnIjzyaU0/qzgHqnlu3gHvl/ju/F/kuRzw3swT6Ugxh0V4NCej6Wdd+j4qee8U8Nqch0jHiu+REt4Ab00/emn/4CQ3Bbw25znSiWIOTngDvDX9+KX9lRLeAe/NfkA6VpwDSvgAfDT96KX9p5PcHPDe7AekE8UcnPAB+Gj68fT+/A9ohZjI
+    </DataArray>
+  </Points>
+
+  <Cells>
+    <DataArray type="Int32" Name="connectivity" format="binary">
+AQAAAPACAADwAgAAMgEAAA==eNoNwwdXDgAAAMAvSipbA9krUdlKSZnZMlIqpBAiyiqjRMosM4RoWaGyEn/O3XsXCAQCQQ51iMGGGOowhxtmhOGOcKSjHeUYxzrecUYaZYzRTnCisU5yslOc5lSnO8NZznS2c4xzrvOMd4HzTTDRhSa5yMUudYnLXG6yK0xxpWmmusp0M1xtpmtc51rXu8EsN7rJzW51i9vc7k53mO0u97jbveaY6z7z3G+B+RZ6wEMetMjDlljsEY9a6jGPe8IyT3rK056x3LNWeM5Kz3vBS160ymqveNmrXrPWGq9b501vWO8tG23wtne8513v+8Bmm3zoI5/42Kc+87ktvvClr2z1tW98a5vvbLfTDrvs9oPv/egne/zsF7/a6zf77PeH3/3pLwf87R8H/edf/wOe8USr
+    </DataArray>
+    <DataArray type="Int32" Name="offsets" format="binary">
+AQAAALwAAAC8AAAAbwAAAA==eNoNwxEXwmAAAMAPgiAIgiAIBkEQBIMgCIJgMBgEwSAIBkEwGARBEATBIAiCYBAE/bjdvXedEELXnn0HDh05NnLi1JlzYxcuXbl2Y2Jq5taduXsPFh49WVp59uLVm3cf1j59+fZj49eff1uAeBGh
+    </DataArray>
+    <DataArray type="UInt8" Name="types" format="binary">
+AQAAAC8AAAAvAAAADAAAAA==eNrj5CQJAAAn1wGo
+    </DataArray>
+  </Cells>
+  <PointData Scalars="scalars">
+    <DataArray type="Float64" Name="level" format="binary">
+AQAAAOAFAADgBQAAHgAAAA==eNpjYCAGfLAfWJrBYZSmJc0xSo/So/QwpAFNQ0DR    </DataArray>
+    <DataArray type="Float64" Name="subdomain" format="binary">
+AQAAAOAFAADgBQAAIgAAAA==eNpjYACBD/sZRukRTBMCH+xHaUpoBodRepQeCBoAveCvGQ==    </DataArray>
+    <DataArray type="Float64" Name="level_subdomain" format="binary">
+AQAAAOAFAADgBQAAKAAAAA==eNpjYKAH+GCPn2ZwIEBTaP5A0wT9N9D0EA/f0fgfpUdpbDQAAW5lwQ==    </DataArray>
+    <DataArray type="Float64" Name="proc_writing" format="binary">
+AQAAAOAFAADgBQAAGAAAAA==eNpjYAADB4ZRepQepUfpUXpY0QCQ5y8B    </DataArray>
+  </PointData>
+ </Piece>
+ </UnstructuredGrid>
+</VTKFile>
+
+

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.