]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add tests for serializing serial triangulation
authorPasquale Africa <pafrica@sissa.it>
Tue, 4 Jul 2023 19:41:32 +0000 (19:41 +0000)
committerPasquale Africa <pafrica@sissa.it>
Fri, 14 Jul 2023 14:50:40 +0000 (14:50 +0000)
tests/grid/save_load_01.cc [new file with mode: 0644]
tests/grid/save_load_01.output [new file with mode: 0644]
tests/grid/serialization_01.cc [new file with mode: 0644]
tests/grid/serialization_01.output [new file with mode: 0644]

diff --git a/tests/grid/save_load_01.cc b/tests/grid/save_load_01.cc
new file mode 100644 (file)
index 0000000..26f83af
--- /dev/null
@@ -0,0 +1,85 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2008 - 2023 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.
+//
+// ---------------------------------------------------------------------
+
+
+
+// Test Triangulation::load()/save().
+// Create a triangulation, save it and load it.
+// The initial and the loaded triangulations must be the same.
+
+#include <deal.II/dofs/dof_handler.h>
+#include <deal.II/dofs/dof_tools.h>
+
+#include <deal.II/fe/fe_q.h>
+
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/tria_description.h>
+
+#include <deal.II/lac/vector.h>
+
+#include "./tests.h"
+
+using namespace dealii;
+
+template <int dim, typename TriangulationType>
+void
+test(TriangulationType &triangulation)
+{
+  const std::string filename = "save_load_" + std::to_string(dim) + "d_out";
+
+  triangulation.save(filename);
+  deallog.push("save");
+  print_statistics(triangulation, false);
+  deallog.pop();
+
+  triangulation.clear();
+
+  triangulation.load(filename);
+  deallog.push("load");
+  print_statistics(triangulation, false);
+  deallog.pop();
+}
+
+
+int
+main(int argc, char **argv)
+{
+  initlog();
+
+  deallog.push("2d");
+  {
+    constexpr int dim = 2;
+
+    Triangulation<dim> triangulation;
+    GridGenerator::hyper_cube(triangulation);
+    triangulation.refine_global(3);
+
+    test<dim>(triangulation);
+  }
+  deallog.pop();
+
+  deallog.push("3d");
+  {
+    constexpr int dim = 3;
+
+    Triangulation<dim> triangulation;
+    GridGenerator::hyper_cube(triangulation);
+    triangulation.refine_global(3);
+
+    test<dim>(triangulation);
+  }
+  deallog.pop();
+}
diff --git a/tests/grid/save_load_01.output b/tests/grid/save_load_01.output
new file mode 100644 (file)
index 0000000..b3c4579
--- /dev/null
@@ -0,0 +1,21 @@
+
+DEAL:2d:save::n_levels:                  4
+DEAL:2d:save::n_cells:                   85
+DEAL:2d:save::n_active_cells:            64
+DEAL:2d:save::has_hanging_nodes:         false
+DEAL:2d:save::
+DEAL:2d:load::n_levels:                  4
+DEAL:2d:load::n_cells:                   85
+DEAL:2d:load::n_active_cells:            64
+DEAL:2d:load::has_hanging_nodes:         false
+DEAL:2d:load::
+DEAL:3d:save::n_levels:                  4
+DEAL:3d:save::n_cells:                   585
+DEAL:3d:save::n_active_cells:            512
+DEAL:3d:save::has_hanging_nodes:         false
+DEAL:3d:save::
+DEAL:3d:load::n_levels:                  4
+DEAL:3d:load::n_cells:                   585
+DEAL:3d:load::n_active_cells:            512
+DEAL:3d:load::has_hanging_nodes:         false
+DEAL:3d:load::
diff --git a/tests/grid/serialization_01.cc b/tests/grid/serialization_01.cc
new file mode 100644 (file)
index 0000000..e1700ef
--- /dev/null
@@ -0,0 +1,122 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2008 - 2023 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.
+//
+// ---------------------------------------------------------------------
+
+
+
+// Test serialization with triangulations.
+
+#include <deal.II/dofs/dof_handler.h>
+
+#include <deal.II/fe/fe_q.h>
+
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/tria.h>
+
+#include <deal.II/lac/la_vector.h>
+
+#include <deal.II/numerics/vector_tools.h>
+
+#include <boost/archive/text_iarchive.hpp>
+#include <boost/archive/text_oarchive.hpp>
+
+#include "./tests.h"
+
+using namespace dealii;
+
+template <int dim>
+class InterpolationFunction : public Function<dim>
+{
+public:
+  InterpolationFunction()
+    : Function<dim>(1)
+  {}
+
+  virtual double
+  value(const Point<dim> &p, const unsigned int component = 0) const
+  {
+    return p.norm();
+  }
+};
+
+template <int dim, typename TriangulationType>
+void
+test(TriangulationType &triangulation)
+{
+  DoFHandler<dim> dof_handler(triangulation);
+  dof_handler.distribute_dofs(FE_Q<dim>(2));
+
+  using VectorType = Vector<double>;
+
+  VectorType vector(dof_handler.n_dofs());
+
+  VectorTools::interpolate(dof_handler, InterpolationFunction<dim>(), vector);
+
+  VectorType vector_loaded(dof_handler.n_dofs());
+
+  print_statistics(triangulation, false);
+
+  std::stringstream stream;
+  {
+    boost::archive::text_oarchive oa(stream);
+    oa << triangulation;
+    oa << vector;
+  }
+
+  triangulation.clear();
+
+  {
+    boost::archive::text_iarchive ia(stream);
+    ia >> triangulation;
+    ia >> vector_loaded;
+  }
+  print_statistics(triangulation, false);
+
+  // Verify that error is 0.
+  VectorType error(vector);
+  error.add(-1, vector_loaded);
+
+  deallog << (error.linfty_norm() < 1e-16 ? "PASSED" : "FAILED") << std::endl;
+}
+
+
+int
+main(int argc, char **argv)
+{
+  initlog();
+
+  deallog.push("2d");
+  {
+    constexpr int dim = 2;
+
+    Triangulation<dim> triangulation;
+    GridGenerator::hyper_cube(triangulation);
+    triangulation.refine_global(3);
+
+    test<dim>(triangulation);
+  }
+  deallog.pop();
+
+  deallog.push("3d");
+  {
+    constexpr int dim = 3;
+
+    Triangulation<dim> triangulation;
+    GridGenerator::hyper_cube(triangulation);
+    triangulation.refine_global(3);
+
+    test<dim>(triangulation);
+  }
+  deallog.pop();
+}
diff --git a/tests/grid/serialization_01.output b/tests/grid/serialization_01.output
new file mode 100644 (file)
index 0000000..117a575
--- /dev/null
@@ -0,0 +1,23 @@
+
+DEAL:2d::n_levels:                  4
+DEAL:2d::n_cells:                   85
+DEAL:2d::n_active_cells:            64
+DEAL:2d::has_hanging_nodes:         false
+DEAL:2d::
+DEAL:2d::n_levels:                  4
+DEAL:2d::n_cells:                   85
+DEAL:2d::n_active_cells:            64
+DEAL:2d::has_hanging_nodes:         false
+DEAL:2d::
+DEAL:2d::PASSED
+DEAL:3d::n_levels:                  4
+DEAL:3d::n_cells:                   585
+DEAL:3d::n_active_cells:            512
+DEAL:3d::has_hanging_nodes:         false
+DEAL:3d::
+DEAL:3d::n_levels:                  4
+DEAL:3d::n_cells:                   585
+DEAL:3d::n_active_cells:            512
+DEAL:3d::has_hanging_nodes:         false
+DEAL:3d::
+DEAL:3d::PASSED

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.