]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Added method get_neighbor_type and EnumClass NeighborType in BoundingBox class. Added... 5117/head
authorGiovanni Alzetta <giovannialzetta@hotmail.it>
Thu, 28 Sep 2017 19:46:02 +0000 (19:46 +0000)
committerGiovanni Alzetta <giovannialzetta@hotmail.it>
Mon, 13 Nov 2017 14:14:22 +0000 (14:14 +0000)
13 files changed:
doc/doxygen/images/bounding_box_predicate.png [new file with mode: 0644]
doc/news/changes/minor/20171020GiovanniAlzetta [new file with mode: 0644]
include/deal.II/base/bounding_box.h
include/deal.II/grid/grid_tools.h
source/base/bounding_box.cc
source/grid/grid_tools.cc
source/grid/grid_tools.inst.in
tests/base/bounding_box_3.cc [new file with mode: 0644]
tests/base/bounding_box_3.output [new file with mode: 0644]
tests/grid/grid_tools_compute_mesh_predicate_bounding_box_1.cc [new file with mode: 0644]
tests/grid/grid_tools_compute_mesh_predicate_bounding_box_1.with_p4est=true.mpirun=2.output [new file with mode: 0644]
tests/grid/grid_tools_compute_mesh_predicate_bounding_box_1.with_p4est=true.mpirun=3.output [new file with mode: 0644]
tests/grid/grid_tools_compute_mesh_predicate_bounding_box_1.with_p4est=true.mpirun=4.output [new file with mode: 0644]

diff --git a/doc/doxygen/images/bounding_box_predicate.png b/doc/doxygen/images/bounding_box_predicate.png
new file mode 100644 (file)
index 0000000..09172d9
Binary files /dev/null and b/doc/doxygen/images/bounding_box_predicate.png differ
diff --git a/doc/news/changes/minor/20171020GiovanniAlzetta b/doc/news/changes/minor/20171020GiovanniAlzetta
new file mode 100644 (file)
index 0000000..675563b
--- /dev/null
@@ -0,0 +1,3 @@
+New: Added BoundingBox::get_neighbor_type() which returns the neighboring relation between two bounding boxes. Added the function GridTools::compute_mesh_predicate_bounding_box() which returns a vector of BoundingBoxes covering the cells for which a given predicate is true.
+<br>
+(Giovanni Alzetta, 2017/10/20)
index 075e53f5a8dc07ebf0f8ec2eafc1fea66738efe7..302e9c5fa578b5a13f7501ef108cc024fb35774b 100644 (file)
 
 DEAL_II_NAMESPACE_OPEN
 
+/**
+ * The enumerator NeighborType describes the neighboring relation between
+ * two bounding boxes.
+ */
+enum class NeighborType
+{
+  /**
+   * not neighbours: the intersection is empty
+   */
+  not_neighbors = 0,
+
+  /**
+   * simple neighbors: the boxes intersect with an intersection of dimension at most spacedim - 2
+   */
+  simple_neighbors = 1,
+
+  /**
+   * attached neighbors: neighbors with an intersection of dimension > spacedim - 2
+   */
+  attached_neighbors = 2,
+
+  /**
+   * mergeable neighbors: neighbors which can be expressed with a single Bounding Box, e.g.
+  *  @code
+  *  .--V--W    .-----V
+  *  |  |  | =  |     |
+  *  V--W--.    V-----.
+  *  @endcode
+  *  or one is inside the other
+   */
+  mergeable_neighbors = 3
+};
+
 /**
  * A class that represents a bounding box in a space with arbitrary dimension
  * <tt>spacedim</tt>.
@@ -76,6 +109,14 @@ public:
    */
   const std::pair<Point<spacedim,Number>,Point<spacedim,Number>> &get_boundary_points () const;
 
+  /**
+   * Check if the current object and @p other_bbox are neighbors, i.e. if the boxes
+   * have dimension spacedim, check if their intersection is non empty.
+   *
+   * Return an enumerator of type NeighborType.
+   */
+  NeighborType get_neighbor_type (const BoundingBox<spacedim,Number> &other_bbox) const;
+
   /**
    * Enlarge the current object so that it contains @p other_bbox .
    * If the current object already contains @p other_bbox then it is not changed
@@ -84,7 +125,7 @@ public:
   void merge_with(const BoundingBox<spacedim,Number> &other_bbox);
 
   /**
-   * Return true if the point is inside the Bounding Box, false otherwise
+   * Return true if the point is inside the Bounding Box, false otherwise.
    */
   bool point_inside (const Point<spacedim, Number> &p) const;
 
index f5d842f74b5f166386c35573131b7feddbec602e..1738c8051681884424eff147dd700682be7ff789 100644 (file)
@@ -17,6 +17,7 @@
 #define dealii_grid_tools_h
 
 
+#include <deal.II/base/bounding_box.h>
 #include <deal.II/base/config.h>
 #include <deal.II/base/bounding_box.h>
 #include <deal.II/base/geometry_info.h>
@@ -1199,7 +1200,6 @@ namespace GridTools
   compute_ghost_cell_layer_within_distance ( const MeshType &mesh,
                                              const double layer_thickness);
 
-
   /**
    * Compute and return a bounding box, defined through a pair of points
    * bottom left and top right, that surrounds a subdomain of the @p mesh.
@@ -1221,6 +1221,55 @@ namespace GridTools
   ( const MeshType                                                                    &mesh,
     const std::function<bool (const typename MeshType::active_cell_iterator &)> &predicate );
 
+  /**
+   * Compute a collection of bounding boxes so that all active cells for which the given predicate is true, are
+   * completely enclosed in at least one of the bounding boxes. Notice the cover is only guaranteed to contain
+   * all these active cells but it's not necessarily exact i.e. it can include a bigger area than their union.
+   *
+   * For each cell at a given refinement level containing active cells for which @p predicate is true,
+   * the function creates a bounding box of its children for which @p predicate is true.
+   *
+   * This results in a cover of all active cells for which @p predicate is true; the parameters
+   * @p allow_merge and @p max_boxes are used to reduce the number of cells at a computational cost and
+   * covering a bigger n-dimensional volume.
+   *
+   * The parameters to control the algorithm are:
+   * - @p predicate : the property of the cells to enclose e.g. IteratorFilters::LocallyOwnedCell .
+   *  The predicate is tested only on active cells.
+   * - @p refinement_level : it defines the level at which the initial bounding box are created. The refinement
+   *  should be set to a coarse refinement level. A bounding box is created for each active cell at coarser
+   *  level than @p refinement_level; if @p refinement_level is higher than the number of levels of the
+   *  triangulation an exception is thrown.
+   * - @p allow_merge : This flag allows for box merging and, by default, is false. The algorithm has a cost of
+   *  O(N^2) where N is the number of the bounding boxes created from the refinement level; for this reason, if
+   *  the flag is set to true, make sure to choose wisely a coarse enough @p refinement_level.
+   * - @p max_boxes : the maximum number of bounding boxes to compute. If more are created the smaller ones are
+   *  merged with neighbors. By default after merging the boxes which can be expressed as a single one no
+   *  more boxes are merged. See the BoundingBox::get_neighbor_type () function for details.
+   *  Notice only neighboring cells are merged (see the @p get_neighbor_type  function in bounding box class): if
+   *  the target number of bounding boxes max_boxes can't be reached by merging neighbors an exception is thrown
+   *
+   * The following image describes an example of the algorithm with @p refinement_level = 2, @p allow_merge = true
+   * and @p max_boxes = 1. The cells with the property predicate are in red, the area of a bounding box is
+   * slightly orange.
+   * @image html bounding_box_predicate.png
+   * - 1. In black we can see the cells of the current level.
+   * - 2. For each cell containing the red area a bounding box is created: by default these are returned.
+   * - 3. Because @p allow_merge = true the number of bounding boxes is reduced while not changing the cover.
+   *  If @p max_boxes was left as default or bigger than 1 these two boxes would be returned.
+   * - 4. Because @p max_boxes = 1 the smallest bounding box is merged to the bigger one.
+   * Notice it is important to choose the parameters wisely. For instance, @p allow_merge = false and
+   * @p refinement_level = 1 returns the very same bounding box but with a fraction of the computational cost.
+   *
+   * This function does not take into account the curvature of cells and thus it is not suited for handling
+   * curved geometry: the mapping is assumed to be linear.
+   */
+  template < class MeshType >
+  std::vector< BoundingBox< MeshType::space_dimension > >
+  compute_mesh_predicate_bounding_box
+  ( const MeshType &mesh,
+    const std::function<bool (const typename MeshType::active_cell_iterator &)> &predicate,
+    const unsigned int &refinement_level = 0, const bool &allow_merge = false, const unsigned int &max_boxes = numbers::invalid_unsigned_int);
 
   /**
    * Return the adjacent cells of all the vertices. If a vertex is also a
index 9390d4b8ccb1db7a396822fa397e66606026bdcc..42843b725e2d5bb66eb0e3cdab37fa1c9141923b 100644 (file)
@@ -12,7 +12,6 @@
 // the top level of the deal.II distribution.
 //
 // ---------------------------------------------------------------------
-
 #include <deal.II/base/bounding_box.h>
 
 DEAL_II_NAMESPACE_OPEN
@@ -25,7 +24,10 @@ bool BoundingBox<spacedim,Number>::point_inside (const Point<spacedim, Number> &
       //Bottom left-top right convention: the point is outside if it's smaller than the
       //first or bigger than the second boundary point
       //The bounding box is defined as a closed set
-      if ( p[i] < this->boundary_points.first[i] || this->boundary_points.second[i] < p[i])
+      if ( std::numeric_limits<Number>::epsilon()*(std::abs( this->boundary_points.first[i] + p[i]) )
+           < this->boundary_points.first[i] - p[i]
+           || std::numeric_limits<Number>::epsilon()*(std::abs( this->boundary_points.second[i] + p[i]) ) <
+           p[i] - this->boundary_points.second[i] )
         return false;
     }
   return true;
@@ -41,6 +43,86 @@ void BoundingBox<spacedim,Number>::merge_with(const BoundingBox<spacedim,Number>
     }
 }
 
+template <int spacedim, typename Number>
+NeighborType BoundingBox < spacedim,Number >::get_neighbor_type (const BoundingBox<spacedim,Number> &other_bbox) const
+{
+  if (spacedim == 1)
+    {
+      // In dimension 1 if the two bounding box are neighbors
+      // we can merge them
+      if ( this->point_inside(other_bbox.boundary_points.first)
+           || this->point_inside(other_bbox.boundary_points.second) )
+        return NeighborType::mergeable_neighbors;
+      return NeighborType::not_neighbors;
+    }
+  else
+    {
+      std::vector< Point<spacedim,Number> > bbox1;
+      bbox1.push_back(this->get_boundary_points().first);
+      bbox1.push_back(this->get_boundary_points().second);
+      std::vector< Point<spacedim,Number> > bbox2;
+      bbox2.push_back(other_bbox.get_boundary_points().first);
+      bbox2.push_back(other_bbox.get_boundary_points().second);
+
+      // Step 1: testing if the boxes are close enough to intersect
+      for ( unsigned int d=0 ; d<spacedim; ++d )
+        if ( bbox1[0][d] * ( 1 - std::numeric_limits<Number>::epsilon() ) > bbox2[1][d]
+             || bbox2[0][d] * (1 - std::numeric_limits<Number>::epsilon() ) > bbox1[1][d] )
+          return NeighborType::not_neighbors;
+
+      // The boxes intersect: we need to understand now how they intersect.
+      // We begin by computing the intersection:
+      std::vector<double> intersect_bbox_min;
+      std::vector<double> intersect_bbox_max;
+      for (unsigned int d=0; d < spacedim ; ++d)
+        {
+          intersect_bbox_min.push_back(std::max(bbox1[0][d],bbox2[0][d]));
+          intersect_bbox_max.push_back(std::min(bbox1[1][d],bbox2[1][d]));
+        }
+
+      // Finding the intersection's dimension
+
+      unsigned int intersect_dim = spacedim;
+      for (unsigned int d=0; d < spacedim; ++d)
+        if ( std::abs( intersect_bbox_min[d] - intersect_bbox_max[d] ) <=
+             std::numeric_limits<Number>::epsilon() * (std::abs( intersect_bbox_min[d]) + std::abs(intersect_bbox_max[d] )) )
+          --intersect_dim;
+
+      if ( intersect_dim == 0 || intersect_dim == spacedim - 2 )
+        return NeighborType::simple_neighbors;
+
+      // Checking the two mergeable cases: first if the boxes are aligned so that they can be merged
+      unsigned int not_align_1 = 0, not_align_2 = 0;
+      bool same_direction = true;
+      for (unsigned int d = 0; d < spacedim; ++d)
+        {
+          if (std::abs(bbox2[0][d] - bbox1[0][d]) >
+              std::numeric_limits<double>::epsilon() * (std::abs( bbox2[0][d]) + std::abs(bbox1[0][d] )) )
+            ++not_align_1;
+          if (std::abs(bbox1[1][d] - bbox2[1][d]) >
+              std::numeric_limits<double>::epsilon() * (std::abs( bbox1[1][d]) + std::abs( bbox2[1][d] )) )
+            ++not_align_2;
+          if (not_align_1 != not_align_2)
+            {
+              same_direction = false;
+              break;
+            }
+        }
+
+      if ( not_align_1 <= 1 && not_align_2 <= 1 && same_direction)
+        return NeighborType::mergeable_neighbors;
+
+      // Second: one box is contained/equal to the other
+      if ( (this->point_inside(bbox2[0]) && this->point_inside(bbox2[1]) )
+           || (other_bbox.point_inside(bbox1[0]) && other_bbox.point_inside(bbox1[1]) ) )
+        return NeighborType::mergeable_neighbors;
+
+      // Degenerate and mergeable cases have been found, it remains:
+      return NeighborType::attached_neighbors;
+    }
+
+}
+
 template <int spacedim, typename Number>
 const std::pair<Point<spacedim,Number>,Point<spacedim,Number>> &BoundingBox<spacedim,Number>::get_boundary_points () const
 {
index 6dd0b32123b5ccfb197e41b9b70ea8af0f4a693e..b80ab571428ffb3797b880564ef6bcddbf6a6e1b 100644 (file)
@@ -55,6 +55,7 @@
 #include <numeric>
 #include <list>
 #include <set>
+#include <tuple>
 
 
 DEAL_II_NAMESPACE_OPEN
@@ -2167,6 +2168,174 @@ next_cell:
 
 
 
+  namespace internal
+  {
+    namespace BoundingBoxPredicate
+    {
+      template < class MeshType >
+      std::tuple< BoundingBox < MeshType::space_dimension >, bool >
+      compute_cell_predicate_bounding_box
+      (const typename MeshType::cell_iterator &parent_cell,
+       const std::function<bool (const typename MeshType::active_cell_iterator &)> &predicate)
+      {
+        bool has_predicate = false; // Start assuming there's no cells with predicate inside
+        std::vector< typename MeshType::active_cell_iterator > active_cells;
+        if (parent_cell->active())
+          active_cells = {parent_cell};
+        else
+          //Finding all active cells descendants of the current one (or the current one if it is active)
+          active_cells = get_active_child_cells < MeshType > (parent_cell);
+
+        const unsigned int spacedim = MeshType::space_dimension;
+
+        // Looking for the first active cell which has the property predicate
+        unsigned int i = 0;
+        while ( i < active_cells.size() && !predicate(active_cells[i]) )
+          ++i;
+
+        // No active cells or no active cells with property
+        if ( active_cells.size() == 0 || i == active_cells.size() )
+          {
+            BoundingBox<spacedim> bbox;
+            return std::make_tuple(bbox, has_predicate);
+          }
+
+        // The two boundary points defining the boundary box
+        Point<spacedim> maxp = active_cells[i]->vertex(0);
+        Point<spacedim> minp = active_cells[i]->vertex(0);
+
+        for (; i < active_cells.size() ; ++i)
+          if ( predicate(active_cells[i]) )
+            for (unsigned int v=0; v<GeometryInfo<spacedim>::vertices_per_cell; ++v)
+              for ( unsigned int d=0; d<spacedim; ++d)
+                {
+                  minp[d] = std::min( minp[d], active_cells[i]->vertex(v)[d]);
+                  maxp[d] = std::max( maxp[d], active_cells[i]->vertex(v)[d]);
+                }
+
+        has_predicate = true;
+        BoundingBox < spacedim > bbox(std::make_pair(minp,maxp));
+        return std::make_tuple(bbox, has_predicate);
+      }
+    }
+  }
+
+
+
+  template < class MeshType >
+  std::vector< BoundingBox<MeshType::space_dimension> >
+  compute_mesh_predicate_bounding_box
+  (const MeshType &mesh,
+   const std::function<bool (const typename MeshType::active_cell_iterator &)> &predicate,
+   const unsigned int &refinement_level,  const bool &allow_merge, const unsigned int &max_boxes)
+  {
+    // Algorithm brief description: begin with creating bounding boxes of all cells at
+    // refinement_level (and coarser levels if there are active cells) which have the predicate
+    // property. These are then merged
+
+    Assert( refinement_level <= mesh.n_levels(),
+            ExcMessage ( "Error: refinement level is higher then total levels in the triangulation!") );
+
+    const unsigned int spacedim = MeshType::space_dimension;
+    std::vector< BoundingBox < spacedim > > bounding_boxes;
+
+    // Creating a bounding box for all active cell on coarser level
+
+    for (unsigned int i=0; i < refinement_level; ++i)
+      for (typename MeshType::cell_iterator cell: mesh.active_cell_iterators_on_level(i))
+        {
+          bool has_predicate = false;
+          BoundingBox < spacedim > bbox;
+          std::tie(bbox, has_predicate) =
+            internal::BoundingBoxPredicate::compute_cell_predicate_bounding_box <MeshType> (cell, predicate);
+          if (has_predicate)
+            bounding_boxes.push_back(bbox);
+        }
+
+    // Creating a Bounding Box for all cells on the chosen refinement_level
+    for (typename MeshType::cell_iterator cell: mesh.cell_iterators_on_level(refinement_level))
+      {
+        bool has_predicate = false;
+        BoundingBox < spacedim > bbox;
+        std::tie(bbox, has_predicate) =
+          internal::BoundingBoxPredicate::compute_cell_predicate_bounding_box <MeshType> (cell, predicate);
+        if (has_predicate)
+          bounding_boxes.push_back(bbox);
+      }
+
+    if ( !allow_merge)
+      // If merging is not requested return the created bounding_boxes
+      return bounding_boxes;
+    else
+      {
+        // Merging part of the algorithm
+        // Part 1: merging neighbors
+        // This array stores the indices of arrays we have already merged
+        std::vector<unsigned int> merged_boxes_idx;
+        bool found_neighbors = true;
+
+        // We merge only nighbors which can be expressed by a single bounding box
+        // e.g. in 1d [0,1] and [1,2] can be described with [0,2] without losing anything
+        while (found_neighbors)
+          {
+            found_neighbors = false;
+            for (unsigned int i=0; i<bounding_boxes.size()-1; ++i)
+              {
+                if ( std::find(merged_boxes_idx.begin(),merged_boxes_idx.end(),i) == merged_boxes_idx.end())
+                  for (unsigned int j=i+1; j<bounding_boxes.size(); ++j)
+                    if ( std::find(merged_boxes_idx.begin(),merged_boxes_idx.end(),j) == merged_boxes_idx.end()
+                         && bounding_boxes[i].get_neighbor_type (bounding_boxes[j]) == NeighborType::mergeable_neighbors  )
+                      {
+                        bounding_boxes[i].merge_with(bounding_boxes[j]);
+                        merged_boxes_idx.push_back(j);
+                        found_neighbors = true;
+                      }
+              }
+          }
+
+        // Copying the merged boxes into merged_b_boxes
+        std::vector< BoundingBox < spacedim > > merged_b_boxes;
+        for (unsigned int i=0; i<bounding_boxes.size(); ++i)
+          if (std::find(merged_boxes_idx.begin(),merged_boxes_idx.end(),i) == merged_boxes_idx.end())
+            merged_b_boxes.push_back(bounding_boxes[i]);
+
+        // Part 2: if there are too many bounding boxes, merging smaller boxes
+        // This has sense only in dimension 2 or greater, since  in dimension 1,
+        // neighboring intervals can always be merged without problems
+        if ( (merged_b_boxes.size() > max_boxes) && (spacedim > 1) )
+          {
+            std::vector<double> volumes;
+            for (unsigned int i=0; i< merged_b_boxes.size(); ++i)
+              volumes.push_back(merged_b_boxes[i].volume());
+
+            while ( merged_b_boxes.size() > max_boxes)
+              {
+                unsigned int min_idx = std::min_element(volumes.begin(),volumes.end()) -
+                                       volumes.begin();
+                volumes.erase(volumes.begin() + min_idx);
+                //Finding a neighbor
+                bool not_removed = true;
+                for (unsigned int i=0; i<merged_b_boxes.size() && not_removed; ++i)
+                  // We merge boxes if we have "attached" or "mergeable" neighbors, even though mergeable should
+                  // be dealt with in Part 1
+                  if ( i != min_idx &&
+                       (merged_b_boxes[i].get_neighbor_type (merged_b_boxes[min_idx]) == NeighborType::attached_neighbors ||
+                        merged_b_boxes[i].get_neighbor_type (merged_b_boxes[min_idx]) == NeighborType::mergeable_neighbors) )
+                    {
+                      merged_b_boxes[i].merge_with(merged_b_boxes[min_idx]);
+                      merged_b_boxes.erase(merged_b_boxes.begin() + min_idx);
+                      not_removed = false;
+                    }
+                Assert( !not_removed,
+                        ExcMessage ( "Error: couldn't merge bounding boxes!") );
+              }
+          }
+        Assert( merged_b_boxes.size() <= max_boxes,
+                ExcMessage ( "Error: couldn't reach target number of bounding boxes!") );
+        return merged_b_boxes;
+      }
+  }
+
   template <int dim, int spacedim>
   std::vector<std::set<typename Triangulation<dim,spacedim>::active_cell_iterator> >
   vertex_to_cell_map(const Triangulation<dim,spacedim> &triangulation)
index c4de24d47308639897b204c5733ffc1b0d7329f5..9d1246976c6cb0982a3c5e834c4fa67cee3c01c1 100644 (file)
@@ -27,6 +27,12 @@ for (X : TRIANGULATIONS; deal_II_dimension : DIMENSIONS ; deal_II_space_dimensio
                                    const std::vector<std::vector<Tensor<1,deal_II_space_dimension> > > &,
                                    const dealii::internal::ActiveCellIterator<deal_II_dimension, deal_II_space_dimension, X>::type &,
                                    const std::vector<bool> &);
+
+    template
+    std::vector < BoundingBox< deal_II_space_dimension > >
+    compute_mesh_predicate_bounding_box< X >
+    (const X &, const std::function<bool (const dealii::internal::ActiveCellIterator<deal_II_dimension, deal_II_space_dimension, X>::type&)> &,
+     const unsigned int &, const bool&, const unsigned int &);
     \}
 
 #endif
@@ -100,6 +106,7 @@ for (X : TRIANGULATION_AND_DOFHANDLERS; deal_II_dimension : DIMENSIONS ; deal_II
     compute_bounding_box (const X &,
                           const std::function<bool (const dealii::internal::ActiveCellIterator<deal_II_dimension, deal_II_space_dimension, X>::type&)> &);
 
+
     template
     std::list<std::pair<X::cell_iterator, X::cell_iterator> >
     get_finest_common_cells (const X &mesh_1,
diff --git a/tests/base/bounding_box_3.cc b/tests/base/bounding_box_3.cc
new file mode 100644 (file)
index 0000000..9baa617
--- /dev/null
@@ -0,0 +1,235 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2017 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.
+//
+// ---------------------------------------------------------------------
+
+
+// test for BoundingBox<unsigned int spacedim> which tests the function
+// get_neighbor_type
+
+#include "../tests.h"
+
+#include <deal.II/base/point.h>
+#include <deal.II/base/bounding_box.h>
+
+template <int spacedim>
+void test_bounding_box()
+{
+
+  std::pair<Point<spacedim>,Point<spacedim>> unit;
+  for (int i=0; i<spacedim; i++)
+    {
+      unit.first[i] = 0.0;
+      unit.second[i] = 1.0;
+    }
+
+  BoundingBox<spacedim> a(unit);
+
+  deallog << "Bounding box boundaries A: " << std::endl;
+  deallog << a.get_boundary_points().first << std::endl;
+  deallog << a.get_boundary_points().second << std::endl;
+  deallog << "Is neighbor of itself: " << (int) a.get_neighbor_type (a) << std::endl;
+
+  // This is a simple neighbor
+  std::pair<Point<spacedim>,Point<spacedim>> second;
+  for (int i=0; i<spacedim; i++)
+    {
+      second.first[i] = 1.0;
+      second.second[i] = 2.0;
+    }
+
+  BoundingBox<spacedim> b(second);
+
+  deallog << "Bounding box boundaries B: " << std::endl;
+  deallog << b.get_boundary_points().first << std::endl;
+  deallog << b.get_boundary_points().second << std::endl;
+
+  deallog << "Is neighbor A with B: " << (int) a.get_neighbor_type (b) << std::endl;
+  deallog << "Is neighbor B with A: " << (int) b.get_neighbor_type (a) << std::endl;
+  deallog << "Is neighbor B with B: " << (int) b.get_neighbor_type (b) << std::endl;
+
+  // Testing the two 3d only cases
+  if (spacedim == 3)
+    {
+      std::pair<Point<spacedim>,Point<spacedim>> second_1;
+      second_1.first[0] = 0.0;
+      second_1.second[0] = 1.0;
+      for (int i=1; i<spacedim; i++)
+        {
+          second_1.first[i] = 1.0;
+          second_1.second[i] = 2.0;
+        }
+      BoundingBox<spacedim> b_1(second_1);
+
+      deallog << "Bounding box boundaries B1: " << std::endl;
+      deallog << b_1.get_boundary_points().first << std::endl;
+      deallog << b_1.get_boundary_points().second << std::endl;
+
+      deallog << "Is neighbor A with B1: " << (int)  a.get_neighbor_type (b_1) << std::endl;
+      deallog << "Is neighbor B1 with A: " << (int) b_1.get_neighbor_type (a) << std::endl;
+      deallog << "Is neighbor B1 with B: " << (int) b_1.get_neighbor_type (b_1) << std::endl;
+
+      //Case of partial edge which is still simple neighbor
+      //Note: in dimension 1 everything is only either not neighbor or mergeable
+      //neighbor
+      second_1.first[0] = 0.2;
+
+      BoundingBox<spacedim> b_2(second_1);
+
+      deallog << "Bounding box boundaries B2: " << std::endl;
+      deallog << b_2.get_boundary_points().first << std::endl;
+      deallog << b_2.get_boundary_points().second << std::endl;
+
+      deallog << "Is neighbor A with B2: " << (int) a.get_neighbor_type (b_2) << std::endl;
+      deallog << "Is neighbor B2 with A: " << (int) b_2.get_neighbor_type (a) << std::endl;
+      deallog << "Is neighbor B2 with B: " << (int) b_2.get_neighbor_type (b_2) << std::endl;
+
+      //Case which is attached_neighbor
+
+      second_1.first[0] = 0.0;
+      second_1.first[2] = 0.8;
+
+      BoundingBox<spacedim> b_3(second_1);
+
+      deallog << "Bounding box boundaries B3: " << std::endl;
+      deallog << b_3.get_boundary_points().first << std::endl;
+      deallog << b_3.get_boundary_points().second << std::endl;
+
+      deallog << "Is neighbor A with B3: " << (int) a.get_neighbor_type (b_3) << std::endl;
+      deallog << "Is neighbor B3 with A: " << (int) b_3.get_neighbor_type (a) << std::endl;
+      deallog << "Is neighbor B3 with B: " << (int) b_3.get_neighbor_type (b_3) << std::endl;
+    }
+
+  // C contains both A and B
+  BoundingBox<spacedim> c(std::make_pair(unit.first,second.second));
+  deallog << "Bounding box boundaries C: " << std::endl;
+  deallog << c.get_boundary_points().first << std::endl;
+  deallog << c.get_boundary_points().second << std::endl;
+
+  deallog << "Is neighbor C with B: " << (int) c.get_neighbor_type (b) << std::endl;
+  deallog << "Is neighbor B with C: " << (int) b.get_neighbor_type (c) << std::endl;
+  deallog << "Is neighbor A with C: " << (int) a.get_neighbor_type (c) << std::endl;
+  deallog << "Is neighbor C with A: " << (int) c.get_neighbor_type (a) << std::endl;
+  deallog << "Is neighbor C with C: " << (int) c.get_neighbor_type (c) << std::endl;
+
+  // D contains A, intersects with B
+  unit.second *= 1.4;
+  BoundingBox<spacedim> d(unit);
+
+  deallog << "Bounding box boundaries D: " << std::endl;
+  deallog << d.get_boundary_points().first << std::endl;
+  deallog << d.get_boundary_points().second << std::endl;
+
+  deallog << "Is neighbor D with A: " << (int) d.get_neighbor_type (a) << std::endl;
+  deallog << "Vice-versa: " << (int) a.get_neighbor_type (d) << std::endl;
+  deallog << "Is neighbor D with B: " << (int) d.get_neighbor_type (b) << std::endl;
+  deallog << "Vice-versa: " << (int) b.get_neighbor_type (d) << std::endl;
+  deallog << "Is neighbor D with C: " << (int) d.get_neighbor_type (c) << std::endl;
+  deallog << "Vice-versa: " << (int) c.get_neighbor_type (d) << std::endl;
+  deallog << "Is neighbor D with D: " << (int) d.get_neighbor_type (d) << std::endl;
+
+  for (int i=0; i<spacedim; i++)
+    {
+      second.first[i] = -10.0;
+      second.second[i] = -8.0;
+    }
+
+  // E is separated from all others
+  BoundingBox<spacedim> e(second);
+
+  deallog << "Bounding box boundaries E: " << std::endl;
+  deallog << e.get_boundary_points().first << std::endl;
+  deallog << e.get_boundary_points().second << std::endl;
+
+  deallog << "Is neighbor E with A: " << (int) e.get_neighbor_type (a) << std::endl;
+  deallog << "Vice-versa: " << (int) a.get_neighbor_type (e) << std::endl;
+  deallog << "Is neighbor E with B: " << (int) e.get_neighbor_type (b) << std::endl;
+  deallog << "Vice-versa: " << (int) b.get_neighbor_type (e) << std::endl;
+  deallog << "Is neighbor E with C: " << (int) e.get_neighbor_type (c) << std::endl;
+  deallog << "Vice-versa: " << (int) c.get_neighbor_type (e) << std::endl;
+  deallog << "Is neighbor E with D: " << (int) e.get_neighbor_type (d) << std::endl;
+  deallog << "Vice-versa: " << (int) d.get_neighbor_type (e) << std::endl;
+  deallog << "Is neighbor E with E: " << (int) e.get_neighbor_type (e) << std::endl;
+
+  // A and F are mergeable because one next to the other
+  for (int i=0; i<spacedim; i++)
+    {
+      second.first[i] = 0.0;
+      second.second[i] = 1.0;
+    }
+  second.first[0] += 1.0;
+  second.second[0] += 1.0;
+
+  BoundingBox<spacedim> f(second);
+
+  deallog << "Bounding box boundaries F: " << std::endl;
+  deallog << f.get_boundary_points().first << std::endl;
+  deallog << f.get_boundary_points().second << std::endl;
+
+  deallog << "Is neighbor F with A: " << (int) f.get_neighbor_type (a) << std::endl;
+  deallog << "Vice-versa: " << (int) a.get_neighbor_type (f) << std::endl;
+  deallog << "Is neighbor F with B: " << (int) f.get_neighbor_type (b) << std::endl;
+  deallog << "Vice-versa: " << (int) b.get_neighbor_type (f) << std::endl;
+  deallog << "Is neighbor F with C: " << (int) f.get_neighbor_type (c) << std::endl;
+  deallog << "Vice-versa: " << (int) c.get_neighbor_type (f) << std::endl;
+  deallog << "Is neighbor F with D: " << (int) f.get_neighbor_type (d) << std::endl;
+  deallog << "Vice-versa: " << (int) d.get_neighbor_type (f) << std::endl;
+  deallog << "Is neighbor F with E: " << (int) f.get_neighbor_type (e) << std::endl;
+  deallog << "Vice-versa: " << (int) e.get_neighbor_type (f) << std::endl;
+  deallog << "Is neighbor F with F: " << (int) f.get_neighbor_type (f) << std::endl;
+
+  deallog << "End test for dimension " << spacedim << std::endl;
+  deallog << std::endl;
+
+  // G is mergeable with A and F and C..
+  second.first[0] -= 0.2;
+  second.second[0] -= 0.2;
+
+  BoundingBox<spacedim> g(second);
+
+  deallog << "Bounding box boundaries G: " << std::endl;
+  deallog << g.get_boundary_points().first << std::endl;
+  deallog << g.get_boundary_points().second << std::endl;
+
+  deallog << "Is neighbor G with A: " << (int) g.get_neighbor_type (a) << std::endl;
+  deallog << "Vice-versa: " << (int) a.get_neighbor_type (g) << std::endl;
+  deallog << "Is neighbor G with B: " << (int) g.get_neighbor_type (b) << std::endl;
+  deallog << "Vice-versa: " << (int) b.get_neighbor_type (g) << std::endl;
+  deallog << "Is neighbor G with C: " << (int) g.get_neighbor_type (c) << std::endl;
+  deallog << "Vice-versa: " << (int) c.get_neighbor_type (g) << std::endl;
+  deallog << "Is neighbor G with D: " << (int) g.get_neighbor_type (d) << std::endl;
+  deallog << "Vice-versa: " << (int) d.get_neighbor_type (g) << std::endl;
+  deallog << "Is neighbor G with E: " << (int) g.get_neighbor_type (e) << std::endl;
+  deallog << "Vice-versa: " << (int) e.get_neighbor_type (g) << std::endl;
+  deallog << "Is neighbor G with F: " << (int) g.get_neighbor_type (f) << std::endl;
+  deallog << "Vice-versa: " << (int) f.get_neighbor_type (g) << std::endl;
+  deallog << "Is neighbor G with G: " << (int) g.get_neighbor_type (g) << std::endl;
+
+  deallog << "End test for dimension " << spacedim << std::endl;
+  deallog << std::endl;
+}
+
+int main()
+{
+  initlog();
+
+  deallog << "Test: Bounding Box class Is neighbor and volume functions" << std::endl;
+  deallog << std::endl << "Test for dimension 1" << std::endl;
+  test_bounding_box<1>();
+
+  deallog << std::endl << "Test for dimension 2" << std::endl;
+  test_bounding_box<2>();
+
+  deallog << std::endl << "Test for dimension 3" << std::endl;
+  test_bounding_box<3>();
+}
diff --git a/tests/base/bounding_box_3.output b/tests/base/bounding_box_3.output
new file mode 100644 (file)
index 0000000..d12f6c4
--- /dev/null
@@ -0,0 +1,248 @@
+
+DEAL::Test: Bounding Box class Is neighbor and volume functions
+DEAL::
+DEAL::Test for dimension 1
+DEAL::Bounding box boundaries A: 
+DEAL::0.00000
+DEAL::1.00000
+DEAL::Is neighbor of itself: 3
+DEAL::Bounding box boundaries B: 
+DEAL::1.00000
+DEAL::2.00000
+DEAL::Is neighbor A with B: 3
+DEAL::Is neighbor B with A: 3
+DEAL::Is neighbor B with B: 3
+DEAL::Bounding box boundaries C: 
+DEAL::0.00000
+DEAL::2.00000
+DEAL::Is neighbor C with B: 3
+DEAL::Is neighbor B with C: 3
+DEAL::Is neighbor A with C: 3
+DEAL::Is neighbor C with A: 3
+DEAL::Is neighbor C with C: 3
+DEAL::Bounding box boundaries D: 
+DEAL::0.00000
+DEAL::1.40000
+DEAL::Is neighbor D with A: 3
+DEAL::Vice-versa: 3
+DEAL::Is neighbor D with B: 3
+DEAL::Vice-versa: 3
+DEAL::Is neighbor D with C: 3
+DEAL::Vice-versa: 3
+DEAL::Is neighbor D with D: 3
+DEAL::Bounding box boundaries E: 
+DEAL::-10.0000
+DEAL::-8.00000
+DEAL::Is neighbor E with A: 0
+DEAL::Vice-versa: 0
+DEAL::Is neighbor E with B: 0
+DEAL::Vice-versa: 0
+DEAL::Is neighbor E with C: 0
+DEAL::Vice-versa: 0
+DEAL::Is neighbor E with D: 0
+DEAL::Vice-versa: 0
+DEAL::Is neighbor E with E: 3
+DEAL::Bounding box boundaries F: 
+DEAL::1.00000
+DEAL::2.00000
+DEAL::Is neighbor F with A: 3
+DEAL::Vice-versa: 3
+DEAL::Is neighbor F with B: 3
+DEAL::Vice-versa: 3
+DEAL::Is neighbor F with C: 3
+DEAL::Vice-versa: 3
+DEAL::Is neighbor F with D: 3
+DEAL::Vice-versa: 3
+DEAL::Is neighbor F with E: 0
+DEAL::Vice-versa: 0
+DEAL::Is neighbor F with F: 3
+DEAL::End test for dimension 1
+DEAL::
+DEAL::Bounding box boundaries G: 
+DEAL::0.800000
+DEAL::1.80000
+DEAL::Is neighbor G with A: 3
+DEAL::Vice-versa: 3
+DEAL::Is neighbor G with B: 3
+DEAL::Vice-versa: 3
+DEAL::Is neighbor G with C: 0
+DEAL::Vice-versa: 3
+DEAL::Is neighbor G with D: 3
+DEAL::Vice-versa: 3
+DEAL::Is neighbor G with E: 0
+DEAL::Vice-versa: 0
+DEAL::Is neighbor G with F: 3
+DEAL::Vice-versa: 3
+DEAL::Is neighbor G with G: 3
+DEAL::End test for dimension 1
+DEAL::
+DEAL::
+DEAL::Test for dimension 2
+DEAL::Bounding box boundaries A: 
+DEAL::0.00000 0.00000
+DEAL::1.00000 1.00000
+DEAL::Is neighbor of itself: 3
+DEAL::Bounding box boundaries B: 
+DEAL::1.00000 1.00000
+DEAL::2.00000 2.00000
+DEAL::Is neighbor A with B: 1
+DEAL::Is neighbor B with A: 1
+DEAL::Is neighbor B with B: 3
+DEAL::Bounding box boundaries C: 
+DEAL::0.00000 0.00000
+DEAL::2.00000 2.00000
+DEAL::Is neighbor C with B: 3
+DEAL::Is neighbor B with C: 3
+DEAL::Is neighbor A with C: 3
+DEAL::Is neighbor C with A: 3
+DEAL::Is neighbor C with C: 3
+DEAL::Bounding box boundaries D: 
+DEAL::0.00000 0.00000
+DEAL::1.40000 1.40000
+DEAL::Is neighbor D with A: 3
+DEAL::Vice-versa: 3
+DEAL::Is neighbor D with B: 2
+DEAL::Vice-versa: 2
+DEAL::Is neighbor D with C: 3
+DEAL::Vice-versa: 3
+DEAL::Is neighbor D with D: 3
+DEAL::Bounding box boundaries E: 
+DEAL::-10.0000 -10.0000
+DEAL::-8.00000 -8.00000
+DEAL::Is neighbor E with A: 0
+DEAL::Vice-versa: 0
+DEAL::Is neighbor E with B: 0
+DEAL::Vice-versa: 0
+DEAL::Is neighbor E with C: 0
+DEAL::Vice-versa: 0
+DEAL::Is neighbor E with D: 0
+DEAL::Vice-versa: 0
+DEAL::Is neighbor E with E: 3
+DEAL::Bounding box boundaries F: 
+DEAL::1.00000 0.00000
+DEAL::2.00000 1.00000
+DEAL::Is neighbor F with A: 3
+DEAL::Vice-versa: 3
+DEAL::Is neighbor F with B: 3
+DEAL::Vice-versa: 3
+DEAL::Is neighbor F with C: 3
+DEAL::Vice-versa: 3
+DEAL::Is neighbor F with D: 2
+DEAL::Vice-versa: 2
+DEAL::Is neighbor F with E: 0
+DEAL::Vice-versa: 0
+DEAL::Is neighbor F with F: 3
+DEAL::End test for dimension 2
+DEAL::
+DEAL::Bounding box boundaries G: 
+DEAL::0.800000 0.00000
+DEAL::1.80000 1.00000
+DEAL::Is neighbor G with A: 3
+DEAL::Vice-versa: 3
+DEAL::Is neighbor G with B: 2
+DEAL::Vice-versa: 2
+DEAL::Is neighbor G with C: 3
+DEAL::Vice-versa: 3
+DEAL::Is neighbor G with D: 2
+DEAL::Vice-versa: 2
+DEAL::Is neighbor G with E: 0
+DEAL::Vice-versa: 0
+DEAL::Is neighbor G with F: 3
+DEAL::Vice-versa: 3
+DEAL::Is neighbor G with G: 3
+DEAL::End test for dimension 2
+DEAL::
+DEAL::
+DEAL::Test for dimension 3
+DEAL::Bounding box boundaries A: 
+DEAL::0.00000 0.00000 0.00000
+DEAL::1.00000 1.00000 1.00000
+DEAL::Is neighbor of itself: 3
+DEAL::Bounding box boundaries B: 
+DEAL::1.00000 1.00000 1.00000
+DEAL::2.00000 2.00000 2.00000
+DEAL::Is neighbor A with B: 1
+DEAL::Is neighbor B with A: 1
+DEAL::Is neighbor B with B: 3
+DEAL::Bounding box boundaries B1: 
+DEAL::0.00000 1.00000 1.00000
+DEAL::1.00000 2.00000 2.00000
+DEAL::Is neighbor A with B1: 1
+DEAL::Is neighbor B1 with A: 1
+DEAL::Is neighbor B1 with B: 3
+DEAL::Bounding box boundaries B2: 
+DEAL::0.200000 1.00000 1.00000
+DEAL::1.00000 2.00000 2.00000
+DEAL::Is neighbor A with B2: 1
+DEAL::Is neighbor B2 with A: 1
+DEAL::Is neighbor B2 with B: 3
+DEAL::Bounding box boundaries B3: 
+DEAL::0.00000 1.00000 0.800000
+DEAL::1.00000 2.00000 2.00000
+DEAL::Is neighbor A with B3: 2
+DEAL::Is neighbor B3 with A: 2
+DEAL::Is neighbor B3 with B: 3
+DEAL::Bounding box boundaries C: 
+DEAL::0.00000 0.00000 0.00000
+DEAL::2.00000 2.00000 2.00000
+DEAL::Is neighbor C with B: 3
+DEAL::Is neighbor B with C: 3
+DEAL::Is neighbor A with C: 3
+DEAL::Is neighbor C with A: 3
+DEAL::Is neighbor C with C: 3
+DEAL::Bounding box boundaries D: 
+DEAL::0.00000 0.00000 0.00000
+DEAL::1.40000 1.40000 1.40000
+DEAL::Is neighbor D with A: 3
+DEAL::Vice-versa: 3
+DEAL::Is neighbor D with B: 2
+DEAL::Vice-versa: 2
+DEAL::Is neighbor D with C: 3
+DEAL::Vice-versa: 3
+DEAL::Is neighbor D with D: 3
+DEAL::Bounding box boundaries E: 
+DEAL::-10.0000 -10.0000 -10.0000
+DEAL::-8.00000 -8.00000 -8.00000
+DEAL::Is neighbor E with A: 0
+DEAL::Vice-versa: 0
+DEAL::Is neighbor E with B: 0
+DEAL::Vice-versa: 0
+DEAL::Is neighbor E with C: 0
+DEAL::Vice-versa: 0
+DEAL::Is neighbor E with D: 0
+DEAL::Vice-versa: 0
+DEAL::Is neighbor E with E: 3
+DEAL::Bounding box boundaries F: 
+DEAL::1.00000 0.00000 0.00000
+DEAL::2.00000 1.00000 1.00000
+DEAL::Is neighbor F with A: 3
+DEAL::Vice-versa: 3
+DEAL::Is neighbor F with B: 1
+DEAL::Vice-versa: 1
+DEAL::Is neighbor F with C: 3
+DEAL::Vice-versa: 3
+DEAL::Is neighbor F with D: 2
+DEAL::Vice-versa: 2
+DEAL::Is neighbor F with E: 0
+DEAL::Vice-versa: 0
+DEAL::Is neighbor F with F: 3
+DEAL::End test for dimension 3
+DEAL::
+DEAL::Bounding box boundaries G: 
+DEAL::0.800000 0.00000 0.00000
+DEAL::1.80000 1.00000 1.00000
+DEAL::Is neighbor G with A: 3
+DEAL::Vice-versa: 3
+DEAL::Is neighbor G with B: 1
+DEAL::Vice-versa: 1
+DEAL::Is neighbor G with C: 3
+DEAL::Vice-versa: 3
+DEAL::Is neighbor G with D: 2
+DEAL::Vice-versa: 2
+DEAL::Is neighbor G with E: 0
+DEAL::Vice-versa: 0
+DEAL::Is neighbor G with F: 3
+DEAL::Vice-versa: 3
+DEAL::Is neighbor G with G: 3
+DEAL::End test for dimension 3
+DEAL::
diff --git a/tests/grid/grid_tools_compute_mesh_predicate_bounding_box_1.cc b/tests/grid/grid_tools_compute_mesh_predicate_bounding_box_1.cc
new file mode 100644 (file)
index 0000000..6b11375
--- /dev/null
@@ -0,0 +1,121 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2017 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.
+//
+// ---------------------------------------------------------------------
+
+
+// test for the function compute_mesh_predicate_bounding_box : as a predicate
+// is_locally_owned is used and, to vary the shapes, various mpi configurations
+// are used on a distributed::parallel::triangulations
+
+
+#include "../tests.h"
+#include <deal.II/base/logstream.h>
+#include <deal.II/base/bounding_box.h>
+#include <deal.II/base/mpi.h>
+#include <deal.II/base/point.h>
+#include <deal.II/base/utilities.h>
+#include <deal.II/distributed/tria.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/grid_tools.h>
+#include <deal.II/grid/manifold_lib.h>
+#include <deal.II/distributed/grid_tools.h>
+#include <deal.II/distributed/grid_refinement.h>
+
+// predicate: test if the cell is locally owned
+template <int spacedim>
+bool
+pred_locally_owned(const typename parallel::distributed::Triangulation<spacedim>::cell_iterator &cell)
+{
+  return cell->is_locally_owned();
+}
+
+template <int dim, int spacedim=dim >
+void test_hypercube(unsigned int ref, unsigned int max_bbox)
+{
+  const MPI_Comm &mpi_communicator = MPI_COMM_WORLD;
+  deallog << "Testing hypercube for spacedim = " << spacedim
+          << " refinement: " << ref << " max number of boxes: " << max_bbox << std::endl;
+
+  parallel::distributed::Triangulation<spacedim> tria(mpi_communicator);
+  GridGenerator::hyper_cube (tria);
+  tria.refine_global(4);
+
+  typedef typename  parallel::distributed::Triangulation<spacedim>::active_cell_iterator cell_iterator;
+
+  std::function<bool (const cell_iterator &)> predicate = pred_locally_owned<spacedim>;
+
+  std::vector< BoundingBox<spacedim> > local_bbox =
+    GridTools::compute_mesh_predicate_bounding_box < parallel::distributed::Triangulation<spacedim> >
+    ( tria, predicate, ref, true, max_bbox);
+
+  if ( local_bbox.size() > 0)
+    {
+      deallog << "Computed Bounding Boxes:" << std::endl;
+      for (BoundingBox<spacedim> b_box: local_bbox)
+        {
+          deallog << b_box.get_boundary_points().first << std::endl;
+          deallog << b_box.get_boundary_points().second << std::endl;
+          deallog << std::endl;
+        }
+    }
+  else
+    deallog << "Has no locally owned children cells" << std::endl;
+
+  //Checking if all the points are inside the bounding boxes
+  bool check = true;
+
+  typename parallel::distributed::Triangulation< dim, spacedim >::active_cell_iterator
+  cell = tria.begin_active();
+  typename parallel::distributed::Triangulation< dim, spacedim >::active_cell_iterator
+  endc = tria.last_active();
+
+  //Looking if every point is at least inside a bounding box
+  for (; cell<endc; ++cell)
+    if (cell->is_locally_owned())
+      for (unsigned int v=0; v<GeometryInfo<spacedim>::vertices_per_cell; ++v)
+        {
+          bool inside_a_box = false;
+          for (unsigned int i=0; i< local_bbox.size(); ++i)
+            {
+              if (local_bbox[i].point_inside(cell->vertex(v)))
+                inside_a_box = true;
+            }
+          if (! inside_a_box)
+            {
+              check = false;
+              deallog << "Point outside " << cell->vertex(v) << std::endl;
+              break;
+            }
+        }
+  deallog << "Bounding Boxes surround locally_owned cells: " << check << std::endl;
+
+  deallog << "Current test END"  << std::endl;
+}
+
+int main (int argc, char *argv[])
+{
+  Utilities::MPI::MPI_InitFinalize mpi_initialization (argc, argv, 1);
+  MPILogInitAll log;
+
+  test_hypercube<2> (0,4);
+  test_hypercube<2> (1,4);
+  test_hypercube<2> (2,4);
+  test_hypercube<2> (2,2);
+  test_hypercube<3> (0,4);
+  test_hypercube<3> (1,4);
+  test_hypercube<3> (2,4);
+  test_hypercube<3> (2,2);
+
+  return 0;
+}
diff --git a/tests/grid/grid_tools_compute_mesh_predicate_bounding_box_1.with_p4est=true.mpirun=2.output b/tests/grid/grid_tools_compute_mesh_predicate_bounding_box_1.with_p4est=true.mpirun=2.output
new file mode 100644 (file)
index 0000000..f146b0e
--- /dev/null
@@ -0,0 +1,115 @@
+
+DEAL:0::Testing hypercube for spacedim = 2 refinement: 0 max number of boxes: 4
+DEAL:0::Computed Bounding Boxes:
+DEAL:0::0.00000 0.00000
+DEAL:0::1.00000 0.500000
+DEAL:0::
+DEAL:0::Bounding Boxes surround locally_owned cells: 1
+DEAL:0::Current test END
+DEAL:0::Testing hypercube for spacedim = 2 refinement: 1 max number of boxes: 4
+DEAL:0::Computed Bounding Boxes:
+DEAL:0::0.00000 0.00000
+DEAL:0::1.00000 0.500000
+DEAL:0::
+DEAL:0::Bounding Boxes surround locally_owned cells: 1
+DEAL:0::Current test END
+DEAL:0::Testing hypercube for spacedim = 2 refinement: 2 max number of boxes: 4
+DEAL:0::Computed Bounding Boxes:
+DEAL:0::0.00000 0.00000
+DEAL:0::1.00000 0.500000
+DEAL:0::
+DEAL:0::Bounding Boxes surround locally_owned cells: 1
+DEAL:0::Current test END
+DEAL:0::Testing hypercube for spacedim = 2 refinement: 2 max number of boxes: 2
+DEAL:0::Computed Bounding Boxes:
+DEAL:0::0.00000 0.00000
+DEAL:0::1.00000 0.500000
+DEAL:0::
+DEAL:0::Bounding Boxes surround locally_owned cells: 1
+DEAL:0::Current test END
+DEAL:0::Testing hypercube for spacedim = 3 refinement: 0 max number of boxes: 4
+DEAL:0::Computed Bounding Boxes:
+DEAL:0::0.00000 0.00000 0.00000
+DEAL:0::1.00000 1.00000 0.500000
+DEAL:0::
+DEAL:0::Bounding Boxes surround locally_owned cells: 1
+DEAL:0::Current test END
+DEAL:0::Testing hypercube for spacedim = 3 refinement: 1 max number of boxes: 4
+DEAL:0::Computed Bounding Boxes:
+DEAL:0::0.00000 0.00000 0.00000
+DEAL:0::1.00000 1.00000 0.500000
+DEAL:0::
+DEAL:0::Bounding Boxes surround locally_owned cells: 1
+DEAL:0::Current test END
+DEAL:0::Testing hypercube for spacedim = 3 refinement: 2 max number of boxes: 4
+DEAL:0::Computed Bounding Boxes:
+DEAL:0::0.00000 0.00000 0.00000
+DEAL:0::1.00000 1.00000 0.500000
+DEAL:0::
+DEAL:0::Bounding Boxes surround locally_owned cells: 1
+DEAL:0::Current test END
+DEAL:0::Testing hypercube for spacedim = 3 refinement: 2 max number of boxes: 2
+DEAL:0::Computed Bounding Boxes:
+DEAL:0::0.00000 0.00000 0.00000
+DEAL:0::1.00000 1.00000 0.500000
+DEAL:0::
+DEAL:0::Bounding Boxes surround locally_owned cells: 1
+DEAL:0::Current test END
+
+DEAL:1::Testing hypercube for spacedim = 2 refinement: 0 max number of boxes: 4
+DEAL:1::Computed Bounding Boxes:
+DEAL:1::0.00000 0.500000
+DEAL:1::1.00000 1.00000
+DEAL:1::
+DEAL:1::Bounding Boxes surround locally_owned cells: 1
+DEAL:1::Current test END
+DEAL:1::Testing hypercube for spacedim = 2 refinement: 1 max number of boxes: 4
+DEAL:1::Computed Bounding Boxes:
+DEAL:1::0.00000 0.500000
+DEAL:1::1.00000 1.00000
+DEAL:1::
+DEAL:1::Bounding Boxes surround locally_owned cells: 1
+DEAL:1::Current test END
+DEAL:1::Testing hypercube for spacedim = 2 refinement: 2 max number of boxes: 4
+DEAL:1::Computed Bounding Boxes:
+DEAL:1::0.00000 0.500000
+DEAL:1::1.00000 1.00000
+DEAL:1::
+DEAL:1::Bounding Boxes surround locally_owned cells: 1
+DEAL:1::Current test END
+DEAL:1::Testing hypercube for spacedim = 2 refinement: 2 max number of boxes: 2
+DEAL:1::Computed Bounding Boxes:
+DEAL:1::0.00000 0.500000
+DEAL:1::1.00000 1.00000
+DEAL:1::
+DEAL:1::Bounding Boxes surround locally_owned cells: 1
+DEAL:1::Current test END
+DEAL:1::Testing hypercube for spacedim = 3 refinement: 0 max number of boxes: 4
+DEAL:1::Computed Bounding Boxes:
+DEAL:1::0.00000 0.00000 0.500000
+DEAL:1::1.00000 1.00000 1.00000
+DEAL:1::
+DEAL:1::Bounding Boxes surround locally_owned cells: 1
+DEAL:1::Current test END
+DEAL:1::Testing hypercube for spacedim = 3 refinement: 1 max number of boxes: 4
+DEAL:1::Computed Bounding Boxes:
+DEAL:1::0.00000 0.00000 0.500000
+DEAL:1::1.00000 1.00000 1.00000
+DEAL:1::
+DEAL:1::Bounding Boxes surround locally_owned cells: 1
+DEAL:1::Current test END
+DEAL:1::Testing hypercube for spacedim = 3 refinement: 2 max number of boxes: 4
+DEAL:1::Computed Bounding Boxes:
+DEAL:1::0.00000 0.00000 0.500000
+DEAL:1::1.00000 1.00000 1.00000
+DEAL:1::
+DEAL:1::Bounding Boxes surround locally_owned cells: 1
+DEAL:1::Current test END
+DEAL:1::Testing hypercube for spacedim = 3 refinement: 2 max number of boxes: 2
+DEAL:1::Computed Bounding Boxes:
+DEAL:1::0.00000 0.00000 0.500000
+DEAL:1::1.00000 1.00000 1.00000
+DEAL:1::
+DEAL:1::Bounding Boxes surround locally_owned cells: 1
+DEAL:1::Current test END
+
diff --git a/tests/grid/grid_tools_compute_mesh_predicate_bounding_box_1.with_p4est=true.mpirun=3.output b/tests/grid/grid_tools_compute_mesh_predicate_bounding_box_1.with_p4est=true.mpirun=3.output
new file mode 100644 (file)
index 0000000..5969a38
--- /dev/null
@@ -0,0 +1,263 @@
+
+DEAL:0::Testing hypercube for spacedim = 2 refinement: 0 max number of boxes: 4
+DEAL:0::Computed Bounding Boxes:
+DEAL:0::0.00000 0.00000
+DEAL:0::0.875000 0.500000
+DEAL:0::
+DEAL:0::Bounding Boxes surround locally_owned cells: 1
+DEAL:0::Current test END
+DEAL:0::Testing hypercube for spacedim = 2 refinement: 1 max number of boxes: 4
+DEAL:0::Computed Bounding Boxes:
+DEAL:0::0.00000 0.00000
+DEAL:0::0.500000 0.500000
+DEAL:0::
+DEAL:0::0.500000 0.00000
+DEAL:0::0.875000 0.250000
+DEAL:0::
+DEAL:0::Bounding Boxes surround locally_owned cells: 1
+DEAL:0::Current test END
+DEAL:0::Testing hypercube for spacedim = 2 refinement: 2 max number of boxes: 4
+DEAL:0::Computed Bounding Boxes:
+DEAL:0::0.00000 0.00000
+DEAL:0::0.750000 0.250000
+DEAL:0::
+DEAL:0::0.00000 0.250000
+DEAL:0::0.500000 0.500000
+DEAL:0::
+DEAL:0::0.750000 0.00000
+DEAL:0::0.875000 0.125000
+DEAL:0::
+DEAL:0::Bounding Boxes surround locally_owned cells: 1
+DEAL:0::Current test END
+DEAL:0::Testing hypercube for spacedim = 2 refinement: 2 max number of boxes: 2
+DEAL:0::Computed Bounding Boxes:
+DEAL:0::0.00000 0.00000
+DEAL:0::0.875000 0.250000
+DEAL:0::
+DEAL:0::0.00000 0.250000
+DEAL:0::0.500000 0.500000
+DEAL:0::
+DEAL:0::Bounding Boxes surround locally_owned cells: 1
+DEAL:0::Current test END
+DEAL:0::Testing hypercube for spacedim = 3 refinement: 0 max number of boxes: 4
+DEAL:0::Computed Bounding Boxes:
+DEAL:0::0.00000 0.00000 0.00000
+DEAL:0::1.00000 1.00000 0.500000
+DEAL:0::
+DEAL:0::Bounding Boxes surround locally_owned cells: 1
+DEAL:0::Current test END
+DEAL:0::Testing hypercube for spacedim = 3 refinement: 1 max number of boxes: 4
+DEAL:0::Computed Bounding Boxes:
+DEAL:0::0.00000 0.00000 0.00000
+DEAL:0::1.00000 0.500000 0.500000
+DEAL:0::
+DEAL:0::0.00000 0.500000 0.00000
+DEAL:0::0.500000 1.00000 0.500000
+DEAL:0::
+DEAL:0::Bounding Boxes surround locally_owned cells: 1
+DEAL:0::Current test END
+DEAL:0::Testing hypercube for spacedim = 3 refinement: 2 max number of boxes: 4
+DEAL:0::Computed Bounding Boxes:
+DEAL:0::0.00000 0.00000 0.00000
+DEAL:0::1.00000 0.500000 0.500000
+DEAL:0::
+DEAL:0::0.00000 0.500000 0.00000
+DEAL:0::0.500000 1.00000 0.250000
+DEAL:0::
+DEAL:0::0.00000 0.500000 0.250000
+DEAL:0::0.250000 0.750000 0.500000
+DEAL:0::
+DEAL:0::0.250000 0.500000 0.250000
+DEAL:0::0.500000 0.750000 0.375000
+DEAL:0::
+DEAL:0::Bounding Boxes surround locally_owned cells: 1
+DEAL:0::Current test END
+DEAL:0::Testing hypercube for spacedim = 3 refinement: 2 max number of boxes: 2
+DEAL:0::Computed Bounding Boxes:
+DEAL:0::0.00000 0.00000 0.00000
+DEAL:0::1.00000 0.750000 0.500000
+DEAL:0::
+DEAL:0::0.00000 0.500000 0.00000
+DEAL:0::0.500000 1.00000 0.250000
+DEAL:0::
+DEAL:0::Bounding Boxes surround locally_owned cells: 1
+DEAL:0::Current test END
+
+DEAL:1::Testing hypercube for spacedim = 2 refinement: 0 max number of boxes: 4
+DEAL:1::Computed Bounding Boxes:
+DEAL:1::0.00000 0.00000
+DEAL:1::1.00000 1.00000
+DEAL:1::
+DEAL:1::Bounding Boxes surround locally_owned cells: 1
+DEAL:1::Current test END
+DEAL:1::Testing hypercube for spacedim = 2 refinement: 1 max number of boxes: 4
+DEAL:1::Computed Bounding Boxes:
+DEAL:1::0.500000 0.00000
+DEAL:1::1.00000 0.500000
+DEAL:1::
+DEAL:1::0.00000 0.500000
+DEAL:1::0.500000 1.00000
+DEAL:1::
+DEAL:1::Bounding Boxes surround locally_owned cells: 1
+DEAL:1::Current test END
+DEAL:1::Testing hypercube for spacedim = 2 refinement: 2 max number of boxes: 4
+DEAL:1::Computed Bounding Boxes:
+DEAL:1::0.750000 0.00000
+DEAL:1::1.00000 0.500000
+DEAL:1::
+DEAL:1::0.500000 0.250000
+DEAL:1::0.750000 0.500000
+DEAL:1::
+DEAL:1::0.00000 0.500000
+DEAL:1::0.500000 0.750000
+DEAL:1::
+DEAL:1::0.00000 0.750000
+DEAL:1::0.250000 1.00000
+DEAL:1::
+DEAL:1::Bounding Boxes surround locally_owned cells: 1
+DEAL:1::Current test END
+DEAL:1::Testing hypercube for spacedim = 2 refinement: 2 max number of boxes: 2
+DEAL:1::Computed Bounding Boxes:
+DEAL:1::0.500000 0.00000
+DEAL:1::1.00000 0.500000
+DEAL:1::
+DEAL:1::0.00000 0.500000
+DEAL:1::0.500000 1.00000
+DEAL:1::
+DEAL:1::Bounding Boxes surround locally_owned cells: 1
+DEAL:1::Current test END
+DEAL:1::Testing hypercube for spacedim = 3 refinement: 0 max number of boxes: 4
+DEAL:1::Computed Bounding Boxes:
+DEAL:1::0.00000 0.00000 0.00000
+DEAL:1::1.00000 1.00000 1.00000
+DEAL:1::
+DEAL:1::Bounding Boxes surround locally_owned cells: 1
+DEAL:1::Current test END
+DEAL:1::Testing hypercube for spacedim = 3 refinement: 1 max number of boxes: 4
+DEAL:1::Computed Bounding Boxes:
+DEAL:1::0.00000 0.500000 0.250000
+DEAL:1::0.500000 1.00000 0.500000
+DEAL:1::
+DEAL:1::0.500000 0.500000 0.00000
+DEAL:1::1.00000 1.00000 0.500000
+DEAL:1::
+DEAL:1::0.00000 0.00000 0.500000
+DEAL:1::0.500000 0.500000 1.00000
+DEAL:1::
+DEAL:1::0.500000 0.00000 0.500000
+DEAL:1::1.00000 0.500000 0.750000
+DEAL:1::
+DEAL:1::Bounding Boxes surround locally_owned cells: 1
+DEAL:1::Current test END
+DEAL:1::Testing hypercube for spacedim = 3 refinement: 2 max number of boxes: 4
+DEAL:1::Computed Bounding Boxes:
+DEAL:1::0.00000 0.500000 0.250000
+DEAL:1::1.00000 1.00000 0.500000
+DEAL:1::
+DEAL:1::0.500000 0.500000 0.00000
+DEAL:1::1.00000 1.00000 0.250000
+DEAL:1::
+DEAL:1::0.00000 0.00000 0.500000
+DEAL:1::1.00000 0.500000 0.750000
+DEAL:1::
+DEAL:1::0.00000 0.00000 0.750000
+DEAL:1::0.500000 0.500000 1.00000
+DEAL:1::
+DEAL:1::Bounding Boxes surround locally_owned cells: 1
+DEAL:1::Current test END
+DEAL:1::Testing hypercube for spacedim = 3 refinement: 2 max number of boxes: 2
+DEAL:1::Computed Bounding Boxes:
+DEAL:1::0.00000 0.500000 0.00000
+DEAL:1::1.00000 1.00000 0.500000
+DEAL:1::
+DEAL:1::0.00000 0.00000 0.500000
+DEAL:1::1.00000 0.500000 1.00000
+DEAL:1::
+DEAL:1::Bounding Boxes surround locally_owned cells: 1
+DEAL:1::Current test END
+
+
+DEAL:2::Testing hypercube for spacedim = 2 refinement: 0 max number of boxes: 4
+DEAL:2::Computed Bounding Boxes:
+DEAL:2::0.125000 0.500000
+DEAL:2::1.00000 1.00000
+DEAL:2::
+DEAL:2::Bounding Boxes surround locally_owned cells: 1
+DEAL:2::Current test END
+DEAL:2::Testing hypercube for spacedim = 2 refinement: 1 max number of boxes: 4
+DEAL:2::Computed Bounding Boxes:
+DEAL:2::0.125000 0.750000
+DEAL:2::0.500000 1.00000
+DEAL:2::
+DEAL:2::0.500000 0.500000
+DEAL:2::1.00000 1.00000
+DEAL:2::
+DEAL:2::Bounding Boxes surround locally_owned cells: 1
+DEAL:2::Current test END
+DEAL:2::Testing hypercube for spacedim = 2 refinement: 2 max number of boxes: 4
+DEAL:2::Computed Bounding Boxes:
+DEAL:2::0.125000 0.875000
+DEAL:2::0.250000 1.00000
+DEAL:2::
+DEAL:2::0.250000 0.750000
+DEAL:2::1.00000 1.00000
+DEAL:2::
+DEAL:2::0.500000 0.500000
+DEAL:2::1.00000 0.750000
+DEAL:2::
+DEAL:2::Bounding Boxes surround locally_owned cells: 1
+DEAL:2::Current test END
+DEAL:2::Testing hypercube for spacedim = 2 refinement: 2 max number of boxes: 2
+DEAL:2::Computed Bounding Boxes:
+DEAL:2::0.125000 0.750000
+DEAL:2::1.00000 1.00000
+DEAL:2::
+DEAL:2::0.500000 0.500000
+DEAL:2::1.00000 0.750000
+DEAL:2::
+DEAL:2::Bounding Boxes surround locally_owned cells: 1
+DEAL:2::Current test END
+DEAL:2::Testing hypercube for spacedim = 3 refinement: 0 max number of boxes: 4
+DEAL:2::Computed Bounding Boxes:
+DEAL:2::0.00000 0.00000 0.500000
+DEAL:2::1.00000 1.00000 1.00000
+DEAL:2::
+DEAL:2::Bounding Boxes surround locally_owned cells: 1
+DEAL:2::Current test END
+DEAL:2::Testing hypercube for spacedim = 3 refinement: 1 max number of boxes: 4
+DEAL:2::Computed Bounding Boxes:
+DEAL:2::0.500000 0.00000 0.500000
+DEAL:2::1.00000 1.00000 1.00000
+DEAL:2::
+DEAL:2::0.00000 0.500000 0.500000
+DEAL:2::0.500000 1.00000 1.00000
+DEAL:2::
+DEAL:2::Bounding Boxes surround locally_owned cells: 1
+DEAL:2::Current test END
+DEAL:2::Testing hypercube for spacedim = 3 refinement: 2 max number of boxes: 4
+DEAL:2::Computed Bounding Boxes:
+DEAL:2::0.500000 0.250000 0.625000
+DEAL:2::0.750000 0.500000 1.00000
+DEAL:2::
+DEAL:2::0.750000 0.250000 0.500000
+DEAL:2::1.00000 0.500000 1.00000
+DEAL:2::
+DEAL:2::0.500000 0.00000 0.750000
+DEAL:2::1.00000 0.250000 1.00000
+DEAL:2::
+DEAL:2::0.00000 0.500000 0.500000
+DEAL:2::1.00000 1.00000 1.00000
+DEAL:2::
+DEAL:2::Bounding Boxes surround locally_owned cells: 1
+DEAL:2::Current test END
+DEAL:2::Testing hypercube for spacedim = 3 refinement: 2 max number of boxes: 2
+DEAL:2::Computed Bounding Boxes:
+DEAL:2::0.500000 0.00000 0.500000
+DEAL:2::1.00000 0.500000 1.00000
+DEAL:2::
+DEAL:2::0.00000 0.500000 0.500000
+DEAL:2::1.00000 1.00000 1.00000
+DEAL:2::
+DEAL:2::Bounding Boxes surround locally_owned cells: 1
+DEAL:2::Current test END
+
diff --git a/tests/grid/grid_tools_compute_mesh_predicate_bounding_box_1.with_p4est=true.mpirun=4.output b/tests/grid/grid_tools_compute_mesh_predicate_bounding_box_1.with_p4est=true.mpirun=4.output
new file mode 100644 (file)
index 0000000..660bd9f
--- /dev/null
@@ -0,0 +1,231 @@
+
+DEAL:0::Testing hypercube for spacedim = 2 refinement: 0 max number of boxes: 4
+DEAL:0::Computed Bounding Boxes:
+DEAL:0::0.00000 0.00000
+DEAL:0::0.500000 0.500000
+DEAL:0::
+DEAL:0::Bounding Boxes surround locally_owned cells: 1
+DEAL:0::Current test END
+DEAL:0::Testing hypercube for spacedim = 2 refinement: 1 max number of boxes: 4
+DEAL:0::Computed Bounding Boxes:
+DEAL:0::0.00000 0.00000
+DEAL:0::0.500000 0.500000
+DEAL:0::
+DEAL:0::Bounding Boxes surround locally_owned cells: 1
+DEAL:0::Current test END
+DEAL:0::Testing hypercube for spacedim = 2 refinement: 2 max number of boxes: 4
+DEAL:0::Computed Bounding Boxes:
+DEAL:0::0.00000 0.00000
+DEAL:0::0.500000 0.500000
+DEAL:0::
+DEAL:0::Bounding Boxes surround locally_owned cells: 1
+DEAL:0::Current test END
+DEAL:0::Testing hypercube for spacedim = 2 refinement: 2 max number of boxes: 2
+DEAL:0::Computed Bounding Boxes:
+DEAL:0::0.00000 0.00000
+DEAL:0::0.500000 0.500000
+DEAL:0::
+DEAL:0::Bounding Boxes surround locally_owned cells: 1
+DEAL:0::Current test END
+DEAL:0::Testing hypercube for spacedim = 3 refinement: 0 max number of boxes: 4
+DEAL:0::Computed Bounding Boxes:
+DEAL:0::0.00000 0.00000 0.00000
+DEAL:0::1.00000 0.500000 0.500000
+DEAL:0::
+DEAL:0::Bounding Boxes surround locally_owned cells: 1
+DEAL:0::Current test END
+DEAL:0::Testing hypercube for spacedim = 3 refinement: 1 max number of boxes: 4
+DEAL:0::Computed Bounding Boxes:
+DEAL:0::0.00000 0.00000 0.00000
+DEAL:0::1.00000 0.500000 0.500000
+DEAL:0::
+DEAL:0::Bounding Boxes surround locally_owned cells: 1
+DEAL:0::Current test END
+DEAL:0::Testing hypercube for spacedim = 3 refinement: 2 max number of boxes: 4
+DEAL:0::Computed Bounding Boxes:
+DEAL:0::0.00000 0.00000 0.00000
+DEAL:0::1.00000 0.500000 0.500000
+DEAL:0::
+DEAL:0::Bounding Boxes surround locally_owned cells: 1
+DEAL:0::Current test END
+DEAL:0::Testing hypercube for spacedim = 3 refinement: 2 max number of boxes: 2
+DEAL:0::Computed Bounding Boxes:
+DEAL:0::0.00000 0.00000 0.00000
+DEAL:0::1.00000 0.500000 0.500000
+DEAL:0::
+DEAL:0::Bounding Boxes surround locally_owned cells: 1
+DEAL:0::Current test END
+
+DEAL:1::Testing hypercube for spacedim = 2 refinement: 0 max number of boxes: 4
+DEAL:1::Computed Bounding Boxes:
+DEAL:1::0.500000 0.00000
+DEAL:1::1.00000 0.500000
+DEAL:1::
+DEAL:1::Bounding Boxes surround locally_owned cells: 1
+DEAL:1::Current test END
+DEAL:1::Testing hypercube for spacedim = 2 refinement: 1 max number of boxes: 4
+DEAL:1::Computed Bounding Boxes:
+DEAL:1::0.500000 0.00000
+DEAL:1::1.00000 0.500000
+DEAL:1::
+DEAL:1::Bounding Boxes surround locally_owned cells: 1
+DEAL:1::Current test END
+DEAL:1::Testing hypercube for spacedim = 2 refinement: 2 max number of boxes: 4
+DEAL:1::Computed Bounding Boxes:
+DEAL:1::0.500000 0.00000
+DEAL:1::1.00000 0.500000
+DEAL:1::
+DEAL:1::Bounding Boxes surround locally_owned cells: 1
+DEAL:1::Current test END
+DEAL:1::Testing hypercube for spacedim = 2 refinement: 2 max number of boxes: 2
+DEAL:1::Computed Bounding Boxes:
+DEAL:1::0.500000 0.00000
+DEAL:1::1.00000 0.500000
+DEAL:1::
+DEAL:1::Bounding Boxes surround locally_owned cells: 1
+DEAL:1::Current test END
+DEAL:1::Testing hypercube for spacedim = 3 refinement: 0 max number of boxes: 4
+DEAL:1::Computed Bounding Boxes:
+DEAL:1::0.00000 0.500000 0.00000
+DEAL:1::1.00000 1.00000 0.500000
+DEAL:1::
+DEAL:1::Bounding Boxes surround locally_owned cells: 1
+DEAL:1::Current test END
+DEAL:1::Testing hypercube for spacedim = 3 refinement: 1 max number of boxes: 4
+DEAL:1::Computed Bounding Boxes:
+DEAL:1::0.00000 0.500000 0.00000
+DEAL:1::1.00000 1.00000 0.500000
+DEAL:1::
+DEAL:1::Bounding Boxes surround locally_owned cells: 1
+DEAL:1::Current test END
+DEAL:1::Testing hypercube for spacedim = 3 refinement: 2 max number of boxes: 4
+DEAL:1::Computed Bounding Boxes:
+DEAL:1::0.00000 0.500000 0.00000
+DEAL:1::1.00000 1.00000 0.500000
+DEAL:1::
+DEAL:1::Bounding Boxes surround locally_owned cells: 1
+DEAL:1::Current test END
+DEAL:1::Testing hypercube for spacedim = 3 refinement: 2 max number of boxes: 2
+DEAL:1::Computed Bounding Boxes:
+DEAL:1::0.00000 0.500000 0.00000
+DEAL:1::1.00000 1.00000 0.500000
+DEAL:1::
+DEAL:1::Bounding Boxes surround locally_owned cells: 1
+DEAL:1::Current test END
+
+
+DEAL:2::Testing hypercube for spacedim = 2 refinement: 0 max number of boxes: 4
+DEAL:2::Computed Bounding Boxes:
+DEAL:2::0.00000 0.500000
+DEAL:2::0.500000 1.00000
+DEAL:2::
+DEAL:2::Bounding Boxes surround locally_owned cells: 1
+DEAL:2::Current test END
+DEAL:2::Testing hypercube for spacedim = 2 refinement: 1 max number of boxes: 4
+DEAL:2::Computed Bounding Boxes:
+DEAL:2::0.00000 0.500000
+DEAL:2::0.500000 1.00000
+DEAL:2::
+DEAL:2::Bounding Boxes surround locally_owned cells: 1
+DEAL:2::Current test END
+DEAL:2::Testing hypercube for spacedim = 2 refinement: 2 max number of boxes: 4
+DEAL:2::Computed Bounding Boxes:
+DEAL:2::0.00000 0.500000
+DEAL:2::0.500000 1.00000
+DEAL:2::
+DEAL:2::Bounding Boxes surround locally_owned cells: 1
+DEAL:2::Current test END
+DEAL:2::Testing hypercube for spacedim = 2 refinement: 2 max number of boxes: 2
+DEAL:2::Computed Bounding Boxes:
+DEAL:2::0.00000 0.500000
+DEAL:2::0.500000 1.00000
+DEAL:2::
+DEAL:2::Bounding Boxes surround locally_owned cells: 1
+DEAL:2::Current test END
+DEAL:2::Testing hypercube for spacedim = 3 refinement: 0 max number of boxes: 4
+DEAL:2::Computed Bounding Boxes:
+DEAL:2::0.00000 0.00000 0.500000
+DEAL:2::1.00000 0.500000 1.00000
+DEAL:2::
+DEAL:2::Bounding Boxes surround locally_owned cells: 1
+DEAL:2::Current test END
+DEAL:2::Testing hypercube for spacedim = 3 refinement: 1 max number of boxes: 4
+DEAL:2::Computed Bounding Boxes:
+DEAL:2::0.00000 0.00000 0.500000
+DEAL:2::1.00000 0.500000 1.00000
+DEAL:2::
+DEAL:2::Bounding Boxes surround locally_owned cells: 1
+DEAL:2::Current test END
+DEAL:2::Testing hypercube for spacedim = 3 refinement: 2 max number of boxes: 4
+DEAL:2::Computed Bounding Boxes:
+DEAL:2::0.00000 0.00000 0.500000
+DEAL:2::1.00000 0.500000 1.00000
+DEAL:2::
+DEAL:2::Bounding Boxes surround locally_owned cells: 1
+DEAL:2::Current test END
+DEAL:2::Testing hypercube for spacedim = 3 refinement: 2 max number of boxes: 2
+DEAL:2::Computed Bounding Boxes:
+DEAL:2::0.00000 0.00000 0.500000
+DEAL:2::1.00000 0.500000 1.00000
+DEAL:2::
+DEAL:2::Bounding Boxes surround locally_owned cells: 1
+DEAL:2::Current test END
+
+
+DEAL:3::Testing hypercube for spacedim = 2 refinement: 0 max number of boxes: 4
+DEAL:3::Computed Bounding Boxes:
+DEAL:3::0.500000 0.500000
+DEAL:3::1.00000 1.00000
+DEAL:3::
+DEAL:3::Bounding Boxes surround locally_owned cells: 1
+DEAL:3::Current test END
+DEAL:3::Testing hypercube for spacedim = 2 refinement: 1 max number of boxes: 4
+DEAL:3::Computed Bounding Boxes:
+DEAL:3::0.500000 0.500000
+DEAL:3::1.00000 1.00000
+DEAL:3::
+DEAL:3::Bounding Boxes surround locally_owned cells: 1
+DEAL:3::Current test END
+DEAL:3::Testing hypercube for spacedim = 2 refinement: 2 max number of boxes: 4
+DEAL:3::Computed Bounding Boxes:
+DEAL:3::0.500000 0.500000
+DEAL:3::1.00000 1.00000
+DEAL:3::
+DEAL:3::Bounding Boxes surround locally_owned cells: 1
+DEAL:3::Current test END
+DEAL:3::Testing hypercube for spacedim = 2 refinement: 2 max number of boxes: 2
+DEAL:3::Computed Bounding Boxes:
+DEAL:3::0.500000 0.500000
+DEAL:3::1.00000 1.00000
+DEAL:3::
+DEAL:3::Bounding Boxes surround locally_owned cells: 1
+DEAL:3::Current test END
+DEAL:3::Testing hypercube for spacedim = 3 refinement: 0 max number of boxes: 4
+DEAL:3::Computed Bounding Boxes:
+DEAL:3::0.00000 0.500000 0.500000
+DEAL:3::1.00000 1.00000 1.00000
+DEAL:3::
+DEAL:3::Bounding Boxes surround locally_owned cells: 1
+DEAL:3::Current test END
+DEAL:3::Testing hypercube for spacedim = 3 refinement: 1 max number of boxes: 4
+DEAL:3::Computed Bounding Boxes:
+DEAL:3::0.00000 0.500000 0.500000
+DEAL:3::1.00000 1.00000 1.00000
+DEAL:3::
+DEAL:3::Bounding Boxes surround locally_owned cells: 1
+DEAL:3::Current test END
+DEAL:3::Testing hypercube for spacedim = 3 refinement: 2 max number of boxes: 4
+DEAL:3::Computed Bounding Boxes:
+DEAL:3::0.00000 0.500000 0.500000
+DEAL:3::1.00000 1.00000 1.00000
+DEAL:3::
+DEAL:3::Bounding Boxes surround locally_owned cells: 1
+DEAL:3::Current test END
+DEAL:3::Testing hypercube for spacedim = 3 refinement: 2 max number of boxes: 2
+DEAL:3::Computed Bounding Boxes:
+DEAL:3::0.00000 0.500000 0.500000
+DEAL:3::1.00000 1.00000 1.00000
+DEAL:3::
+DEAL:3::Bounding Boxes surround locally_owned cells: 1
+DEAL:3::Current test END
+

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.