]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Added possibility to attach data to bounding boxes. 10282/head
authorLuca Heltai <luca.heltai@sissa.it>
Wed, 27 May 2020 12:57:40 +0000 (14:57 +0200)
committerLuca Heltai <luca.heltai@sissa.it>
Wed, 27 May 2020 12:57:40 +0000 (14:57 +0200)
include/deal.II/base/bounding_box_data_out.h
tests/base/bounding_box_data_out_03.cc [new file with mode: 0644]
tests/base/bounding_box_data_out_03.output [new file with mode: 0644]

index f2c8ffaeaa52f5e9203ba25622f66e62d914b8ec..a0299bb92ab7938b78d59231f395ec944f49aa32 100644 (file)
@@ -45,13 +45,22 @@ template <int dim>
 class BoundingBoxDataOut : public DataOutInterface<dim, dim>
 {
 public:
+  /**
+   * Constructor.
+   */
   BoundingBoxDataOut() = default;
 
+  /**
+   * Destructor.
+   */
   ~BoundingBoxDataOut() = default;
 
   /**
    * Generate patches from a range of objects that can be converted by boost to
    * a collection of BoundingBox objects.
+   *
+   * You could pass to this function iterators to BoundingBox objects, or to
+   * pairs or tuples in which the first element is a BoundingBox object.
    */
   template <class ConvertibleToBoundingBoxIterator>
   void
@@ -61,11 +70,29 @@ public:
   /**
    * Generate patches from a container of objects that can be converted by boost
    * to a collection of BoundingBox objects.
+   *
+   * You could pass to this function a container of BoundingBox objects, or a
+   * container of pairs or tuples in which the first element is a BoundingBox
+   * object.
    */
   template <class Container>
   void
   build_patches(const Container &boxes);
 
+  /**
+   * Attach data to each output patch that was generated by build_patches().
+   *
+   * The @p datasets parameter is expected to have the same size of the
+   * container you used in the call to build_patches(), and each entry should
+   * have the same size of the @p dataset_names argument.
+   *
+   * @param[in] datasets The actual data to attach to each patch
+   * @param[in] dataset_names The name of each component of the dataset
+   */
+  void
+  add_datasets(const std::vector<std::vector<double>> &datasets,
+               const std::vector<std::string> &        dataset_names);
+
 protected:
   // Copy doc
   virtual const std::vector<dealii::DataOutBase::Patch<dim, dim>> &
@@ -116,9 +143,10 @@ BoundingBoxDataOut<dim>::build_patches(
       boost::geometry::convert(getter(*value), box);
       for (unsigned int v = 0; v < GeometryInfo<dim>::vertices_per_cell; ++v)
         {
-          patches[i].vertices[v]    = box.vertex(v);
-          patches[i].patch_index    = i;
-          patches[i].n_subdivisions = 1;
+          patches[i].vertices[v]          = box.vertex(v);
+          patches[i].patch_index          = i;
+          patches[i].n_subdivisions       = 1;
+          patches[i].points_are_available = false;
         }
       ++i;
     }
@@ -136,6 +164,27 @@ BoundingBoxDataOut<dim>::build_patches(const Container &boxes)
 
 
 
+template <int dim>
+void
+BoundingBoxDataOut<dim>::add_datasets(
+  const std::vector<std::vector<double>> &datasets,
+  const std::vector<std::string> &        names)
+{
+  AssertDimension(datasets.size(), patches.size());
+  dataset_names = names;
+  for (unsigned int i = 0; i < datasets.size(); ++i)
+    {
+      AssertDimension(datasets[i].size(), names.size());
+      patches[i].data.reinit(names.size(),
+                             GeometryInfo<dim>::vertices_per_cell);
+      for (unsigned int j = 0; j < names.size(); ++j)
+        for (unsigned int k = 0; k < GeometryInfo<dim>::vertices_per_cell; ++k)
+          patches[i].data(j, k) = datasets[i][j];
+    }
+}
+
+
+
 template <int dim>
 const std::vector<DataOutBase::Patch<dim, dim>> &
 BoundingBoxDataOut<dim>::get_patches() const
diff --git a/tests/base/bounding_box_data_out_03.cc b/tests/base/bounding_box_data_out_03.cc
new file mode 100644 (file)
index 0000000..7fcbc3c
--- /dev/null
@@ -0,0 +1,86 @@
+// ---------------------------------------------------------------------
+//
+// 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.
+//
+// ---------------------------------------------------------------------
+
+// Test BoundingBoxDataOut with an rtree of std::pair<BoundingBox<dim>, int>,
+// and output also as data the second entry of the pair.
+
+#include <deal.II/base/bounding_box_data_out.h>
+#include <deal.II/base/geometry_info.h>
+
+#include <deal.II/numerics/rtree.h>
+
+#include "../tests.h"
+
+using namespace dealii;
+
+template <int dim>
+void
+test()
+{
+  const unsigned int N = 3;
+  // Use a non-trivial object, convertible to BoundingBox by boost:
+  std::vector<std::pair<BoundingBox<dim>, unsigned int>> boxes(N);
+  auto unit = create_unit_bounding_box<dim>();
+
+  Point<dim> ones;
+  for (unsigned int i = 0; i < dim; ++i)
+    ones[i] = 1;
+
+  unsigned int i = 0;
+  for (auto &box : boxes)
+    {
+      const auto c = random_point<dim>();
+      const auto d = random_value();
+      box          = std::make_pair(BoundingBox<dim>({Point<dim>(c - d * ones),
+                                             Point<dim>(c + d * ones)}),
+                           i++);
+    }
+
+  const auto tree = pack_rtree(boxes);
+
+  std::vector<std::vector<double>> datasets(boxes.size(),
+                                            std::vector<double>(1));
+
+  {
+    unsigned i = 0;
+    for (auto p : tree)
+      datasets[i++][0] = p.second;
+  }
+
+  std::string fname = "boxes_" + std::to_string(dim) + ".vtk";
+  {
+    std::ofstream            ofile(fname);
+    BoundingBoxDataOut<dim>  data_out;
+    DataOutBase::VtkFlags    flags;
+    std::vector<std::string> names = {"id"};
+    flags.print_date_and_time      = false;
+    data_out.set_flags(flags);
+    data_out.build_patches(tree);
+    data_out.add_datasets(datasets, names);
+    data_out.write_vtk(ofile);
+  }
+  cat_file(fname.c_str());
+}
+
+
+int
+main()
+{
+  initlog();
+
+  test<1>();
+  test<2>();
+  test<3>();
+}
diff --git a/tests/base/bounding_box_data_out_03.output b/tests/base/bounding_box_data_out_03.output
new file mode 100644 (file)
index 0000000..146cf50
--- /dev/null
@@ -0,0 +1,100 @@
+
+# vtk DataFile Version 3.0
+#This file was generated by the deal.II library.
+ASCII
+DATASET UNSTRUCTURED_GRID
+
+POINTS 6 double
+0.445805 0 0
+1.23457 0 0
+-0.0153408 0 0
+1.58154 0 0
+0.714096 0 0
+1.1092 0 0
+
+CELLS 3 9
+2      0       1
+2      2       3
+2      4       5
+
+CELL_TYPES 3
+ 3 3 3
+POINT_DATA 6
+SCALARS id double 1
+LOOKUP_TABLE default
+0 0 1 1 2 2 
+
+# vtk DataFile Version 3.0
+#This file was generated by the deal.II library.
+ASCII
+DATASET UNSTRUCTURED_GRID
+
+POINTS 12 double
+0.057448 0.490455 0
+0.612997 0.490455 0
+0.057448 1.046 0
+0.612997 1.046 0
+-0.074901 -0.151474 0
+1.18284 -0.151474 0
+-0.074901 1.10627 0
+1.18284 1.10627 0
+-0.587445 -0.438829 0
+1.31701 -0.438829 0
+-0.587445 1.46563 0
+1.31701 1.46563 0
+
+CELLS 3 15
+4      0       1       3       2
+4      4       5       7       6
+4      8       9       11      10
+
+CELL_TYPES 3
+ 9 9 9
+POINT_DATA 12
+SCALARS id double 1
+LOOKUP_TABLE default
+0 0 0 0 1 1 1 1 2 2 2 2 
+
+# vtk DataFile Version 3.0
+#This file was generated by the deal.II library.
+ASCII
+DATASET UNSTRUCTURED_GRID
+
+POINTS 24 double
+0.774593 0.494109 0.575694
+1.0578 0.494109 0.575694
+0.774593 0.777314 0.575694
+1.0578 0.777314 0.575694
+0.774593 0.494109 0.858899
+1.0578 0.494109 0.858899
+0.774593 0.777314 0.858899
+1.0578 0.777314 0.858899
+0.469737 -0.120931 0.105655
+0.7442 -0.120931 0.105655
+0.469737 0.153532 0.105655
+0.7442 0.153532 0.105655
+0.469737 -0.120931 0.380118
+0.7442 -0.120931 0.380118
+0.469737 0.153532 0.380118
+0.7442 0.153532 0.380118
+0.674386 0.0268886 0.271154
+0.933967 0.0268886 0.271154
+0.674386 0.28647 0.271154
+0.933967 0.28647 0.271154
+0.674386 0.0268886 0.530735
+0.933967 0.0268886 0.530735
+0.674386 0.28647 0.530735
+0.933967 0.28647 0.530735
+
+CELLS 3 27
+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
+
+CELL_TYPES 3
+ 12 12 12
+POINT_DATA 24
+SCALARS id double 1
+LOOKUP_TABLE default
+0 0 0 0 0 0 0 0 1 1 1 1 1 1 1 1 2 2 2 2 2 2 2 2 
+

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.