]> https://gitweb.dealii.org/ - dealii.git/commitdiff
GridTools::compute_bounding_box(), compute_active_cells_within_skin(), compute_ghost_... 4093/head
authorVishal Boddu <vishal.boddu@fau.de>
Sat, 29 Apr 2017 17:16:16 +0000 (19:16 +0200)
committerVishal Boddu <vishal.boddu@fau.de>
Sat, 29 Apr 2017 17:16:16 +0000 (19:16 +0200)
doc/doxygen/images/active_cell_layer_within_distance.png [new file with mode: 0644]
doc/news/changes/minor/20170323VishalBoddu [new file with mode: 0644]
include/deal.II/grid/grid_tools.h
source/grid/grid_tools.cc
source/grid/grid_tools.inst.in
tests/grid/grid_tools_active_cell_layer_within_distance_01.cc [new file with mode: 0644]
tests/grid/grid_tools_active_cell_layer_within_distance_01.output [new file with mode: 0644]
tests/grid/grid_tools_active_cell_layer_within_distance_02.cc [new file with mode: 0644]
tests/grid/grid_tools_active_cell_layer_within_distance_02.output [new file with mode: 0644]
tests/grid/grid_tools_bounding_box_01.cc [new file with mode: 0644]
tests/grid/grid_tools_bounding_box_01.output [new file with mode: 0644]

diff --git a/doc/doxygen/images/active_cell_layer_within_distance.png b/doc/doxygen/images/active_cell_layer_within_distance.png
new file mode 100644 (file)
index 0000000..a73e0dd
Binary files /dev/null and b/doc/doxygen/images/active_cell_layer_within_distance.png differ
diff --git a/doc/news/changes/minor/20170323VishalBoddu b/doc/news/changes/minor/20170323VishalBoddu
new file mode 100644 (file)
index 0000000..ce1b584
--- /dev/null
@@ -0,0 +1,10 @@
+New: GridTools::compute_bounding_box() computes a bounding box of a subdomain
+whose active cells conform to a given predicate. 
+GridTools::compute_active_cell_layer_within_distance() computes a collection
+of active cells that are within a given distance from  
+predicate subdomain.
+GridTools::compute_ghost_cell_layer_within_distance() computes a collection
+of active cells that are within a given distance from 
+locally owned active cells.
+<br>
+(Vishal Boddu, Denis Davydov 2017/03/23)
index ba7889c48354741cd00d57806a648e9f18d144ee..d1db8be1d5fc51475fab45abd8e8f5ec9f68a408 100644 (file)
@@ -813,6 +813,110 @@ namespace GridTools
   std::vector<typename MeshType::active_cell_iterator>
   compute_ghost_cell_halo_layer (const MeshType &mesh);
 
+  /**
+   * Extract and return the set of active cells within a geometric distance of
+   * @p layer_thickness around a subdomain (set of active cells) in the @p mesh.
+   * Here, the "subdomain" consists of exactly all of
+   * those cells for which the @p predicate returns @p true.
+   *
+   * The function first computes the cells that form the 'surface' of the
+   * subdomain that consists of all of the active cells for which the predicate
+   * is true. Using compute_bounding_box(), a bounding box is
+   * computed for this subdomain and extended by @p layer_thickness. These
+   * cells are called interior subdomain boundary cells.
+   * The active cells with all of their vertices outside the extended
+   * bounding box are ignored.
+   * The cells that are inside the extended bounding box are then checked for
+   * their proximity to the interior subdomain boundary cells. This implies
+   * checking the distance between a pair of arbitrarily oriented cells,
+   * which is not trivial in general. To simplify this, the algorithm checks
+   * the distance between the two enclosing spheres of the cells.
+   * This will definitely result in slightly more cells being marked but
+   * also greatly simplifies the arithmetic complexity of the algorithm.
+   *
+   * @image html active_cell_layer_within_distance.png
+   * The image shows a mesh generated by subdivided_hyper_rectangle(). The cells
+   * are marked using three different colors. If the grey colored cells in the
+   * image are the cells for which the predicate is true, then the function
+   * compute_active_cell_layer_within_distance() will return a set of cell
+   * iterators corresponding to the cells colored in red.
+   * The red colored cells are the active cells that are within a given
+   * distance to the grey colored cells.
+   *
+   * @tparam MeshType A type that satisfies the requirements of the
+   * @ref ConceptMeshType "MeshType concept".
+   * @param mesh A mesh (i.e. objects of type Triangulation, DoFHandler,
+   * or hp::DoFHandler).
+   * @param predicate A function  (or object of a type with an operator())
+   * defining the subdomain around which the halo layer is to be extracted. It
+   * is a function that takes in an active cell and returns a boolean.
+   * @param layer_thickness specifies the geometric distance within
+   * which the function searches for active cells from the predicate domain.
+   * If the minimal distance between the enclosing sphere of the an
+   * active cell and the enclosing sphere of any of the cells for which
+   * the @p predicate returns @p true is less than @p layer_thickness,
+   * then the active cell is an \a active_cell_wthin_distance.
+   * @return A list of active cells within a given geometric distance
+   * @p layer_thickness from the set of active cells for which the @p predicate
+   * returns @p true.
+   *
+   * See compute_active_cell_halo_layer().
+   *
+   * @author Vishal Boddu, Denis Davydov, 2017
+   */
+  template <class MeshType>
+  std::vector<typename MeshType::active_cell_iterator>
+  compute_active_cell_layer_within_distance
+  (const MeshType                                                                    &mesh,
+   const std::function<bool (const typename MeshType::active_cell_iterator &)> &predicate,
+   const double                                                                       layer_thickness);
+
+  /**
+   * Extract and return a set of ghost cells which are within a
+   * @p layer_thickness around all locally owned cells.
+   * This is most relevant for parallel::shared::Triangulation
+   * where it will return a subset of all ghost cells on a process, but for
+   * parallel::distributed::Triangulation this will return all the ghost cells.
+   * All the cells for the parallel::shared::Triangulation class that
+   * are not owned by the current processor can be considered as ghost cells;
+   * in particular, they do not only form a single layer of cells around the
+   * locally owned ones.
+   *
+   * @tparam MeshType A type that satisfies the requirements of the
+   * @ref ConceptMeshType "MeshType concept".
+   * @param mesh A mesh (i.e. objects of type Triangulation, DoFHandler,
+   * or hp::DoFHandler).
+   * @param layer_thickness specifies the geometric distance within
+   * which the function searches for active cells from the locally owned cells.
+   * @return A subset of ghost cells within a given geometric distance of @p
+   * layer_thickness from the locally owned cells of a current process.
+   *
+   * Also see compute_ghost_cell_halo_layer() and
+   * compute_active_cell_layer_within_distance().
+   *
+   * @author Vishal Boddu, Denis Davydov, 2017
+   */
+  template <int spacedim, class MeshType>
+  std::vector<typename MeshType::active_cell_iterator>
+  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.
+   * Here, the "subdomain" consists of exactly all of those
+   * active cells for which the @p predicate returns @p true.
+   *
+   * For a description of how predicate works,
+   * see compute_active_cell_halo_layer()
+   */
+  template<class MeshType>
+  std::pair< Point<MeshType::space_dimension>, Point<MeshType::space_dimension> >
+  compute_bounding_box
+  ( const MeshType                                                                    &mesh,
+    const std::function<bool (const typename MeshType::active_cell_iterator &)> &predicate );
+
 
   /**
    * Return the adjacent cells of all the vertices. If a vertex is also a
index c3e19f32abe089d0822dc8ff2f595d1aa4cd2560..9c2a764e7a04bf1569df1ed5be608d4e27a240d5 100644 (file)
@@ -1716,6 +1716,208 @@ next_cell:
 
 
 
+  template <class MeshType>
+  std::vector<typename MeshType::active_cell_iterator>
+  compute_active_cell_layer_within_distance
+  (const MeshType                                                                    &mesh,
+   const std::function<bool (const typename MeshType::active_cell_iterator &)> &predicate,
+   const double                                                                       layer_thickness)
+  {
+    std::vector<typename MeshType::active_cell_iterator> subdomain_boundary_cells, active_cell_layer_within_distance;
+    std::vector<bool> vertices_outside_subdomain ( mesh.get_triangulation().n_vertices(),
+                                                   false);
+
+    const unsigned int spacedim = MeshType::space_dimension;
+
+    unsigned int n_non_predicate_cells = 0; // Number of non predicate cells
+
+    // Find the layer of cells for which predicate is true and that
+    // are on the boundary with other cells. These are
+    // subdomain boundary cells.
+
+    // Find the cells for which the predicate is false
+    // These are the cells which are around the predicate subdomain
+    for ( typename MeshType::active_cell_iterator
+          cell = mesh.begin_active();
+          cell != mesh.end(); ++cell)
+      if ( !predicate(cell)) // Negation of predicate --> Not Part of subdomain
+        {
+          for (unsigned int v=0; v<GeometryInfo<MeshType::dimension>::vertices_per_cell; ++v)
+            vertices_outside_subdomain[cell->vertex_index(v)] = true;
+          n_non_predicate_cells++;
+        }
+
+    // If all the active cells conform to the predicate
+    // or if none of the active cells conform to the predicate
+    // there is no active cell layer around the predicate
+    // subdomain (within any distance)
+    if ( n_non_predicate_cells == 0  || n_non_predicate_cells == mesh.get_triangulation().n_active_cells() )
+      return std::vector<typename MeshType::active_cell_iterator>();
+
+    // Find the cells that conform to the predicate
+    // but share a vertex with the cell not in the predicate subdomain
+    for ( typename MeshType::active_cell_iterator
+          cell = mesh.begin_active();
+          cell != mesh.end(); ++cell)
+      if ( predicate(cell)) // True predicate --> Potential boundary cell of the subdomain
+        for (unsigned int v=0; v<GeometryInfo<MeshType::dimension>::vertices_per_cell; ++v)
+          if (vertices_outside_subdomain[cell->vertex_index(v)] == true)
+            {
+              subdomain_boundary_cells.push_back(cell);
+              break; // No need to go through remaining vertices
+            }
+
+    // To cheaply filter out some cells located far away from the predicate subdomain,
+    // get the bounding box of the predicate subdomain.
+    std::pair< Point<spacedim>, Point<spacedim> > bounding_box = compute_bounding_box( mesh,
+        predicate );
+
+    // DOUBLE_EPSILON to compare really close double values
+    const double &DOUBLE_EPSILON = 100.*std::numeric_limits<double>::epsilon();
+
+    // Add layer_thickness to the bounding box
+    for ( unsigned int d=0; d<spacedim; ++d)
+      {
+        bounding_box.first[d]  -= (layer_thickness+DOUBLE_EPSILON); // minp
+        bounding_box.second[d] += (layer_thickness+DOUBLE_EPSILON); // maxp
+      }
+
+    std::vector<Point<spacedim> > subdomain_boundary_cells_centers; // cache all the subdomain boundary cells centers here
+    std::vector<double> subdomain_boundary_cells_radii; // cache all the subdomain boundary cells radii
+    subdomain_boundary_cells_centers.reserve (subdomain_boundary_cells.size());
+    subdomain_boundary_cells_radii.reserve (subdomain_boundary_cells.size());
+    // compute cell radius for each boundary cell of the predicate subdomain
+    for ( typename std::vector<typename MeshType::active_cell_iterator>::const_iterator
+          subdomain_boundary_cell_iterator  = subdomain_boundary_cells.begin();
+          subdomain_boundary_cell_iterator != subdomain_boundary_cells.end(); ++subdomain_boundary_cell_iterator )
+      {
+        const std::pair<Point<spacedim>, double> &
+        subdomain_boundary_cell_enclosing_ball = (*subdomain_boundary_cell_iterator)->enclosing_ball();
+
+        subdomain_boundary_cells_centers.push_back( subdomain_boundary_cell_enclosing_ball.first);
+        subdomain_boundary_cells_radii.push_back( subdomain_boundary_cell_enclosing_ball.second);
+      }
+    AssertThrow( subdomain_boundary_cells_radii.size() == subdomain_boundary_cells_centers.size(),
+                 ExcInternalError());
+
+    // Find the cells that are within layer_thickness of predicate subdomain boundary
+    // distance but are inside the extended bounding box.
+    // Most cells might be outside the extended bounding box, so we could skip them.
+    // Those cells that are inside the extended bounding box but are not part of the
+    // predicate subdomain are possible candidates to be within the distance to the
+    // boundary cells of the predicate subdomain.
+    for ( typename MeshType::active_cell_iterator
+          cell = mesh.begin_active();
+          cell != mesh.end(); ++cell)
+      {
+        // Ignore all the cells that are in the predicate subdomain
+        if ( predicate(cell))
+          continue;
+
+        const std::pair<Point<spacedim>, double> &cell_enclosing_ball
+          = cell->enclosing_ball();
+
+        const Point<spacedim> &cell_enclosing_ball_center = cell_enclosing_ball.first;
+        const double &cell_enclosing_ball_radius = cell_enclosing_ball.second;
+
+        bool cell_inside = true; // reset for each cell
+
+        for (unsigned int d = 0; d < spacedim; ++d)
+          cell_inside &= (cell_enclosing_ball_center[d] + cell_enclosing_ball_radius > bounding_box.first[d])
+                         && (cell_enclosing_ball_center[d] - cell_enclosing_ball_radius < bounding_box.second[d]);
+        // cell_inside is true if its enclosing ball intersects the extended bounding box
+
+        // Ignore all the cells that are outside the extended bounding box
+        if (cell_inside)
+          for (unsigned int i =0; i< subdomain_boundary_cells_radii.size(); ++i)
+            if ( cell_enclosing_ball_center.distance_square(subdomain_boundary_cells_centers[i])
+                 <  Utilities::fixed_power<2>( cell_enclosing_ball_radius +
+                                               subdomain_boundary_cells_radii[i] +
+                                               layer_thickness + DOUBLE_EPSILON ))
+              {
+                active_cell_layer_within_distance.push_back(cell);
+                break; // Exit the loop checking all the remaining subdomain boundary cells
+              }
+
+      }
+    return active_cell_layer_within_distance;
+  }
+
+
+
+  template <class MeshType>
+  std::vector<typename MeshType::active_cell_iterator>
+  compute_ghost_cell_layer_within_distance ( const MeshType &mesh, const double layer_thickness)
+  {
+    IteratorFilters::LocallyOwnedCell locally_owned_cell_predicate;
+    std::function<bool (const typename MeshType::active_cell_iterator &)> predicate (locally_owned_cell_predicate);
+
+    const std::vector<typename MeshType::active_cell_iterator>
+    ghost_cell_layer_within_distance = compute_active_cell_layer_within_distance (mesh, predicate, layer_thickness);
+
+    // Check that we never return locally owned or artificial cells
+    // What is left should only be the ghost cells
+    Assert(contains_locally_owned_cells<MeshType>(ghost_cell_layer_within_distance) == false,
+           ExcMessage("Ghost cells within layer_thickness contains locally owned cells."));
+    Assert(contains_artificial_cells<MeshType>(ghost_cell_layer_within_distance) == false,
+           ExcMessage("Ghost cells within layer_thickness contains artificial cells."
+                      "The function compute_ghost_cell_layer_within_distance "
+                      "is probably called while using parallel::distributed::Triangulation. "
+                      "In such case please refer to the description of this function."));
+
+    return ghost_cell_layer_within_distance;
+  }
+
+
+
+  template< class MeshType>
+  std::pair< Point<MeshType::space_dimension>, Point<MeshType::space_dimension> >
+  compute_bounding_box
+  ( const MeshType                                                                    &mesh,
+    const std::function<bool (const typename MeshType::active_cell_iterator &)> &predicate )
+  {
+    std::vector<bool> locally_active_vertices_on_subdomain (mesh.get_triangulation().n_vertices(),
+                                                            false);
+
+    const unsigned int spacedim = MeshType::space_dimension;
+
+    // Two extreme points can define the bounding box
+    // around the active cells that conform to the given predicate.
+    Point<MeshType::space_dimension> maxp, minp;
+
+    // initialize minp and maxp with the first predicate cell center
+    for ( typename MeshType::active_cell_iterator
+          cell = mesh.begin_active();
+          cell != mesh.end(); ++cell)
+      if ( predicate(cell))
+        {
+          minp = cell->center();
+          maxp = cell->center();
+          break;
+        }
+
+    // Run through all the cells to check if it belongs to predicate domain,
+    // if it belongs to the predicate domain, extend the bounding box.
+    for ( typename MeshType::active_cell_iterator
+          cell = mesh.begin_active();
+          cell != mesh.end(); ++cell)
+      if (predicate(cell)) // True predicate --> Part of subdomain
+        for (unsigned int v=0; v<GeometryInfo<MeshType::dimension>::vertices_per_cell; ++v)
+          if (locally_active_vertices_on_subdomain[cell->vertex_index(v)] == false)
+            {
+              locally_active_vertices_on_subdomain[cell->vertex_index(v)] = true;
+              for ( unsigned int d=0; d<spacedim; ++d)
+                {
+                  minp[d] = std::min( minp[d], cell->vertex(v)[d]);
+                  maxp[d] = std::max( maxp[d], cell->vertex(v)[d]);
+                }
+            }
+
+    return std::make_pair( minp, maxp );
+  }
+
+
+
   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 0c70f351634695440b23b158c3ef62ec3f89b675..4b5c9a308e3c71723b24e3ef703621db8a18600a 100644 (file)
@@ -57,6 +57,24 @@ for (X : TRIANGULATION_AND_DOFHANDLERS; deal_II_dimension : DIMENSIONS ; deal_II
     std::vector<dealii::internal::ActiveCellIterator<deal_II_dimension, deal_II_space_dimension, X>::type>
     compute_ghost_cell_halo_layer (const X &);
 
+
+    template
+    std::vector<dealii::internal::ActiveCellIterator<deal_II_dimension, deal_II_space_dimension, X>::type>
+    compute_active_cell_layer_within_distance (const X &,
+            const std::function<bool (const dealii::internal::ActiveCellIterator<deal_II_dimension, deal_II_space_dimension, X>::type&)> &,
+            const double);
+
+
+    template
+    std::vector<dealii::internal::ActiveCellIterator<deal_II_dimension, deal_II_space_dimension, X>::type>
+    compute_ghost_cell_layer_within_distance (const X &, const double);
+
+
+    template
+    std::pair< Point<X::space_dimension>, Point<X::space_dimension> >
+    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/grid/grid_tools_active_cell_layer_within_distance_01.cc b/tests/grid/grid_tools_active_cell_layer_within_distance_01.cc
new file mode 100644 (file)
index 0000000..abfed26
--- /dev/null
@@ -0,0 +1,130 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2001 - 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.
+//
+// ---------------------------------------------------------------------
+
+
+// This test checks functionality of compute_active_cell_layer_within_distance using
+// a mesh generated by hyper_cube with a layer_thickness of zero.
+// The output is identical to that of grid_tools_halo_layer_01.cc
+// The numerals in the blessed file are exact copy.
+
+
+#include "../tests.h"
+#include <deal.II/base/logstream.h>
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/grid_tools.h>
+#include <deal.II/grid/grid_out.h>
+
+template<int dim>
+bool
+pred_mat_id(const typename Triangulation<dim>::active_cell_iterator &cell)
+{
+  return cell->material_id() == 2;
+}
+
+template <int dim>
+void
+write_mat_id_to_file (const Triangulation<dim> &tria)
+{
+  int count = 0;
+  typename Triangulation<dim>::active_cell_iterator
+  cell = tria.begin_active(),
+  endc = tria.end();
+  for (; cell != endc; ++cell, ++count)
+    {
+      deallog
+          << count << " "
+          << static_cast<int>(cell->material_id())
+          << std::endl;
+    }
+  deallog << std::endl;
+}
+
+
+template <int dim>
+void test ()
+{
+  deallog << "dim = " << dim << std::endl;
+
+  Triangulation<dim> tria;
+  GridGenerator::hyper_cube(tria);
+  tria.refine_global(2);
+
+  typedef typename Triangulation<dim>::active_cell_iterator cell_iterator;
+
+  // Mark a small block at the corner of the hypercube
+  cell_iterator
+  cell = tria.begin_active(),
+  endc = tria.end();
+  for (; cell != endc; ++cell)
+    {
+      bool mark = true;
+      for (unsigned int d=0; d < dim; ++d)
+        if (cell->center()[d] > 0.5)
+          {
+            mark = false;
+            break;
+          }
+
+      if (mark == true)
+        cell->set_material_id(2);
+      else
+        cell->set_material_id(1);
+    }
+
+  deallog << "Grid without skin:" << std::endl;
+  write_mat_id_to_file(tria);
+  // Write to file to visually check result
+  {
+    const std::string filename = "grid_no_halo_" + Utilities::int_to_string(dim) + "d.vtk";
+    std::ofstream f(filename.c_str());
+    GridOut().write_vtk (tria, f);
+  }
+
+  std::function<bool (const cell_iterator &)> predicate = pred_mat_id<dim>;
+
+  // Compute a halo layer around material id 2 and set it to material id 3
+  const std::vector<cell_iterator> active_halo_layer
+    = GridTools::compute_active_cell_layer_within_distance(tria, predicate, 0.); // General predicate
+
+  AssertThrow(active_halo_layer.size() > 0, ExcMessage("No halo layer found."));
+  for (typename std::vector<cell_iterator>::const_iterator
+       it = active_halo_layer.begin();
+       it != active_halo_layer.end(); ++it)
+    {
+      (*it)->set_material_id(3);
+    }
+
+  deallog << "Grid with skin:" << std::endl;
+  write_mat_id_to_file(tria);
+  // Write to file to visually check result
+  {
+    const std::string filename = "grid_with_halo_" + Utilities::int_to_string(dim) + "d.vtk";
+    std::ofstream f(filename.c_str());
+    GridOut().write_vtk (tria, f);
+  }
+}
+
+
+int main ()
+{
+  initlog();
+
+  test<1> ();
+  test<2> ();
+  test<3> ();
+
+  return 0;
+}
diff --git a/tests/grid/grid_tools_active_cell_layer_within_distance_01.output b/tests/grid/grid_tools_active_cell_layer_within_distance_01.output
new file mode 100644 (file)
index 0000000..ce99417
--- /dev/null
@@ -0,0 +1,184 @@
+
+DEAL::dim = 1
+DEAL::Grid without skin:
+DEAL::0 2
+DEAL::1 2
+DEAL::2 1
+DEAL::3 1
+DEAL::
+DEAL::Grid with skin:
+DEAL::0 2
+DEAL::1 2
+DEAL::2 3
+DEAL::3 1
+DEAL::
+DEAL::dim = 2
+DEAL::Grid without skin:
+DEAL::0 2
+DEAL::1 2
+DEAL::2 2
+DEAL::3 2
+DEAL::4 1
+DEAL::5 1
+DEAL::6 1
+DEAL::7 1
+DEAL::8 1
+DEAL::9 1
+DEAL::10 1
+DEAL::11 1
+DEAL::12 1
+DEAL::13 1
+DEAL::14 1
+DEAL::15 1
+DEAL::
+DEAL::Grid with skin:
+DEAL::0 2
+DEAL::1 2
+DEAL::2 2
+DEAL::3 2
+DEAL::4 3
+DEAL::5 1
+DEAL::6 3
+DEAL::7 1
+DEAL::8 3
+DEAL::9 3
+DEAL::10 1
+DEAL::11 1
+DEAL::12 3
+DEAL::13 1
+DEAL::14 1
+DEAL::15 1
+DEAL::
+DEAL::dim = 3
+DEAL::Grid without skin:
+DEAL::0 2
+DEAL::1 2
+DEAL::2 2
+DEAL::3 2
+DEAL::4 2
+DEAL::5 2
+DEAL::6 2
+DEAL::7 2
+DEAL::8 1
+DEAL::9 1
+DEAL::10 1
+DEAL::11 1
+DEAL::12 1
+DEAL::13 1
+DEAL::14 1
+DEAL::15 1
+DEAL::16 1
+DEAL::17 1
+DEAL::18 1
+DEAL::19 1
+DEAL::20 1
+DEAL::21 1
+DEAL::22 1
+DEAL::23 1
+DEAL::24 1
+DEAL::25 1
+DEAL::26 1
+DEAL::27 1
+DEAL::28 1
+DEAL::29 1
+DEAL::30 1
+DEAL::31 1
+DEAL::32 1
+DEAL::33 1
+DEAL::34 1
+DEAL::35 1
+DEAL::36 1
+DEAL::37 1
+DEAL::38 1
+DEAL::39 1
+DEAL::40 1
+DEAL::41 1
+DEAL::42 1
+DEAL::43 1
+DEAL::44 1
+DEAL::45 1
+DEAL::46 1
+DEAL::47 1
+DEAL::48 1
+DEAL::49 1
+DEAL::50 1
+DEAL::51 1
+DEAL::52 1
+DEAL::53 1
+DEAL::54 1
+DEAL::55 1
+DEAL::56 1
+DEAL::57 1
+DEAL::58 1
+DEAL::59 1
+DEAL::60 1
+DEAL::61 1
+DEAL::62 1
+DEAL::63 1
+DEAL::
+DEAL::Grid with skin:
+DEAL::0 2
+DEAL::1 2
+DEAL::2 2
+DEAL::3 2
+DEAL::4 2
+DEAL::5 2
+DEAL::6 2
+DEAL::7 2
+DEAL::8 3
+DEAL::9 1
+DEAL::10 3
+DEAL::11 1
+DEAL::12 3
+DEAL::13 1
+DEAL::14 3
+DEAL::15 1
+DEAL::16 3
+DEAL::17 3
+DEAL::18 1
+DEAL::19 1
+DEAL::20 3
+DEAL::21 3
+DEAL::22 1
+DEAL::23 1
+DEAL::24 3
+DEAL::25 1
+DEAL::26 1
+DEAL::27 1
+DEAL::28 3
+DEAL::29 1
+DEAL::30 1
+DEAL::31 1
+DEAL::32 3
+DEAL::33 3
+DEAL::34 3
+DEAL::35 3
+DEAL::36 1
+DEAL::37 1
+DEAL::38 1
+DEAL::39 1
+DEAL::40 3
+DEAL::41 1
+DEAL::42 3
+DEAL::43 1
+DEAL::44 1
+DEAL::45 1
+DEAL::46 1
+DEAL::47 1
+DEAL::48 3
+DEAL::49 3
+DEAL::50 1
+DEAL::51 1
+DEAL::52 1
+DEAL::53 1
+DEAL::54 1
+DEAL::55 1
+DEAL::56 3
+DEAL::57 1
+DEAL::58 1
+DEAL::59 1
+DEAL::60 1
+DEAL::61 1
+DEAL::62 1
+DEAL::63 1
+DEAL::
diff --git a/tests/grid/grid_tools_active_cell_layer_within_distance_02.cc b/tests/grid/grid_tools_active_cell_layer_within_distance_02.cc
new file mode 100644 (file)
index 0000000..c937957
--- /dev/null
@@ -0,0 +1,170 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2001 - 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.
+//
+// ---------------------------------------------------------------------
+
+// This test checks functionality of compute_active_cell_layer_within_distance using
+// a mesh generated by subdivided_hyper_rectangle. The motivation here is to
+// check predicate skin using a slightly complex mesh than a refined hyper_cube,
+// that could still be visually verified.
+// The output is verified solely by visualizing the resulting active cells
+// within skin distance to predicate subdomain using vtk files
+
+
+// write VTK files for visual inspection
+//#define WRITE_VTK
+
+
+#include "../tests.h"
+#include <deal.II/base/logstream.h>
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/grid_tools.h>
+#include <deal.II/grid/grid_out.h>
+
+#include <vector>
+#include <cmath>
+
+template<int dim>
+bool
+pred_mat_id(const typename Triangulation<dim>::active_cell_iterator &cell)
+{
+  return cell->material_id() == 2;
+}
+
+template <int dim>
+void
+write_mat_id_to_file (const Triangulation<dim> &tria)
+{
+  int count = 0;
+  typename Triangulation<dim>::active_cell_iterator
+  cell = tria.begin_active(),
+  endc = tria.end();
+  for (; cell != endc; ++cell, ++count)
+    {
+      deallog
+          << count << " "
+          << static_cast<int>(cell->material_id())
+          << std::endl;
+    }
+  deallog << std::endl;
+}
+
+
+template <int dim>
+void test ()
+{
+  deallog << "dim = " << dim << std::endl;
+
+
+  std::vector<double> step_sizes_i; // step sizes in i direction (similar in all directions)
+  std::vector< std::vector< double> > step_sizes; // step sizes as input to the subdivided_hyper_rectangle
+
+
+  unsigned int n_steps = 25;
+  double size =0.; // size of domain in i direction
+  for ( unsigned int j = 1; j < n_steps; ++j)
+    {
+      step_sizes_i.push_back( (1.+double(j))/ double(n_steps) );
+      size += (1.+double(j))/ double(n_steps);
+    }
+
+  for (unsigned int d=0; d < dim; ++d)
+    step_sizes.push_back( step_sizes_i );
+
+  const Point<dim> bottom_left;
+  const Point<dim> upper_right = dim == 1 ?
+                                 Point<dim> (size):
+                                 dim == 2 ?
+                                 Point<dim> (size, size):
+                                 Point<dim> (size, size, size);
+
+  Triangulation<dim> tria;
+  GridGenerator::subdivided_hyper_rectangle(tria,
+                                            step_sizes,
+                                            bottom_left,
+                                            upper_right,
+                                            true);
+
+  typedef typename Triangulation<dim>::active_cell_iterator cell_iterator;
+
+  // Mark a small block at the corner of the hypercube
+  cell_iterator
+  cell = tria.begin_active(),
+  endc = tria.end();
+  for (; cell != endc; ++cell)
+    {
+      bool mark = true;
+      for (unsigned int d=0; d < dim; ++d)
+        mark &= cell->center()[d] > 5. && cell->center()[d] < 9.;
+
+      if (mark == true)
+        cell->set_material_id(2);
+      else
+        cell->set_material_id(1);
+
+    }
+
+  deallog << "Grid without skin:" << std::endl;
+  write_mat_id_to_file(tria);
+#ifdef WRITE_VTK
+  // Write to file to visually check result
+  {
+    const std::string filename = "grid_no_skin_" + Utilities::int_to_string(dim) + "d.vtk";
+    std::ofstream f(filename.c_str());
+    GridOut().write_vtk (tria, f);
+  }
+#endif
+
+  std::function<bool (const cell_iterator &)> predicate = pred_mat_id<dim>;
+
+  // Compute a halo layer around material id 2 and set it to material id 3
+  const std::vector<cell_iterator> cells_within_skin
+    = GridTools::compute_active_cell_layer_within_distance(tria, predicate, 1.25); // General predicate
+
+  // This test should fail if skin_thickness is 1./3. (which will accumulate an extra layer of cells)
+
+  AssertThrow(cells_within_skin.size() > 0, ExcMessage("No skin cells found."));
+  for (typename std::vector<cell_iterator>::const_iterator
+       it = cells_within_skin.begin();
+       it != cells_within_skin.end(); ++it)
+    {
+      (*it)->set_material_id(3);
+    }
+
+  deallog << "Grid with skin:" << std::endl;
+  write_mat_id_to_file(tria);
+
+#ifdef WRITE_VTK
+  // Write to file to visually check result
+  {
+    const std::string filename = "grid_with_skin_" + Utilities::int_to_string(dim) + "d.vtk";
+    std::ofstream f(filename.c_str());
+    GridOut().write_vtk (tria, f);
+  }
+#endif
+
+}
+
+
+int main ()
+{
+  initlog();
+
+  test<1> ();
+  test<2> ();    // visually checked for 2D case
+  // test<3> (); // the output is too long to include in the blessed file
+  // visually checked for 3D case
+
+  return 0;
+}
diff --git a/tests/grid/grid_tools_active_cell_layer_within_distance_02.output b/tests/grid/grid_tools_active_cell_layer_within_distance_02.output
new file mode 100644 (file)
index 0000000..18f7879
--- /dev/null
@@ -0,0 +1,1211 @@
+
+DEAL::dim = 1
+DEAL::Grid without skin:
+DEAL::0 1
+DEAL::1 1
+DEAL::2 1
+DEAL::3 1
+DEAL::4 1
+DEAL::5 1
+DEAL::6 1
+DEAL::7 1
+DEAL::8 1
+DEAL::9 1
+DEAL::10 1
+DEAL::11 1
+DEAL::12 1
+DEAL::13 1
+DEAL::14 2
+DEAL::15 2
+DEAL::16 2
+DEAL::17 2
+DEAL::18 2
+DEAL::19 2
+DEAL::20 1
+DEAL::21 1
+DEAL::22 1
+DEAL::23 1
+DEAL::
+DEAL::Grid with skin:
+DEAL::0 1
+DEAL::1 1
+DEAL::2 1
+DEAL::3 1
+DEAL::4 1
+DEAL::5 1
+DEAL::6 1
+DEAL::7 1
+DEAL::8 1
+DEAL::9 1
+DEAL::10 1
+DEAL::11 3
+DEAL::12 3
+DEAL::13 3
+DEAL::14 2
+DEAL::15 2
+DEAL::16 2
+DEAL::17 2
+DEAL::18 2
+DEAL::19 2
+DEAL::20 3
+DEAL::21 3
+DEAL::22 1
+DEAL::23 1
+DEAL::
+DEAL::dim = 2
+DEAL::Grid without skin:
+DEAL::0 1
+DEAL::1 1
+DEAL::2 1
+DEAL::3 1
+DEAL::4 1
+DEAL::5 1
+DEAL::6 1
+DEAL::7 1
+DEAL::8 1
+DEAL::9 1
+DEAL::10 1
+DEAL::11 1
+DEAL::12 1
+DEAL::13 1
+DEAL::14 1
+DEAL::15 1
+DEAL::16 1
+DEAL::17 1
+DEAL::18 1
+DEAL::19 1
+DEAL::20 1
+DEAL::21 1
+DEAL::22 1
+DEAL::23 1
+DEAL::24 1
+DEAL::25 1
+DEAL::26 1
+DEAL::27 1
+DEAL::28 1
+DEAL::29 1
+DEAL::30 1
+DEAL::31 1
+DEAL::32 1
+DEAL::33 1
+DEAL::34 1
+DEAL::35 1
+DEAL::36 1
+DEAL::37 1
+DEAL::38 1
+DEAL::39 1
+DEAL::40 1
+DEAL::41 1
+DEAL::42 1
+DEAL::43 1
+DEAL::44 1
+DEAL::45 1
+DEAL::46 1
+DEAL::47 1
+DEAL::48 1
+DEAL::49 1
+DEAL::50 1
+DEAL::51 1
+DEAL::52 1
+DEAL::53 1
+DEAL::54 1
+DEAL::55 1
+DEAL::56 1
+DEAL::57 1
+DEAL::58 1
+DEAL::59 1
+DEAL::60 1
+DEAL::61 1
+DEAL::62 1
+DEAL::63 1
+DEAL::64 1
+DEAL::65 1
+DEAL::66 1
+DEAL::67 1
+DEAL::68 1
+DEAL::69 1
+DEAL::70 1
+DEAL::71 1
+DEAL::72 1
+DEAL::73 1
+DEAL::74 1
+DEAL::75 1
+DEAL::76 1
+DEAL::77 1
+DEAL::78 1
+DEAL::79 1
+DEAL::80 1
+DEAL::81 1
+DEAL::82 1
+DEAL::83 1
+DEAL::84 1
+DEAL::85 1
+DEAL::86 1
+DEAL::87 1
+DEAL::88 1
+DEAL::89 1
+DEAL::90 1
+DEAL::91 1
+DEAL::92 1
+DEAL::93 1
+DEAL::94 1
+DEAL::95 1
+DEAL::96 1
+DEAL::97 1
+DEAL::98 1
+DEAL::99 1
+DEAL::100 1
+DEAL::101 1
+DEAL::102 1
+DEAL::103 1
+DEAL::104 1
+DEAL::105 1
+DEAL::106 1
+DEAL::107 1
+DEAL::108 1
+DEAL::109 1
+DEAL::110 1
+DEAL::111 1
+DEAL::112 1
+DEAL::113 1
+DEAL::114 1
+DEAL::115 1
+DEAL::116 1
+DEAL::117 1
+DEAL::118 1
+DEAL::119 1
+DEAL::120 1
+DEAL::121 1
+DEAL::122 1
+DEAL::123 1
+DEAL::124 1
+DEAL::125 1
+DEAL::126 1
+DEAL::127 1
+DEAL::128 1
+DEAL::129 1
+DEAL::130 1
+DEAL::131 1
+DEAL::132 1
+DEAL::133 1
+DEAL::134 1
+DEAL::135 1
+DEAL::136 1
+DEAL::137 1
+DEAL::138 1
+DEAL::139 1
+DEAL::140 1
+DEAL::141 1
+DEAL::142 1
+DEAL::143 1
+DEAL::144 1
+DEAL::145 1
+DEAL::146 1
+DEAL::147 1
+DEAL::148 1
+DEAL::149 1
+DEAL::150 1
+DEAL::151 1
+DEAL::152 1
+DEAL::153 1
+DEAL::154 1
+DEAL::155 1
+DEAL::156 1
+DEAL::157 1
+DEAL::158 1
+DEAL::159 1
+DEAL::160 1
+DEAL::161 1
+DEAL::162 1
+DEAL::163 1
+DEAL::164 1
+DEAL::165 1
+DEAL::166 1
+DEAL::167 1
+DEAL::168 1
+DEAL::169 1
+DEAL::170 1
+DEAL::171 1
+DEAL::172 1
+DEAL::173 1
+DEAL::174 1
+DEAL::175 1
+DEAL::176 1
+DEAL::177 1
+DEAL::178 1
+DEAL::179 1
+DEAL::180 1
+DEAL::181 1
+DEAL::182 1
+DEAL::183 1
+DEAL::184 1
+DEAL::185 1
+DEAL::186 1
+DEAL::187 1
+DEAL::188 1
+DEAL::189 1
+DEAL::190 1
+DEAL::191 1
+DEAL::192 1
+DEAL::193 1
+DEAL::194 1
+DEAL::195 1
+DEAL::196 1
+DEAL::197 1
+DEAL::198 1
+DEAL::199 1
+DEAL::200 1
+DEAL::201 1
+DEAL::202 1
+DEAL::203 1
+DEAL::204 1
+DEAL::205 1
+DEAL::206 1
+DEAL::207 1
+DEAL::208 1
+DEAL::209 1
+DEAL::210 1
+DEAL::211 1
+DEAL::212 1
+DEAL::213 1
+DEAL::214 1
+DEAL::215 1
+DEAL::216 1
+DEAL::217 1
+DEAL::218 1
+DEAL::219 1
+DEAL::220 1
+DEAL::221 1
+DEAL::222 1
+DEAL::223 1
+DEAL::224 1
+DEAL::225 1
+DEAL::226 1
+DEAL::227 1
+DEAL::228 1
+DEAL::229 1
+DEAL::230 1
+DEAL::231 1
+DEAL::232 1
+DEAL::233 1
+DEAL::234 1
+DEAL::235 1
+DEAL::236 1
+DEAL::237 1
+DEAL::238 1
+DEAL::239 1
+DEAL::240 1
+DEAL::241 1
+DEAL::242 1
+DEAL::243 1
+DEAL::244 1
+DEAL::245 1
+DEAL::246 1
+DEAL::247 1
+DEAL::248 1
+DEAL::249 1
+DEAL::250 1
+DEAL::251 1
+DEAL::252 1
+DEAL::253 1
+DEAL::254 1
+DEAL::255 1
+DEAL::256 1
+DEAL::257 1
+DEAL::258 1
+DEAL::259 1
+DEAL::260 1
+DEAL::261 1
+DEAL::262 1
+DEAL::263 1
+DEAL::264 1
+DEAL::265 1
+DEAL::266 1
+DEAL::267 1
+DEAL::268 1
+DEAL::269 1
+DEAL::270 1
+DEAL::271 1
+DEAL::272 1
+DEAL::273 1
+DEAL::274 1
+DEAL::275 1
+DEAL::276 1
+DEAL::277 1
+DEAL::278 1
+DEAL::279 1
+DEAL::280 1
+DEAL::281 1
+DEAL::282 1
+DEAL::283 1
+DEAL::284 1
+DEAL::285 1
+DEAL::286 1
+DEAL::287 1
+DEAL::288 1
+DEAL::289 1
+DEAL::290 1
+DEAL::291 1
+DEAL::292 1
+DEAL::293 1
+DEAL::294 1
+DEAL::295 1
+DEAL::296 1
+DEAL::297 1
+DEAL::298 1
+DEAL::299 1
+DEAL::300 1
+DEAL::301 1
+DEAL::302 1
+DEAL::303 1
+DEAL::304 1
+DEAL::305 1
+DEAL::306 1
+DEAL::307 1
+DEAL::308 1
+DEAL::309 1
+DEAL::310 1
+DEAL::311 1
+DEAL::312 1
+DEAL::313 1
+DEAL::314 1
+DEAL::315 1
+DEAL::316 1
+DEAL::317 1
+DEAL::318 1
+DEAL::319 1
+DEAL::320 1
+DEAL::321 1
+DEAL::322 1
+DEAL::323 1
+DEAL::324 1
+DEAL::325 1
+DEAL::326 1
+DEAL::327 1
+DEAL::328 1
+DEAL::329 1
+DEAL::330 1
+DEAL::331 1
+DEAL::332 1
+DEAL::333 1
+DEAL::334 1
+DEAL::335 1
+DEAL::336 1
+DEAL::337 1
+DEAL::338 1
+DEAL::339 1
+DEAL::340 1
+DEAL::341 1
+DEAL::342 1
+DEAL::343 1
+DEAL::344 1
+DEAL::345 1
+DEAL::346 1
+DEAL::347 1
+DEAL::348 1
+DEAL::349 1
+DEAL::350 2
+DEAL::351 2
+DEAL::352 2
+DEAL::353 2
+DEAL::354 2
+DEAL::355 2
+DEAL::356 1
+DEAL::357 1
+DEAL::358 1
+DEAL::359 1
+DEAL::360 1
+DEAL::361 1
+DEAL::362 1
+DEAL::363 1
+DEAL::364 1
+DEAL::365 1
+DEAL::366 1
+DEAL::367 1
+DEAL::368 1
+DEAL::369 1
+DEAL::370 1
+DEAL::371 1
+DEAL::372 1
+DEAL::373 1
+DEAL::374 2
+DEAL::375 2
+DEAL::376 2
+DEAL::377 2
+DEAL::378 2
+DEAL::379 2
+DEAL::380 1
+DEAL::381 1
+DEAL::382 1
+DEAL::383 1
+DEAL::384 1
+DEAL::385 1
+DEAL::386 1
+DEAL::387 1
+DEAL::388 1
+DEAL::389 1
+DEAL::390 1
+DEAL::391 1
+DEAL::392 1
+DEAL::393 1
+DEAL::394 1
+DEAL::395 1
+DEAL::396 1
+DEAL::397 1
+DEAL::398 2
+DEAL::399 2
+DEAL::400 2
+DEAL::401 2
+DEAL::402 2
+DEAL::403 2
+DEAL::404 1
+DEAL::405 1
+DEAL::406 1
+DEAL::407 1
+DEAL::408 1
+DEAL::409 1
+DEAL::410 1
+DEAL::411 1
+DEAL::412 1
+DEAL::413 1
+DEAL::414 1
+DEAL::415 1
+DEAL::416 1
+DEAL::417 1
+DEAL::418 1
+DEAL::419 1
+DEAL::420 1
+DEAL::421 1
+DEAL::422 2
+DEAL::423 2
+DEAL::424 2
+DEAL::425 2
+DEAL::426 2
+DEAL::427 2
+DEAL::428 1
+DEAL::429 1
+DEAL::430 1
+DEAL::431 1
+DEAL::432 1
+DEAL::433 1
+DEAL::434 1
+DEAL::435 1
+DEAL::436 1
+DEAL::437 1
+DEAL::438 1
+DEAL::439 1
+DEAL::440 1
+DEAL::441 1
+DEAL::442 1
+DEAL::443 1
+DEAL::444 1
+DEAL::445 1
+DEAL::446 2
+DEAL::447 2
+DEAL::448 2
+DEAL::449 2
+DEAL::450 2
+DEAL::451 2
+DEAL::452 1
+DEAL::453 1
+DEAL::454 1
+DEAL::455 1
+DEAL::456 1
+DEAL::457 1
+DEAL::458 1
+DEAL::459 1
+DEAL::460 1
+DEAL::461 1
+DEAL::462 1
+DEAL::463 1
+DEAL::464 1
+DEAL::465 1
+DEAL::466 1
+DEAL::467 1
+DEAL::468 1
+DEAL::469 1
+DEAL::470 2
+DEAL::471 2
+DEAL::472 2
+DEAL::473 2
+DEAL::474 2
+DEAL::475 2
+DEAL::476 1
+DEAL::477 1
+DEAL::478 1
+DEAL::479 1
+DEAL::480 1
+DEAL::481 1
+DEAL::482 1
+DEAL::483 1
+DEAL::484 1
+DEAL::485 1
+DEAL::486 1
+DEAL::487 1
+DEAL::488 1
+DEAL::489 1
+DEAL::490 1
+DEAL::491 1
+DEAL::492 1
+DEAL::493 1
+DEAL::494 1
+DEAL::495 1
+DEAL::496 1
+DEAL::497 1
+DEAL::498 1
+DEAL::499 1
+DEAL::500 1
+DEAL::501 1
+DEAL::502 1
+DEAL::503 1
+DEAL::504 1
+DEAL::505 1
+DEAL::506 1
+DEAL::507 1
+DEAL::508 1
+DEAL::509 1
+DEAL::510 1
+DEAL::511 1
+DEAL::512 1
+DEAL::513 1
+DEAL::514 1
+DEAL::515 1
+DEAL::516 1
+DEAL::517 1
+DEAL::518 1
+DEAL::519 1
+DEAL::520 1
+DEAL::521 1
+DEAL::522 1
+DEAL::523 1
+DEAL::524 1
+DEAL::525 1
+DEAL::526 1
+DEAL::527 1
+DEAL::528 1
+DEAL::529 1
+DEAL::530 1
+DEAL::531 1
+DEAL::532 1
+DEAL::533 1
+DEAL::534 1
+DEAL::535 1
+DEAL::536 1
+DEAL::537 1
+DEAL::538 1
+DEAL::539 1
+DEAL::540 1
+DEAL::541 1
+DEAL::542 1
+DEAL::543 1
+DEAL::544 1
+DEAL::545 1
+DEAL::546 1
+DEAL::547 1
+DEAL::548 1
+DEAL::549 1
+DEAL::550 1
+DEAL::551 1
+DEAL::552 1
+DEAL::553 1
+DEAL::554 1
+DEAL::555 1
+DEAL::556 1
+DEAL::557 1
+DEAL::558 1
+DEAL::559 1
+DEAL::560 1
+DEAL::561 1
+DEAL::562 1
+DEAL::563 1
+DEAL::564 1
+DEAL::565 1
+DEAL::566 1
+DEAL::567 1
+DEAL::568 1
+DEAL::569 1
+DEAL::570 1
+DEAL::571 1
+DEAL::572 1
+DEAL::573 1
+DEAL::574 1
+DEAL::575 1
+DEAL::
+DEAL::Grid with skin:
+DEAL::0 1
+DEAL::1 1
+DEAL::2 1
+DEAL::3 1
+DEAL::4 1
+DEAL::5 1
+DEAL::6 1
+DEAL::7 1
+DEAL::8 1
+DEAL::9 1
+DEAL::10 1
+DEAL::11 1
+DEAL::12 1
+DEAL::13 1
+DEAL::14 1
+DEAL::15 1
+DEAL::16 1
+DEAL::17 1
+DEAL::18 1
+DEAL::19 1
+DEAL::20 1
+DEAL::21 1
+DEAL::22 1
+DEAL::23 1
+DEAL::24 1
+DEAL::25 1
+DEAL::26 1
+DEAL::27 1
+DEAL::28 1
+DEAL::29 1
+DEAL::30 1
+DEAL::31 1
+DEAL::32 1
+DEAL::33 1
+DEAL::34 1
+DEAL::35 1
+DEAL::36 1
+DEAL::37 1
+DEAL::38 1
+DEAL::39 1
+DEAL::40 1
+DEAL::41 1
+DEAL::42 1
+DEAL::43 1
+DEAL::44 1
+DEAL::45 1
+DEAL::46 1
+DEAL::47 1
+DEAL::48 1
+DEAL::49 1
+DEAL::50 1
+DEAL::51 1
+DEAL::52 1
+DEAL::53 1
+DEAL::54 1
+DEAL::55 1
+DEAL::56 1
+DEAL::57 1
+DEAL::58 1
+DEAL::59 1
+DEAL::60 1
+DEAL::61 1
+DEAL::62 1
+DEAL::63 1
+DEAL::64 1
+DEAL::65 1
+DEAL::66 1
+DEAL::67 1
+DEAL::68 1
+DEAL::69 1
+DEAL::70 1
+DEAL::71 1
+DEAL::72 1
+DEAL::73 1
+DEAL::74 1
+DEAL::75 1
+DEAL::76 1
+DEAL::77 1
+DEAL::78 1
+DEAL::79 1
+DEAL::80 1
+DEAL::81 1
+DEAL::82 1
+DEAL::83 1
+DEAL::84 1
+DEAL::85 1
+DEAL::86 1
+DEAL::87 1
+DEAL::88 1
+DEAL::89 1
+DEAL::90 1
+DEAL::91 1
+DEAL::92 1
+DEAL::93 1
+DEAL::94 1
+DEAL::95 1
+DEAL::96 1
+DEAL::97 1
+DEAL::98 1
+DEAL::99 1
+DEAL::100 1
+DEAL::101 1
+DEAL::102 1
+DEAL::103 1
+DEAL::104 1
+DEAL::105 1
+DEAL::106 1
+DEAL::107 1
+DEAL::108 1
+DEAL::109 1
+DEAL::110 1
+DEAL::111 1
+DEAL::112 1
+DEAL::113 1
+DEAL::114 1
+DEAL::115 1
+DEAL::116 1
+DEAL::117 1
+DEAL::118 1
+DEAL::119 1
+DEAL::120 1
+DEAL::121 1
+DEAL::122 1
+DEAL::123 1
+DEAL::124 1
+DEAL::125 1
+DEAL::126 1
+DEAL::127 1
+DEAL::128 1
+DEAL::129 1
+DEAL::130 1
+DEAL::131 1
+DEAL::132 1
+DEAL::133 1
+DEAL::134 1
+DEAL::135 1
+DEAL::136 1
+DEAL::137 1
+DEAL::138 1
+DEAL::139 1
+DEAL::140 1
+DEAL::141 1
+DEAL::142 1
+DEAL::143 1
+DEAL::144 1
+DEAL::145 1
+DEAL::146 1
+DEAL::147 1
+DEAL::148 1
+DEAL::149 1
+DEAL::150 1
+DEAL::151 1
+DEAL::152 1
+DEAL::153 1
+DEAL::154 1
+DEAL::155 1
+DEAL::156 1
+DEAL::157 1
+DEAL::158 1
+DEAL::159 1
+DEAL::160 1
+DEAL::161 1
+DEAL::162 1
+DEAL::163 1
+DEAL::164 1
+DEAL::165 1
+DEAL::166 1
+DEAL::167 1
+DEAL::168 1
+DEAL::169 1
+DEAL::170 1
+DEAL::171 1
+DEAL::172 1
+DEAL::173 1
+DEAL::174 1
+DEAL::175 1
+DEAL::176 1
+DEAL::177 1
+DEAL::178 1
+DEAL::179 1
+DEAL::180 1
+DEAL::181 1
+DEAL::182 1
+DEAL::183 1
+DEAL::184 1
+DEAL::185 1
+DEAL::186 1
+DEAL::187 1
+DEAL::188 1
+DEAL::189 1
+DEAL::190 1
+DEAL::191 1
+DEAL::192 1
+DEAL::193 1
+DEAL::194 1
+DEAL::195 1
+DEAL::196 1
+DEAL::197 1
+DEAL::198 1
+DEAL::199 1
+DEAL::200 1
+DEAL::201 1
+DEAL::202 1
+DEAL::203 1
+DEAL::204 1
+DEAL::205 1
+DEAL::206 1
+DEAL::207 1
+DEAL::208 1
+DEAL::209 1
+DEAL::210 1
+DEAL::211 1
+DEAL::212 1
+DEAL::213 1
+DEAL::214 1
+DEAL::215 1
+DEAL::216 1
+DEAL::217 1
+DEAL::218 1
+DEAL::219 1
+DEAL::220 1
+DEAL::221 1
+DEAL::222 1
+DEAL::223 1
+DEAL::224 1
+DEAL::225 1
+DEAL::226 1
+DEAL::227 1
+DEAL::228 1
+DEAL::229 1
+DEAL::230 1
+DEAL::231 1
+DEAL::232 1
+DEAL::233 1
+DEAL::234 1
+DEAL::235 1
+DEAL::236 1
+DEAL::237 1
+DEAL::238 1
+DEAL::239 1
+DEAL::240 1
+DEAL::241 1
+DEAL::242 1
+DEAL::243 1
+DEAL::244 1
+DEAL::245 1
+DEAL::246 1
+DEAL::247 1
+DEAL::248 1
+DEAL::249 1
+DEAL::250 1
+DEAL::251 1
+DEAL::252 1
+DEAL::253 1
+DEAL::254 1
+DEAL::255 1
+DEAL::256 1
+DEAL::257 1
+DEAL::258 1
+DEAL::259 1
+DEAL::260 1
+DEAL::261 1
+DEAL::262 1
+DEAL::263 1
+DEAL::264 1
+DEAL::265 1
+DEAL::266 1
+DEAL::267 1
+DEAL::268 1
+DEAL::269 1
+DEAL::270 1
+DEAL::271 1
+DEAL::272 1
+DEAL::273 1
+DEAL::274 1
+DEAL::275 1
+DEAL::276 1
+DEAL::277 3
+DEAL::278 3
+DEAL::279 3
+DEAL::280 3
+DEAL::281 3
+DEAL::282 3
+DEAL::283 3
+DEAL::284 3
+DEAL::285 1
+DEAL::286 1
+DEAL::287 1
+DEAL::288 1
+DEAL::289 1
+DEAL::290 1
+DEAL::291 1
+DEAL::292 1
+DEAL::293 1
+DEAL::294 1
+DEAL::295 1
+DEAL::296 1
+DEAL::297 1
+DEAL::298 1
+DEAL::299 1
+DEAL::300 3
+DEAL::301 3
+DEAL::302 3
+DEAL::303 3
+DEAL::304 3
+DEAL::305 3
+DEAL::306 3
+DEAL::307 3
+DEAL::308 3
+DEAL::309 3
+DEAL::310 1
+DEAL::311 1
+DEAL::312 1
+DEAL::313 1
+DEAL::314 1
+DEAL::315 1
+DEAL::316 1
+DEAL::317 1
+DEAL::318 1
+DEAL::319 1
+DEAL::320 1
+DEAL::321 1
+DEAL::322 1
+DEAL::323 3
+DEAL::324 3
+DEAL::325 3
+DEAL::326 3
+DEAL::327 3
+DEAL::328 3
+DEAL::329 3
+DEAL::330 3
+DEAL::331 3
+DEAL::332 3
+DEAL::333 3
+DEAL::334 1
+DEAL::335 1
+DEAL::336 1
+DEAL::337 1
+DEAL::338 1
+DEAL::339 1
+DEAL::340 1
+DEAL::341 1
+DEAL::342 1
+DEAL::343 1
+DEAL::344 1
+DEAL::345 1
+DEAL::346 1
+DEAL::347 3
+DEAL::348 3
+DEAL::349 3
+DEAL::350 2
+DEAL::351 2
+DEAL::352 2
+DEAL::353 2
+DEAL::354 2
+DEAL::355 2
+DEAL::356 3
+DEAL::357 3
+DEAL::358 1
+DEAL::359 1
+DEAL::360 1
+DEAL::361 1
+DEAL::362 1
+DEAL::363 1
+DEAL::364 1
+DEAL::365 1
+DEAL::366 1
+DEAL::367 1
+DEAL::368 1
+DEAL::369 1
+DEAL::370 1
+DEAL::371 3
+DEAL::372 3
+DEAL::373 3
+DEAL::374 2
+DEAL::375 2
+DEAL::376 2
+DEAL::377 2
+DEAL::378 2
+DEAL::379 2
+DEAL::380 3
+DEAL::381 3
+DEAL::382 1
+DEAL::383 1
+DEAL::384 1
+DEAL::385 1
+DEAL::386 1
+DEAL::387 1
+DEAL::388 1
+DEAL::389 1
+DEAL::390 1
+DEAL::391 1
+DEAL::392 1
+DEAL::393 1
+DEAL::394 1
+DEAL::395 3
+DEAL::396 3
+DEAL::397 3
+DEAL::398 2
+DEAL::399 2
+DEAL::400 2
+DEAL::401 2
+DEAL::402 2
+DEAL::403 2
+DEAL::404 3
+DEAL::405 3
+DEAL::406 1
+DEAL::407 1
+DEAL::408 1
+DEAL::409 1
+DEAL::410 1
+DEAL::411 1
+DEAL::412 1
+DEAL::413 1
+DEAL::414 1
+DEAL::415 1
+DEAL::416 1
+DEAL::417 1
+DEAL::418 1
+DEAL::419 3
+DEAL::420 3
+DEAL::421 3
+DEAL::422 2
+DEAL::423 2
+DEAL::424 2
+DEAL::425 2
+DEAL::426 2
+DEAL::427 2
+DEAL::428 3
+DEAL::429 3
+DEAL::430 1
+DEAL::431 1
+DEAL::432 1
+DEAL::433 1
+DEAL::434 1
+DEAL::435 1
+DEAL::436 1
+DEAL::437 1
+DEAL::438 1
+DEAL::439 1
+DEAL::440 1
+DEAL::441 1
+DEAL::442 1
+DEAL::443 3
+DEAL::444 3
+DEAL::445 3
+DEAL::446 2
+DEAL::447 2
+DEAL::448 2
+DEAL::449 2
+DEAL::450 2
+DEAL::451 2
+DEAL::452 3
+DEAL::453 3
+DEAL::454 1
+DEAL::455 1
+DEAL::456 1
+DEAL::457 1
+DEAL::458 1
+DEAL::459 1
+DEAL::460 1
+DEAL::461 1
+DEAL::462 1
+DEAL::463 1
+DEAL::464 1
+DEAL::465 1
+DEAL::466 1
+DEAL::467 3
+DEAL::468 3
+DEAL::469 3
+DEAL::470 2
+DEAL::471 2
+DEAL::472 2
+DEAL::473 2
+DEAL::474 2
+DEAL::475 2
+DEAL::476 3
+DEAL::477 3
+DEAL::478 1
+DEAL::479 1
+DEAL::480 1
+DEAL::481 1
+DEAL::482 1
+DEAL::483 1
+DEAL::484 1
+DEAL::485 1
+DEAL::486 1
+DEAL::487 1
+DEAL::488 1
+DEAL::489 1
+DEAL::490 1
+DEAL::491 3
+DEAL::492 3
+DEAL::493 3
+DEAL::494 3
+DEAL::495 3
+DEAL::496 3
+DEAL::497 3
+DEAL::498 3
+DEAL::499 3
+DEAL::500 3
+DEAL::501 3
+DEAL::502 1
+DEAL::503 1
+DEAL::504 1
+DEAL::505 1
+DEAL::506 1
+DEAL::507 1
+DEAL::508 1
+DEAL::509 1
+DEAL::510 1
+DEAL::511 1
+DEAL::512 1
+DEAL::513 1
+DEAL::514 1
+DEAL::515 1
+DEAL::516 3
+DEAL::517 3
+DEAL::518 3
+DEAL::519 3
+DEAL::520 3
+DEAL::521 3
+DEAL::522 3
+DEAL::523 3
+DEAL::524 3
+DEAL::525 3
+DEAL::526 1
+DEAL::527 1
+DEAL::528 1
+DEAL::529 1
+DEAL::530 1
+DEAL::531 1
+DEAL::532 1
+DEAL::533 1
+DEAL::534 1
+DEAL::535 1
+DEAL::536 1
+DEAL::537 1
+DEAL::538 1
+DEAL::539 1
+DEAL::540 1
+DEAL::541 1
+DEAL::542 1
+DEAL::543 1
+DEAL::544 1
+DEAL::545 1
+DEAL::546 1
+DEAL::547 1
+DEAL::548 1
+DEAL::549 1
+DEAL::550 1
+DEAL::551 1
+DEAL::552 1
+DEAL::553 1
+DEAL::554 1
+DEAL::555 1
+DEAL::556 1
+DEAL::557 1
+DEAL::558 1
+DEAL::559 1
+DEAL::560 1
+DEAL::561 1
+DEAL::562 1
+DEAL::563 1
+DEAL::564 1
+DEAL::565 1
+DEAL::566 1
+DEAL::567 1
+DEAL::568 1
+DEAL::569 1
+DEAL::570 1
+DEAL::571 1
+DEAL::572 1
+DEAL::573 1
+DEAL::574 1
+DEAL::575 1
+DEAL::
diff --git a/tests/grid/grid_tools_bounding_box_01.cc b/tests/grid/grid_tools_bounding_box_01.cc
new file mode 100644 (file)
index 0000000..3ea31e6
--- /dev/null
@@ -0,0 +1,83 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2001 - 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.
+//
+// ---------------------------------------------------------------------
+
+// Short test to validate compute_bounding_box()
+
+#include "../tests.h"
+#include <deal.II/base/logstream.h>
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/grid_tools.h>
+#include <deal.II/grid/grid_out.h>
+
+template<int dim>
+bool
+pred_mat_id(const typename Triangulation<dim>::active_cell_iterator &cell)
+{
+  return cell->material_id() == 2;
+}
+
+template <int dim>
+void test ()
+{
+  deallog << "dim = " << dim << std::endl;
+
+  Triangulation<dim> tria;
+  GridGenerator::hyper_cube(tria);
+  tria.refine_global(2);
+
+  typedef typename Triangulation<dim>::active_cell_iterator cell_iterator;
+
+  // Mark a small block at the corner of the hypercube
+  cell_iterator
+  cell = tria.begin_active(),
+  endc = tria.end();
+  for (; cell != endc; ++cell)
+    {
+      bool mark = true;
+      for (unsigned int d=0; d < dim; ++d)
+        if (cell->center()[d] > 0.33 )
+          {
+            mark = false;
+            break;
+          }
+
+      if (mark == true)
+        cell->set_material_id(2);
+      else
+        cell->set_material_id(1);
+    }
+
+  std::function<bool (const cell_iterator &)> predicate = pred_mat_id<dim>;
+
+  // Find bounding box that surrounds cells with material id 2
+  auto bounding_box
+    = GridTools::compute_bounding_box(tria, predicate); // General predicate
+
+  deallog << bounding_box.first << " " << bounding_box.second << std::endl;
+
+}
+
+
+int main ()
+{
+  initlog();
+
+  test<1> ();
+  test<2> ();
+  test<3> ();
+
+  return 0;
+}
diff --git a/tests/grid/grid_tools_bounding_box_01.output b/tests/grid/grid_tools_bounding_box_01.output
new file mode 100644 (file)
index 0000000..e045db7
--- /dev/null
@@ -0,0 +1,7 @@
+
+DEAL::dim = 1
+DEAL::0.00000 0.250000
+DEAL::dim = 2
+DEAL::0.00000 0.00000 0.250000 0.250000
+DEAL::dim = 3
+DEAL::0.00000 0.00000 0.00000 0.250000 0.250000 0.250000

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.