]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Write and read VTU format that contains serialized triangulation. 9376/head
authorLuca Heltai <luca.heltai@sissa.it>
Fri, 17 Jan 2020 17:23:15 +0000 (18:23 +0100)
committerLuca Heltai <luca.heltai@sissa.it>
Wed, 8 Apr 2020 09:00:25 +0000 (11:00 +0200)
doc/news/changes/minor/20200117LucaHeltaiNicolaGiuliani_bis [new file with mode: 0644]
include/deal.II/grid/grid_in.h
include/deal.II/grid/grid_out.h
source/grid/grid_in.cc
source/grid/grid_out.cc
tests/grid/grid_in_out_vtu_01.cc [new file with mode: 0644]
tests/grid/grid_in_out_vtu_01.output [new file with mode: 0644]

diff --git a/doc/news/changes/minor/20200117LucaHeltaiNicolaGiuliani_bis b/doc/news/changes/minor/20200117LucaHeltaiNicolaGiuliani_bis
new file mode 100644 (file)
index 0000000..4fb5d3e
--- /dev/null
@@ -0,0 +1,4 @@
+New: GridOut::write_vtu() and GridIn::read_vtu() allow to save and restore a locally refined triangulation,
+using a xml section in the vtu file that is ignored by vtu readers.
+<br>
+(Luca Heltai, Nicola Giuliani, 2020/01/17)
index 04dadbf15a9d5d7548cea03cd1b61a700d5fab67..84e2f356340aaee5aa7b6191cede748aab1907a7 100644 (file)
@@ -333,6 +333,8 @@ public:
     tecplot,
     /// Use read_vtk()
     vtk,
+    /// Use read_vtu()
+    vtu,
     /// Use read_assimp()
     assimp,
   };
@@ -404,6 +406,26 @@ public:
   void
   read_vtk(std::istream &in);
 
+  /**
+   * Read grid data from a unstructured vtu file, saved by deal.II using
+   * GridOut::write_vtu(), with the flag
+   * GridOutFlags::Vtu::serialize_triangulation set to true.
+   *
+   * Notice that this function does not support reading in arbitrary vtu files,
+   * but only files that were written by deal.II itself, using the function
+   * GridOut::write_vtu and setting GridOutFlags::Vtu::serialize_triangulation
+   * to true.
+   *
+   * When this flag is set to true, the generated vtu file contains the
+   * triangulation in a xml section which is ignored by vtu general vtu readers.
+   * If this section is absent, an exception is thrown.
+   *
+   * @author Luca Heltai, Nicola Giuliani, 2020
+   */
+  void
+  read_vtu(std::istream &in);
+
+
   /**
    * Read grid data from an unv file as generated by the Salome mesh
    * generator. Numerical data is ignored.
index 22c8472c725c483c5f8e4168530d6f7361975337..3dc7ef51689d102afcded1b367cde940b687b83a 100644 (file)
@@ -921,12 +921,23 @@ namespace GridOutFlags
 
   /**
    * Flags for grid output in Vtu format. These flags are the same as those
-   * declared in DataOutBase::VtkFlags.
+   * declared in DataOutBase::VtkFlags, with the addition of a flag that
+   * determines if you want to add a entry in the vtu file (which is really
+   * a xml file) containing the entire serialization of the triangulation.
    *
    * @ingroup output
    */
   struct Vtu : public DataOutBase::VtkFlags
-  {};
+  {
+    Vtu(const bool serialize_triangulation = false)
+      : serialize_triangulation(serialize_triangulation)
+    {}
+
+    /**
+     * Add to the vtu file also the serialized triangulation.
+     */
+    bool serialize_triangulation;
+  };
 } // namespace GridOutFlags
 
 
index 105ce03f79a11d95fe998532a406f7b6f8e0fad2..bbff7456d8cd274694b0a4387ecd00da6976a45a 100644 (file)
 #include <deal.II/grid/grid_tools.h>
 #include <deal.II/grid/tria.h>
 
+#include <boost/archive/binary_iarchive.hpp>
 #include <boost/io/ios_state.hpp>
+#include <boost/property_tree/ptree.hpp>
+#include <boost/property_tree/xml_parser.hpp>
+#include <boost/serialization/serialization.hpp>
 
 #include <algorithm>
 #include <cctype>
@@ -495,6 +499,34 @@ GridIn<dim, spacedim>::read_vtk(std::istream &in)
 }
 
 
+
+template <int dim, int spacedim>
+void
+GridIn<dim, spacedim>::read_vtu(std::istream &in)
+{
+  namespace pt = boost::property_tree;
+  pt::ptree tree;
+  pt::read_xml(in, tree);
+  auto section = tree.get_optional<std::string>("VTKFile.dealiiData");
+
+  AssertThrow(section,
+              ExcMessage(
+                "While reading a VTU file, failed to find dealiiData section. "
+                "Notice that we can only read grid files in .vtu format that "
+                "were created by the deal.II library, using a call to "
+                "GridOut::write_vtu(), where the flag "
+                "GridOutFlags::Vtu::serialize_triangulation is set to true."));
+
+  const auto decoded =
+    Utilities::decode_base64({section->begin(), section->end() - 1});
+  const auto string_archive =
+    Utilities::decompress({decoded.begin(), decoded.end()});
+  std::istringstream              in_stream(string_archive);
+  boost::archive::binary_iarchive ia(in_stream);
+  tria->load(ia, 0);
+}
+
+
 template <int dim, int spacedim>
 void
 GridIn<dim, spacedim>::read_unv(std::istream &in)
@@ -3376,6 +3408,10 @@ GridIn<dim, spacedim>::read(std::istream &in, Format format)
         read_vtk(in);
         return;
 
+      case vtu:
+        read_vtu(in);
+        return;
+
       case unv:
         read_unv(in);
         return;
@@ -3430,6 +3466,8 @@ GridIn<dim, spacedim>::default_suffix(const Format format)
         return ".msh";
       case vtk:
         return ".vtk";
+      case vtu:
+        return ".vtu";
       case unv:
         return ".unv";
       case ucd:
@@ -3467,6 +3505,9 @@ GridIn<dim, spacedim>::parse_format(const std::string &format_name)
   if (format_name == "vtk")
     return vtk;
 
+  if (format_name == "vtu")
+    return vtu;
+
   // This is also the typical extension of Abaqus input files.
   if (format_name == "inp")
     return ucd;
@@ -3512,7 +3553,7 @@ template <int dim, int spacedim>
 std::string
 GridIn<dim, spacedim>::get_format_names()
 {
-  return "dbmesh|msh|unv|vtk|ucd|abaqus|xda|netcdf|tecplot|assimp";
+  return "dbmesh|msh|unv|vtk|vtu|ucd|abaqus|xda|netcdf|tecplot|assimp";
 }
 
 
index e2f173a744aabcfb394a907bafb319704bdd4b15..f03760be66026dec3f452cbc7132a32a5190759f 100644 (file)
@@ -29,6 +29,9 @@
 
 #include <deal.II/numerics/data_out.h>
 
+#include <boost/algorithm/string.hpp>
+#include <boost/archive/binary_oarchive.hpp>
+
 #include <algorithm>
 #include <cmath>
 #include <cstring>
@@ -38,7 +41,6 @@
 #include <list>
 #include <set>
 
-
 DEAL_II_NAMESPACE_OPEN
 
 
@@ -3369,7 +3371,9 @@ GridOut::write_vtu(const Triangulation<dim, spacedim> &tria,
   std::vector<DataOutBase::Patch<dim, spacedim>> patches;
   patches.reserve(tria.n_active_cells());
   generate_triangulation_patches(patches, tria.begin_active(), tria.end());
-  DataOutBase::write_vtu(
+
+  DataOutBase::write_vtu_header(out, vtu_flags);
+  DataOutBase::write_vtu_main(
     patches,
     triangulation_patch_data_names(),
     std::vector<
@@ -3379,7 +3383,22 @@ GridOut::write_vtu(const Triangulation<dim, spacedim> &tria,
                  DataComponentInterpretation::DataComponentInterpretation>>(),
     vtu_flags,
     out);
+  if (vtu_flags.serialize_triangulation)
+    {
+      out << " </UnstructuredGrid>\n";
+      out << "<dealiiData  encoding=\"base64\">";
+      std::stringstream               outstring;
+      boost::archive::binary_oarchive ia(outstring);
+      tria.save(ia, 0);
+      const auto compressed = Utilities::compress(outstring.str());
+      out << Utilities::encode_base64({compressed.begin(), compressed.end()});
+      out << "\n</dealiiData>\n";
+      out << "</VTKFile>\n";
+    }
+  else
+    DataOutBase::write_vtu_footer(out);
 
+  out << std::flush;
   AssertThrow(out, ExcIO());
 }
 
diff --git a/tests/grid/grid_in_out_vtu_01.cc b/tests/grid/grid_in_out_vtu_01.cc
new file mode 100644 (file)
index 0000000..483829d
--- /dev/null
@@ -0,0 +1,81 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2020 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.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+
+// read and write a 2d file in the VTU format
+
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/grid_in.h>
+#include <deal.II/grid/grid_out.h>
+#include <deal.II/grid/tria.h>
+
+#include <boost/archive/text_oarchive.hpp>
+
+#include <string>
+
+#include "../tests.h"
+
+int
+main()
+{
+  initlog();
+
+  Triangulation<2> tria;
+  Triangulation<2> tria2;
+  GridGenerator::hyper_ball(tria);
+
+  for (const auto cell : tria.active_cell_iterators())
+    {
+      if (cell->center()[0] < 0)
+        cell->set_refine_flag();
+    }
+  tria.execute_coarsening_and_refinement();
+
+  for (const auto cell : tria.active_cell_iterators())
+    {
+      if (cell->center()[1] < 0)
+        cell->set_refine_flag();
+    }
+  tria.execute_coarsening_and_refinement();
+
+  {
+    std::ofstream     out("grid.vtu");
+    GridOut           go;
+    GridOutFlags::Vtu vtu_flags(true);
+    go.set_flags(vtu_flags);
+    go.write_vtu(tria, out);
+  }
+  {
+    std::ifstream in("grid.vtu");
+    GridIn<2>     gi;
+    gi.attach_triangulation(tria2);
+    gi.read_vtu(in);
+  }
+  std::stringstream stream1;
+  std::stringstream stream2;
+  {
+    boost::archive::text_oarchive oa(stream1);
+    tria.save(oa, 0);
+  }
+  {
+    boost::archive::text_oarchive oa(stream2);
+    tria2.save(oa, 0);
+  }
+  if (stream1.str() == stream2.str())
+    deallog << "OK" << std::endl;
+  else
+    deallog << "Not OK" << std::endl;
+  std::remove("grid.vtu");
+}
\ No newline at end of file
diff --git a/tests/grid/grid_in_out_vtu_01.output b/tests/grid/grid_in_out_vtu_01.output
new file mode 100644 (file)
index 0000000..0fd8fc1
--- /dev/null
@@ -0,0 +1,2 @@
+
+DEAL::OK

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.