]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Implementation of general cell halo layer function inside GridTools. 1396/head
authorJean-Paul Pelteret <jppelteret@gmail.com>
Thu, 20 Aug 2015 15:24:43 +0000 (17:24 +0200)
committerJean-Paul Pelteret <jppelteret@gmail.com>
Wed, 26 Aug 2015 09:28:52 +0000 (11:28 +0200)
Added filtered iterators that work on material id and active FE index, with the
option of only extracting locally owned cells.
Introduced a GridTools function to extract the halo layer that is composed
of a subset of ghost cells (triangulation type dependent).
Multiple tests to check for output based on a general predicate and the
implemented IteratorFilters.
The output of GridTools::compute_ghost_cell_halo_layer is tested against
a distributed triangulation, ensuring that we return all of the ghost cells
on each processor.

18 files changed:
doc/news/changes.h
include/deal.II/grid/filtered_iterator.h
include/deal.II/grid/grid_tools.h
source/grid/grid_tools.cc
source/grid/grid_tools.inst.in
tests/deal.II/grid_tools_halo_layer_01.cc [new file with mode: 0644]
tests/deal.II/grid_tools_halo_layer_01.output [new file with mode: 0644]
tests/deal.II/grid_tools_halo_layer_02.cc [new file with mode: 0644]
tests/deal.II/grid_tools_halo_layer_02.output [new file with mode: 0644]
tests/deal.II/grid_tools_halo_layer_03.cc [new file with mode: 0644]
tests/deal.II/grid_tools_halo_layer_03.output [new file with mode: 0644]
tests/deal.II/grid_tools_halo_layer_04.cc [new file with mode: 0644]
tests/deal.II/grid_tools_halo_layer_04.output [new file with mode: 0644]
tests/deal.II/grid_tools_halo_layer_05.cc [new file with mode: 0644]
tests/deal.II/grid_tools_halo_layer_05.output [new file with mode: 0644]
tests/deal.II/grid_tools_halo_layer_ghost_cells.cc [new file with mode: 0644]
tests/deal.II/grid_tools_halo_layer_ghost_cells.mpirun=2.output [new file with mode: 0644]
tests/deal.II/grid_tools_halo_layer_ghost_cells.mpirun=3.output [new file with mode: 0644]

index de656a9d3e0ff9cb8fc48b0762ad0e522093a48b..2e6bdc88bd52163bdf116344f8e7b6de7896c736 100644 (file)
@@ -191,6 +191,14 @@ inconvenience this causes.
 
 
 <ol>
+  <li> New: There are now a collection of functions named GridTools::compute_active_cell_halo_layer()
+  that determine which cells form a layer around a specified subdomain. There is also a function
+  GridTools::compute_ghost_cell_halo_layer() that returns the smallest layer of ghost cells around
+  all locally relevant cells.
+  <br>
+  (Jean-Paul Pelteret, Denis Davydov, Wolfgang Bangerth, 2015/08/21)
+  </li>
+
   <li> Documentation: How to set up a testsuite in a user project is now
   properly documented.
   <br>
index 7b4c804056451d5d654ca394c928ab09c3943534..39ddaae81319ec97954f91a294ff4464892bc577 100644 (file)
@@ -135,7 +135,7 @@ namespace IteratorFilters
    * Filter for iterators that evaluates to true if either the iterator is
    * past the end or the subdomain id of the object pointed to is equal to a
    * value given to the constructor, assuming that the iterator allows
-   * querying for a subdomain id).
+   * querying for a subdomain id.
    *
    * @ingroup Iterators
    */
@@ -204,6 +204,99 @@ namespace IteratorFilters
     template <class Iterator>
     bool operator () (const Iterator &i) const;
   };
+
+
+  /**
+   * Filter for iterators that evaluates to true if the iterator of the object
+   * pointed to is equal to a value or set of values given to the constructor,
+   * assuming that the iterator allows querying for a material id.
+   *
+   * @author Jean-Paul Pelteret, Denis Davydov, 2015
+   *
+   * @ingroup Iterators
+   */
+  class MaterialIdEqualTo
+  {
+  public:
+    /**
+     * Constructor. Store the material id which iterators shall have to be
+     * evaluated to true and state if the iterator must be locally owned.
+     */
+    MaterialIdEqualTo (const types::material_id material_id,
+                       const bool only_locally_owned = false);
+
+    /**
+     * Constructor. Store a collection of material ids which iterators
+     * shall have to be evaluated to true and state if the iterator must be
+     * locally owned.
+     */
+    MaterialIdEqualTo (const std::set<types::material_id> material_ids,
+                       const bool only_locally_owned = false);
+
+    /**
+     * Evaluation operator. Returns true if the material id of the object
+     * pointed to is equal within the stored set of value allowable values
+     * and, if required, if the cell is locally owned.
+     */
+    template <class Iterator>
+    bool operator () (const Iterator &i) const;
+
+  protected:
+    /**
+     * Stored value to compare the material id with.
+     */
+    const std::set<types::material_id> material_ids;
+    /**
+     * Flag stating whether only locally owned cells must return true.
+     */
+    const bool only_locally_owned;
+  };
+
+  /**
+   * Filter for iterators that evaluates to true if the iterator of the object
+   * pointed to is equal to a value or set of values given to the constructor,
+   * assuming that the iterator allows querying for an active FE index.
+   *
+   * @author Jean-Paul Pelteret, Denis Davydov, 2015
+   *
+   * @ingroup Iterators
+   */
+  class ActiveFEIndexEqualTo
+  {
+  public:
+    /**
+     * Constructor. Store the active FE index which iterators shall have to be
+     * evaluated to true and state if the iterator must be locally owned.
+     */
+    ActiveFEIndexEqualTo (const unsigned int active_fe_index,
+                          const bool only_locally_owned = false);
+
+    /**
+     * Constructor. Store a collection of active FE indices which iterators
+     * shall have to be evaluated to true and state if the iterator must be
+     * locally owned.
+     */
+    ActiveFEIndexEqualTo (const std::set<unsigned int> active_fe_indices,
+                          const bool only_locally_owned = false);
+
+    /**
+     * Evaluation operator. Returns true if the active FE index of the object
+     * pointed to is equal within the stored set of value allowable values
+     * and, if required, if the cell is locally owned.
+     */
+    template <class Iterator>
+    bool operator () (const Iterator &i) const;
+
+  protected:
+    /**
+     * Stored value to compare the material id with.
+     */
+    const std::set<unsigned int> active_fe_indices;
+    /**
+     * Flag stating whether only locally owned cells must return true.
+     */
+    const bool only_locally_owned;
+  };
 }
 
 
@@ -1043,6 +1136,73 @@ namespace IteratorFilters
   {
     return (i->is_locally_owned_on_level());
   }
+
+
+
+// ---------------- IteratorFilters::MaterialIdEqualTo ---------
+  inline
+  MaterialIdEqualTo::MaterialIdEqualTo (const types::material_id material_id,
+                                        const bool only_locally_owned)
+    :
+    material_ids ({material_id}),
+               only_locally_owned (only_locally_owned)
+  {}
+
+
+
+  inline
+  MaterialIdEqualTo::MaterialIdEqualTo (const std::set<types::material_id> material_ids,
+                                        const bool only_locally_owned)
+    :
+    material_ids (material_ids),
+    only_locally_owned (only_locally_owned)
+  {}
+
+
+
+  template <class Iterator>
+  inline
+  bool
+  MaterialIdEqualTo::operator () (const Iterator &i) const
+  {
+    return only_locally_owned == true ?
+           (material_ids.find(i->material_id()) != material_ids.end() && i->is_locally_owned()):
+           material_ids.find(i->material_id()) != material_ids.end();
+  }
+
+
+
+// ---------------- IteratorFilters::ActiveFEIndexEqualTo ---------
+  inline
+  ActiveFEIndexEqualTo::ActiveFEIndexEqualTo (const unsigned int active_fe_index,
+                                              const bool only_locally_owned)
+    :
+    active_fe_indices ({active_fe_index}),
+                    only_locally_owned (only_locally_owned)
+  {}
+
+
+
+  inline
+  ActiveFEIndexEqualTo::ActiveFEIndexEqualTo (const std::set<unsigned int> active_fe_indices,
+                                              const bool only_locally_owned)
+    :
+    active_fe_indices (active_fe_indices),
+    only_locally_owned (only_locally_owned)
+  {}
+
+
+
+  template <class Iterator>
+  inline
+  bool
+  ActiveFEIndexEqualTo::operator () (const Iterator &i) const
+  {
+    return only_locally_owned == true ?
+           (active_fe_indices.find(i->active_fe_index()) != active_fe_indices.end() && i->is_locally_owned()):
+           active_fe_indices.find(i->active_fe_index()) != active_fe_indices.end();
+  }
+
 }
 
 
index daf233effd87e727b0244d506ed094029330b165..a3155cf56d8acfb4f7be591e0d3f895780000dfc 100644 (file)
@@ -556,6 +556,80 @@ namespace GridTools
   get_active_neighbors (const typename Container::active_cell_iterator        &cell,
                         std::vector<typename Container::active_cell_iterator> &active_neighbors);
 
+  /**
+   * Extract and return the active cell layer around a subdomain (set of active
+   * cells) in the @p container (i.e. those that share a common set of vertices
+   * with the subdomain but are not a part of it).
+   * Here, the "subdomain" consists of exactly all of those cells for which the
+   * @p predicate returns @p true.
+   *
+   * An example of a custom predicate is one that checks for a given material id
+   * @code
+   * template<int dim>
+   * bool
+   * pred_mat_id(const typename Triangulation<dim>::active_cell_iterator & cell)
+   * {
+   *   return cell->material_id() ==  1;
+   * }
+   * @endcode
+   * and we can then extract the layer of cells around this material with the
+   * following call:
+   * @code
+   * GridTools::compute_active_cell_halo_layer(tria, pred_mat_id<dim>);
+   * @endcode
+   *
+   * Predicates that are frequently useful can be found in namespace IteratorFilters.
+   * For example, it is possible to extracting a layer based on material id
+   * @code
+   * GridTools::compute_active_cell_halo_layer(tria,
+   *                                           IteratorFilters::MaterialIdEqualTo(1, true));
+   * @endcode
+   * or based on a set of active FE indices for an hp::DoFHandler
+   * @code
+   * GridTools::compute_active_cell_halo_layer(hp_dof_handler,
+   *                                           IteratorFilters::ActiveFEIndexEqualTo({1,2}, true));
+   * @endcode
+   * Note that in the last two examples we ensure that the predicate returns
+   * true only for locally owned cells. This means that the halo layer will
+   * not contain any artificial cells.
+   *
+   * @tparam Container A type that satisfies the requirements of a mesh
+   * container (see @ref GlossMeshAsAContainer).
+   * @param[in] container A mesh container (i.e. objects of type Triangulation,
+   * DoFHandler, or hp::DoFHandler).
+   * @param[in] 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.
+   * @return A list of active cells sharing at least one common vertex with the
+   * predicated subdomain.
+   *
+   * @author Jean-Paul Pelteret, Denis Davydov, Wolfgang Bangerth, 2015
+   */
+  template <class Container>
+  std::vector<typename Container::active_cell_iterator>
+  compute_active_cell_halo_layer (const Container                                                                    &container,
+                                  const std_cxx11::function<bool (const typename Container::active_cell_iterator &)> &predicate);
+
+  /**
+   * Extract and return ghost cells which are the active cell layer
+   * 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 processor, but for parallel::distributed::Triangulation
+   * this will return all the ghost cells.
+   *
+   * @tparam Container A type that satisfies the requirements of a mesh
+   * container (see @ref GlossMeshAsAContainer).
+   * @param[in] container A mesh container (i.e. objects of type Triangulation,
+   * DoFHandler, or hp::DoFHandler).
+   * @return A list of ghost cells
+   *
+   * @author Jean-Paul Pelteret, Denis Davydov, Wolfgang Bangerth, 2015
+   */
+  template <class Container>
+  std::vector<typename Container::active_cell_iterator>
+  compute_ghost_cell_halo_layer (const Container &container);
+
+
   /*@}*/
   /**
    * @name Partitions and subdomains of triangulations
@@ -1280,6 +1354,7 @@ namespace GridTools
 
 namespace GridTools
 {
+
   template <int dim, typename Predicate, int spacedim>
   void transform (const Predicate    &predicate,
                   Triangulation<dim, spacedim> &triangulation)
@@ -1429,7 +1504,6 @@ namespace GridTools
 
 
 
-
 // declaration of explicit specializations
 
   template <>
index 7a04c594b9015358092f3d7436d6929bce01f299..48cc49d8dbff8de2737f00164fe7471833e1d328 100644 (file)
@@ -26,6 +26,7 @@
 #include <deal.II/lac/constraint_matrix.h>
 #include <deal.II/lac/sparsity_pattern.h>
 #include <deal.II/lac/sparsity_tools.h>
+#include <deal.II/grid/filtered_iterator.h>
 #include <deal.II/grid/tria.h>
 #include <deal.II/distributed/tria.h>
 #include <deal.II/grid/tria_accessor.h>
@@ -1501,6 +1502,95 @@ next_cell:
   }
 
 
+  namespace
+  {
+
+    template<class Container>
+    bool
+    contains_locally_owned_cells (const std::vector<typename Container::active_cell_iterator> &cells)
+    {
+      for (typename std::vector<typename Container::active_cell_iterator>::const_iterator
+           it = cells.begin(); it != cells.end(); ++it)
+        {
+          if ((*it)->is_locally_owned())
+            return true;
+        }
+      return false;
+    }
+
+    template<class Container>
+    bool
+    contains_artificial_cells (const std::vector<typename Container::active_cell_iterator> &cells)
+    {
+      for (typename std::vector<typename Container::active_cell_iterator>::const_iterator
+           it = cells.begin(); it != cells.end(); ++it)
+        {
+          if ((*it)->is_artificial())
+            return true;
+        }
+      return false;
+    }
+
+  }
+
+
+
+  template <class Container>
+  std::vector<typename Container::active_cell_iterator>
+  compute_active_cell_halo_layer (const Container                                                                    &container,
+                                  const std_cxx11::function<bool (const typename Container::active_cell_iterator &)> &predicate)
+  {
+    std::vector<typename Container::active_cell_iterator> active_halo_layer;
+    std::vector<bool> locally_active_vertices_on_subdomain (get_tria(container).n_vertices(),
+                                                            false);
+
+    // Find the cells for which the predicate is true
+    // These are the cells around which we wish to construct
+    // the halo layer
+    for (typename Container::active_cell_iterator
+         cell = container.begin_active();
+         cell != container.end(); ++cell)
+      if (predicate(cell)) // True predicate --> Part of subdomain
+        for (unsigned int v=0; v<GeometryInfo<Container::dimension>::vertices_per_cell; ++v)
+          locally_active_vertices_on_subdomain[cell->vertex_index(v)] = true;
+
+    // Find the cells that do not conform to the predicate
+    // but share a vertex with the selected subdomain
+    // These comprise the halo layer
+    for (typename Container::active_cell_iterator
+         cell = container.begin_active();
+         cell != container.end(); ++cell)
+      if (!predicate(cell)) // False predicate --> Potential halo cell
+        for (unsigned int v=0; v<GeometryInfo<Container::dimension>::vertices_per_cell; ++v)
+          if (locally_active_vertices_on_subdomain[cell->vertex_index(v)] == true)
+            {
+              active_halo_layer.push_back(cell);
+              break;
+            }
+
+    return active_halo_layer;
+  }
+
+
+
+  template <class Container>
+  std::vector<typename Container::active_cell_iterator>
+  compute_ghost_cell_halo_layer (const Container &container)
+  {
+    const std::vector<typename Container::active_cell_iterator>
+    active_halo_layer = compute_active_cell_halo_layer (container, IteratorFilters::LocallyOwnedCell());
+
+    // Check that we never return locally owned or artificial cells
+    // What is left should only be the ghost cells
+    Assert(contains_locally_owned_cells<Container>(active_halo_layer) == false,
+           ExcMessage("Halo layer contains locally owned cells"));
+    Assert(contains_artificial_cells<Container>(active_halo_layer) == false,
+           ExcMessage("Halo layer contains artificial cells"));
+
+    return active_halo_layer;
+  }
+
+
 
   template <int dim, int spacedim>
   void
index c1c5fe585c539ecada7c3360565abfe87e702532..7123da054c818095704cdf2705d925d44c2c3bc6 100644 (file)
@@ -40,6 +40,15 @@ for (X : TRIANGULATION_AND_DOFHANDLERS; deal_II_dimension : DIMENSIONS ; deal_II
                                    const X &,
                                    const Point<deal_II_space_dimension> &);
 
+  template
+    std::vector<dealii::internal::ActiveCellIterator<deal_II_dimension, deal_II_space_dimension, X>::type>
+    compute_active_cell_halo_layer (const X &,
+                                    const std_cxx11::function<bool (const dealii::internal::ActiveCellIterator<deal_II_dimension, deal_II_space_dimension, X>::type&)> &);
+
+  template
+    std::vector<dealii::internal::ActiveCellIterator<deal_II_dimension, deal_II_space_dimension, X>::type>
+    compute_ghost_cell_halo_layer (const X &);
+
   template
     std::list<std::pair<X::cell_iterator, X::cell_iterator> >
     get_finest_common_cells (const X &mesh_1,
@@ -84,6 +93,16 @@ for (deal_II_dimension : DIMENSIONS ; deal_II_space_dimension : SPACE_DIMENSIONS
                                    const parallel::distributed::Triangulation<deal_II_dimension,deal_II_space_dimension> &,
                                    const Point<deal_II_space_dimension> &);
 
+  template
+    std::vector<dealii::internal::ActiveCellIterator<deal_II_dimension, deal_II_space_dimension, parallel::distributed::Triangulation<deal_II_dimension, deal_II_space_dimension>>::type>
+    compute_active_cell_halo_layer (const parallel::distributed::Triangulation<deal_II_dimension, deal_II_space_dimension> &,
+                                    const std_cxx11::function<bool (const dealii::internal::ActiveCellIterator<deal_II_dimension, deal_II_space_dimension, 
+                                                                                                               parallel::distributed::Triangulation<deal_II_dimension, deal_II_space_dimension>>::type&)> &);
+
+  template
+    std::vector<dealii::internal::ActiveCellIterator<deal_II_dimension, deal_II_space_dimension, parallel::distributed::Triangulation<deal_II_dimension, deal_II_space_dimension>>::type>
+    compute_ghost_cell_halo_layer (const parallel::distributed::Triangulation<deal_II_dimension, deal_II_space_dimension> &);
+
   template
     std::list<std::pair<parallel::distributed::Triangulation<deal_II_dimension,deal_II_space_dimension>::cell_iterator, parallel::distributed::Triangulation<deal_II_dimension,deal_II_space_dimension>::cell_iterator> >
     get_finest_common_cells (const parallel::distributed::Triangulation<deal_II_dimension,deal_II_space_dimension> &mesh_1,
@@ -190,8 +209,6 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension :  SPACE_DIMENSIONS
        const hp::DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
        const Point<deal_II_space_dimension> &);
 
-
-
     template
       void get_subdomain_association (const Triangulation<deal_II_dimension, deal_II_space_dimension>  &,
                                       std::vector<types::subdomain_id> &);
diff --git a/tests/deal.II/grid_tools_halo_layer_01.cc b/tests/deal.II/grid_tools_halo_layer_01.cc
new file mode 100644 (file)
index 0000000..911899c
--- /dev/null
@@ -0,0 +1,122 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2001 - 2014 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.
+//
+// ---------------------------------------------------------------------
+
+
+
+#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 halo:" << std::endl;
+  write_mat_id_to_file(tria);
+  // Write to file to visually check result
+  {
+    const std::string filename = "grid_no_halo_" + std::to_string(dim) + "d.vtk";
+    std::ofstream f(filename);
+    GridOut().write_vtk (tria, f);
+  }
+
+  // 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_halo_layer(tria, pred_mat_id<dim>); // 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 halo:" << std::endl;
+  write_mat_id_to_file(tria);
+  // Write to file to visually check result
+  {
+    const std::string filename = "grid_with_halo_" + std::to_string(dim) + "d.vtk";
+    std::ofstream f(filename);
+    GridOut().write_vtk (tria, f);
+  }
+}
+
+
+int main ()
+{
+  initlog();
+  deallog.depth_console(0);
+
+  test<2> ();
+  test<3> ();
+
+  return 0;
+}
diff --git a/tests/deal.II/grid_tools_halo_layer_01.output b/tests/deal.II/grid_tools_halo_layer_01.output
new file mode 100644 (file)
index 0000000..e91b7c2
--- /dev/null
@@ -0,0 +1,171 @@
+
+DEAL::dim = 2
+DEAL::Grid without halo:
+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 halo:
+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 halo:
+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 halo:
+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/deal.II/grid_tools_halo_layer_02.cc b/tests/deal.II/grid_tools_halo_layer_02.cc
new file mode 100644 (file)
index 0000000..8cf6744
--- /dev/null
@@ -0,0 +1,116 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2001 - 2014 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.
+//
+// ---------------------------------------------------------------------
+
+
+
+#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 <deal.II/grid/filtered_iterator.h>
+
+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 halo:" << std::endl;
+  write_mat_id_to_file(tria);
+  // Write to file to visually check result
+  {
+    const std::string filename = "grid_no_halo_" + std::to_string(dim) + "d.vtk";
+    std::ofstream f(filename);
+    GridOut().write_vtk (tria, f);
+  }
+
+  // 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_halo_layer(tria, IteratorFilters::MaterialIdEqualTo(2, true));
+  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 halo:" << std::endl;
+  write_mat_id_to_file(tria);
+  // Write to file to visually check result
+  {
+    const std::string filename = "grid_with_halo_" + std::to_string(dim) + "d.vtk";
+    std::ofstream f(filename);
+    GridOut().write_vtk (tria, f);
+  }
+}
+
+
+int main ()
+{
+  initlog();
+  deallog.depth_console(0);
+
+  test<2> ();
+  test<3> ();
+
+  return 0;
+}
diff --git a/tests/deal.II/grid_tools_halo_layer_02.output b/tests/deal.II/grid_tools_halo_layer_02.output
new file mode 100644 (file)
index 0000000..e91b7c2
--- /dev/null
@@ -0,0 +1,171 @@
+
+DEAL::dim = 2
+DEAL::Grid without halo:
+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 halo:
+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 halo:
+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 halo:
+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/deal.II/grid_tools_halo_layer_03.cc b/tests/deal.II/grid_tools_halo_layer_03.cc
new file mode 100644 (file)
index 0000000..c2de9ae
--- /dev/null
@@ -0,0 +1,121 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2001 - 2014 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.
+//
+// ---------------------------------------------------------------------
+
+
+
+#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 <deal.II/grid/filtered_iterator.h>
+
+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)
+    {
+      if (cell->center()[0] < 0.25)
+        cell->set_material_id(2);
+      else
+        cell->set_material_id(3);
+    }
+    else
+      cell->set_material_id(1);
+  }
+
+  deallog << "Grid without halo:" << std::endl;
+  write_mat_id_to_file(tria);
+  // Write to file to visually check result
+  {
+    const std::string filename = "grid_no_halo_" + std::to_string(dim) + "d.vtk";
+    std::ofstream f(filename);
+    GridOut().write_vtk (tria, f);
+  }
+
+  // Compute a halo layer around material ids 2 and 3 and set it to material id 4
+  const std::vector<cell_iterator> active_halo_layer
+    = GridTools::compute_active_cell_halo_layer(tria, IteratorFilters::MaterialIdEqualTo({2,3}, true));
+  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(4);
+  }
+
+  deallog << "Grid with halo:" << std::endl;
+  write_mat_id_to_file(tria);
+  // Write to file to visually check result
+  {
+    const std::string filename = "grid_with_halo_" + std::to_string(dim) + "d.vtk";
+    std::ofstream f(filename);
+    GridOut().write_vtk (tria, f);
+  }
+}
+
+
+int main ()
+{
+  initlog();
+  deallog.depth_console(0);
+
+  test<2> ();
+  test<3> ();
+
+  return 0;
+}
diff --git a/tests/deal.II/grid_tools_halo_layer_03.output b/tests/deal.II/grid_tools_halo_layer_03.output
new file mode 100644 (file)
index 0000000..b8b799e
--- /dev/null
@@ -0,0 +1,171 @@
+
+DEAL::dim = 2
+DEAL::Grid without halo:
+DEAL::0 2
+DEAL::1 3
+DEAL::2 2
+DEAL::3 3
+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 halo:
+DEAL::0 2
+DEAL::1 3
+DEAL::2 2
+DEAL::3 3
+DEAL::4 4
+DEAL::5 1
+DEAL::6 4
+DEAL::7 1
+DEAL::8 4
+DEAL::9 4
+DEAL::10 1
+DEAL::11 1
+DEAL::12 4
+DEAL::13 1
+DEAL::14 1
+DEAL::15 1
+DEAL::
+DEAL::dim = 3
+DEAL::Grid without halo:
+DEAL::0 2
+DEAL::1 3
+DEAL::2 2
+DEAL::3 3
+DEAL::4 2
+DEAL::5 3
+DEAL::6 2
+DEAL::7 3
+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 halo:
+DEAL::0 2
+DEAL::1 3
+DEAL::2 2
+DEAL::3 3
+DEAL::4 2
+DEAL::5 3
+DEAL::6 2
+DEAL::7 3
+DEAL::8 4
+DEAL::9 1
+DEAL::10 4
+DEAL::11 1
+DEAL::12 4
+DEAL::13 1
+DEAL::14 4
+DEAL::15 1
+DEAL::16 4
+DEAL::17 4
+DEAL::18 1
+DEAL::19 1
+DEAL::20 4
+DEAL::21 4
+DEAL::22 1
+DEAL::23 1
+DEAL::24 4
+DEAL::25 1
+DEAL::26 1
+DEAL::27 1
+DEAL::28 4
+DEAL::29 1
+DEAL::30 1
+DEAL::31 1
+DEAL::32 4
+DEAL::33 4
+DEAL::34 4
+DEAL::35 4
+DEAL::36 1
+DEAL::37 1
+DEAL::38 1
+DEAL::39 1
+DEAL::40 4
+DEAL::41 1
+DEAL::42 4
+DEAL::43 1
+DEAL::44 1
+DEAL::45 1
+DEAL::46 1
+DEAL::47 1
+DEAL::48 4
+DEAL::49 4
+DEAL::50 1
+DEAL::51 1
+DEAL::52 1
+DEAL::53 1
+DEAL::54 1
+DEAL::55 1
+DEAL::56 4
+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/deal.II/grid_tools_halo_layer_04.cc b/tests/deal.II/grid_tools_halo_layer_04.cc
new file mode 100644 (file)
index 0000000..9db8a25
--- /dev/null
@@ -0,0 +1,146 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2001 - 2014 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.
+//
+// ---------------------------------------------------------------------
+
+
+
+#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/filtered_iterator.h>
+#include <deal.II/numerics/data_out.h>
+#include <deal.II/hp/dof_handler.h>
+
+template <int dim>
+void
+write_active_fe_index_to_file (const hp::DoFHandler<dim> & dof_handler)
+{
+  int count = 0;
+  typename hp::DoFHandler<dim>::active_cell_iterator
+  cell = dof_handler.begin_active(),
+  endc = dof_handler.end();
+  for (; cell != endc; ++cell, ++count)
+  {
+    deallog
+      << count << " "
+      << cell->active_fe_index()
+      << std::endl;
+  }
+  deallog << std::endl;
+}
+
+template <int dim>
+void
+write_vtk (const hp::DoFHandler<dim> & dof_handler, const std::string filename)
+{
+  Vector<double> active_fe_index (dof_handler.get_tria().n_active_cells());
+  int count = 0;
+  typename hp::DoFHandler<dim>::active_cell_iterator
+  cell = dof_handler.begin_active(),
+  endc = dof_handler.end();
+  for (; cell != endc; ++cell, ++count)
+  {
+    active_fe_index[count] = cell->active_fe_index();
+  }
+
+  const std::vector<DataComponentInterpretation::DataComponentInterpretation>
+  data_component_interpretation
+  (1, DataComponentInterpretation::component_is_scalar);
+  const std::vector<std::string> data_names (1, "active_fe_index");
+
+  DataOut<dim,hp::DoFHandler<dim> > data_out;
+  data_out.attach_dof_handler (dof_handler);
+  data_out.add_data_vector (active_fe_index, data_names,
+                            DataOut<dim,hp::DoFHandler<dim> >::type_cell_data,
+                            data_component_interpretation);
+  data_out.build_patches ();
+
+  std::ofstream output (filename.c_str());
+  data_out.write_vtk (output);
+}
+
+
+template <int dim>
+void test ()
+{
+  deallog << "dim = " << dim << std::endl;
+
+  Triangulation<dim> tria;
+  GridGenerator::hyper_cube(tria);
+  tria.refine_global(2);
+  hp::DoFHandler<dim> dof_handler (tria);
+
+  typedef typename hp::DoFHandler<dim>::active_cell_iterator cell_iterator;
+
+  // Mark a small block at the corner of the hypercube
+  cell_iterator
+  cell = dof_handler.begin_active(),
+  endc = dof_handler.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_active_fe_index(2);
+    else
+      cell->set_active_fe_index(1);
+  }
+
+  deallog << "Grid without halo:" << std::endl;
+  write_active_fe_index_to_file(dof_handler);
+  // Write to file to visually check result
+  {
+    const std::string filename = "grid_no_halo_" + std::to_string(dim) + "d.vtk";
+    write_vtk(dof_handler, filename);
+  }
+
+  // Compute a halo layer around active fe index 2 and set it to active fe index 3
+  std::vector<cell_iterator> active_halo_layer
+    = GridTools::compute_active_cell_halo_layer(dof_handler, IteratorFilters::ActiveFEIndexEqualTo(2, true));
+  AssertThrow(active_halo_layer.size() > 0, ExcMessage("No halo layer found."));
+  for (typename std::vector<cell_iterator>::iterator
+       it = active_halo_layer.begin();
+       it != active_halo_layer.end(); ++it)
+  {
+    (*it)->set_active_fe_index(3);
+  }
+
+  deallog << "Grid with halo:" << std::endl;
+  write_active_fe_index_to_file(dof_handler);
+  // Write to file to visually check result
+  {
+    const std::string filename = "grid_with_halo_" + std::to_string(dim) + "d.vtk";
+    write_vtk(dof_handler, filename);
+  }
+}
+
+
+int main ()
+{
+  initlog();
+  deallog.depth_console(0);
+
+  test<2> ();
+  test<3> ();
+
+  return 0;
+}
diff --git a/tests/deal.II/grid_tools_halo_layer_04.output b/tests/deal.II/grid_tools_halo_layer_04.output
new file mode 100644 (file)
index 0000000..e91b7c2
--- /dev/null
@@ -0,0 +1,171 @@
+
+DEAL::dim = 2
+DEAL::Grid without halo:
+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 halo:
+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 halo:
+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 halo:
+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/deal.II/grid_tools_halo_layer_05.cc b/tests/deal.II/grid_tools_halo_layer_05.cc
new file mode 100644 (file)
index 0000000..8edbba3
--- /dev/null
@@ -0,0 +1,152 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2001 - 2014 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.
+//
+// ---------------------------------------------------------------------
+
+
+
+#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/filtered_iterator.h>
+#include <deal.II/numerics/data_out.h>
+#include <deal.II/hp/dof_handler.h>
+
+template <int dim>
+void
+write_active_fe_index_to_file (const hp::DoFHandler<dim> & dof_handler)
+{
+  int count = 0;
+  typename hp::DoFHandler<dim>::active_cell_iterator
+  cell = dof_handler.begin_active(),
+  endc = dof_handler.end();
+  for (; cell != endc; ++cell, ++count)
+  {
+    deallog
+      << count << " "
+      << cell->active_fe_index()
+      << std::endl;
+  }
+  deallog << std::endl;
+}
+
+template <int dim>
+void
+write_vtk (const hp::DoFHandler<dim> & dof_handler, const std::string filename)
+{
+  Vector<double> active_fe_index (dof_handler.get_tria().n_active_cells());
+  int count = 0;
+  typename hp::DoFHandler<dim>::active_cell_iterator
+  cell = dof_handler.begin_active(),
+  endc = dof_handler.end();
+  for (; cell != endc; ++cell, ++count)
+  {
+    active_fe_index[count] = cell->active_fe_index();
+  }
+
+  const std::vector<DataComponentInterpretation::DataComponentInterpretation>
+  data_component_interpretation
+  (1, DataComponentInterpretation::component_is_scalar);
+  const std::vector<std::string> data_names (1, "active_fe_index");
+
+  DataOut<dim,hp::DoFHandler<dim> > data_out;
+  data_out.attach_dof_handler (dof_handler);
+  data_out.add_data_vector (active_fe_index, data_names,
+                            DataOut<dim,hp::DoFHandler<dim> >::type_cell_data,
+                            data_component_interpretation);
+  data_out.build_patches ();
+
+  std::ofstream output (filename.c_str());
+  data_out.write_vtk (output);
+}
+
+
+template <int dim>
+void test ()
+{
+  deallog << "dim = " << dim << std::endl;
+
+  Triangulation<dim> tria;
+  GridGenerator::hyper_cube(tria);
+  tria.refine_global(2);
+  hp::DoFHandler<dim> dof_handler (tria);
+
+  typedef typename hp::DoFHandler<dim>::active_cell_iterator cell_iterator;
+
+  // Mark a small block at the corner of the hypercube
+  cell_iterator
+  cell = dof_handler.begin_active(),
+  endc = dof_handler.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)
+    {
+
+      if (cell->center()[0] < 0.25)
+        cell->set_active_fe_index(2);
+      else
+        cell->set_active_fe_index(3);
+    }
+    else
+      cell->set_active_fe_index(1);
+  }
+
+  deallog << "Grid without halo:" << std::endl;
+  write_active_fe_index_to_file(dof_handler);
+  // Write to file to visually check result
+  {
+    const std::string filename = "grid_no_halo_" + std::to_string(dim) + "d.vtk";
+    write_vtk(dof_handler, filename);
+  }
+
+  // Compute a halo layer around active fe indices 2,3 and set it to active fe index 4
+  std::vector<cell_iterator> active_halo_layer
+    = GridTools::compute_active_cell_halo_layer(dof_handler, IteratorFilters::ActiveFEIndexEqualTo({2,3}, true));
+  AssertThrow(active_halo_layer.size() > 0, ExcMessage("No halo layer found."));
+  for (typename std::vector<cell_iterator>::iterator
+       it = active_halo_layer.begin();
+       it != active_halo_layer.end(); ++it)
+  {
+    (*it)->set_active_fe_index(4);
+  }
+
+  deallog << "Grid with halo:" << std::endl;
+  write_active_fe_index_to_file(dof_handler);
+  // Write to file to visually check result
+  {
+    const std::string filename = "grid_with_halo_" + std::to_string(dim) + "d.vtk";
+    write_vtk(dof_handler, filename);
+  }
+}
+
+
+int main ()
+{
+  initlog();
+  deallog.depth_console(0);
+
+  test<2> ();
+  test<3> ();
+
+  return 0;
+}
diff --git a/tests/deal.II/grid_tools_halo_layer_05.output b/tests/deal.II/grid_tools_halo_layer_05.output
new file mode 100644 (file)
index 0000000..b8b799e
--- /dev/null
@@ -0,0 +1,171 @@
+
+DEAL::dim = 2
+DEAL::Grid without halo:
+DEAL::0 2
+DEAL::1 3
+DEAL::2 2
+DEAL::3 3
+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 halo:
+DEAL::0 2
+DEAL::1 3
+DEAL::2 2
+DEAL::3 3
+DEAL::4 4
+DEAL::5 1
+DEAL::6 4
+DEAL::7 1
+DEAL::8 4
+DEAL::9 4
+DEAL::10 1
+DEAL::11 1
+DEAL::12 4
+DEAL::13 1
+DEAL::14 1
+DEAL::15 1
+DEAL::
+DEAL::dim = 3
+DEAL::Grid without halo:
+DEAL::0 2
+DEAL::1 3
+DEAL::2 2
+DEAL::3 3
+DEAL::4 2
+DEAL::5 3
+DEAL::6 2
+DEAL::7 3
+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 halo:
+DEAL::0 2
+DEAL::1 3
+DEAL::2 2
+DEAL::3 3
+DEAL::4 2
+DEAL::5 3
+DEAL::6 2
+DEAL::7 3
+DEAL::8 4
+DEAL::9 1
+DEAL::10 4
+DEAL::11 1
+DEAL::12 4
+DEAL::13 1
+DEAL::14 4
+DEAL::15 1
+DEAL::16 4
+DEAL::17 4
+DEAL::18 1
+DEAL::19 1
+DEAL::20 4
+DEAL::21 4
+DEAL::22 1
+DEAL::23 1
+DEAL::24 4
+DEAL::25 1
+DEAL::26 1
+DEAL::27 1
+DEAL::28 4
+DEAL::29 1
+DEAL::30 1
+DEAL::31 1
+DEAL::32 4
+DEAL::33 4
+DEAL::34 4
+DEAL::35 4
+DEAL::36 1
+DEAL::37 1
+DEAL::38 1
+DEAL::39 1
+DEAL::40 4
+DEAL::41 1
+DEAL::42 4
+DEAL::43 1
+DEAL::44 1
+DEAL::45 1
+DEAL::46 1
+DEAL::47 1
+DEAL::48 4
+DEAL::49 4
+DEAL::50 1
+DEAL::51 1
+DEAL::52 1
+DEAL::53 1
+DEAL::54 1
+DEAL::55 1
+DEAL::56 4
+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/deal.II/grid_tools_halo_layer_ghost_cells.cc b/tests/deal.II/grid_tools_halo_layer_ghost_cells.cc
new file mode 100644 (file)
index 0000000..f8fc22f
--- /dev/null
@@ -0,0 +1,102 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2001 - 2014 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.
+//
+// ---------------------------------------------------------------------
+
+
+
+#include "../tests.h"
+#include <deal.II/base/logstream.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/grid_out.h>
+#include <deal.II/grid/cell_id.h>
+
+#include <algorithm>
+
+template <int dim>
+void test ()
+{
+  const MPI_Comm &mpi_communicator = MPI_COMM_WORLD;
+  deallog << "dim = " << dim << std::endl;
+
+  parallel::distributed::Triangulation<dim> tria (mpi_communicator);
+  GridGenerator::hyper_cube(tria);
+  tria.refine_global(2);
+
+//  // Write to file for visual inspection
+//  const std::string filename = "grid_" + std::to_string(dim) + "d";
+//  tria.write_mesh_vtk(filename.c_str());
+
+  typedef typename parallel::distributed::Triangulation<dim>::active_cell_iterator cell_iterator;
+
+  // Mark a small block at the corner of the hypercube
+  std::vector<cell_iterator> ghost_cells_tria;
+  cell_iterator
+  cell = tria.begin_active(),
+  endc = tria.end();
+  for (; cell != endc; ++cell)
+  {
+    if (cell->is_ghost() == true)
+      ghost_cells_tria.push_back(cell);
+  }
+  std::sort(ghost_cells_tria.begin(),
+            ghost_cells_tria.end());
+
+  // Compute a halo layer around the locally owned cells.
+  // These should all be ghost cells
+  std::vector<cell_iterator> ghost_cell_halo_layer
+    = GridTools::compute_ghost_cell_halo_layer(tria);
+  AssertThrow(ghost_cell_halo_layer.size() > 0, ExcMessage("Ghost cell halo layer found."));
+  AssertThrow(ghost_cell_halo_layer.size() == ghost_cells_tria.size(), ExcMessage("Ghost cell halo layer wrong size."));
+  std::sort(ghost_cell_halo_layer.begin(),
+            ghost_cell_halo_layer.end());
+
+  for (unsigned int proc = 0; proc < Utilities::MPI::n_mpi_processes(mpi_communicator); ++proc)
+  {
+    if (proc == Utilities::MPI::this_mpi_process(mpi_communicator))
+    {
+      for (typename std::vector<cell_iterator>::const_iterator
+           it_1 = ghost_cells_tria.begin(),
+           it_2 = ghost_cell_halo_layer.begin();
+           it_1 != ghost_cells_tria.end() &&
+           it_2 != ghost_cell_halo_layer.end();
+           ++it_1, ++it_2)
+      {
+          const cell_iterator & cell_1 = *it_1;
+          const cell_iterator & cell_2 = *it_2;
+          AssertThrow(cell_1->is_ghost() == true, ExcMessage("Cell is not a ghost cell!"));
+          AssertThrow(cell_2->is_ghost() == true, ExcMessage("Halo cell is not a ghost cell!"));
+          deallog
+            << "Ghost " << cell_1->level() << " " << cell_1->index() << " " << cell_1->id() << " " << cell_1->id().to_string() << " "
+            << "Halo "  << cell_2->level() << " " << cell_2->index() << " " << cell_2->id() << " " << cell_2->id().to_string()
+            << std::endl;
+          AssertThrow(cell_2 == cell_1, ExcMessage("Halo cell is not identical to ghost cell."));
+      }
+    }
+  }
+}
+
+
+int main (int argc, char *argv[])
+{
+  Utilities::MPI::MPI_InitFinalize mpi_initialization (argc, argv, 1);
+  MPILogInitAll log;
+  deallog.depth_console(0);
+
+  test<2> ();
+  test<3> ();
+
+  return 0;
+}
diff --git a/tests/deal.II/grid_tools_halo_layer_ghost_cells.mpirun=2.output b/tests/deal.II/grid_tools_halo_layer_ghost_cells.mpirun=2.output
new file mode 100644 (file)
index 0000000..a50adba
--- /dev/null
@@ -0,0 +1,47 @@
+
+DEAL:0::dim = 2
+DEAL:0::Ghost 2 8 0_2:20 0_2:20 Halo 2 8 0_2:20 0_2:20
+DEAL:0::Ghost 2 9 0_2:21 0_2:21 Halo 2 9 0_2:21 0_2:21
+DEAL:0::Ghost 2 12 0_2:30 0_2:30 Halo 2 12 0_2:30 0_2:30
+DEAL:0::Ghost 2 13 0_2:31 0_2:31 Halo 2 13 0_2:31 0_2:31
+DEAL:0::dim = 3
+DEAL:0::Ghost 2 32 0_2:40 0_2:40 Halo 2 32 0_2:40 0_2:40
+DEAL:0::Ghost 2 33 0_2:41 0_2:41 Halo 2 33 0_2:41 0_2:41
+DEAL:0::Ghost 2 34 0_2:42 0_2:42 Halo 2 34 0_2:42 0_2:42
+DEAL:0::Ghost 2 35 0_2:43 0_2:43 Halo 2 35 0_2:43 0_2:43
+DEAL:0::Ghost 2 40 0_2:50 0_2:50 Halo 2 40 0_2:50 0_2:50
+DEAL:0::Ghost 2 41 0_2:51 0_2:51 Halo 2 41 0_2:51 0_2:51
+DEAL:0::Ghost 2 42 0_2:52 0_2:52 Halo 2 42 0_2:52 0_2:52
+DEAL:0::Ghost 2 43 0_2:53 0_2:53 Halo 2 43 0_2:53 0_2:53
+DEAL:0::Ghost 2 48 0_2:60 0_2:60 Halo 2 48 0_2:60 0_2:60
+DEAL:0::Ghost 2 49 0_2:61 0_2:61 Halo 2 49 0_2:61 0_2:61
+DEAL:0::Ghost 2 50 0_2:62 0_2:62 Halo 2 50 0_2:62 0_2:62
+DEAL:0::Ghost 2 51 0_2:63 0_2:63 Halo 2 51 0_2:63 0_2:63
+DEAL:0::Ghost 2 56 0_2:70 0_2:70 Halo 2 56 0_2:70 0_2:70
+DEAL:0::Ghost 2 57 0_2:71 0_2:71 Halo 2 57 0_2:71 0_2:71
+DEAL:0::Ghost 2 58 0_2:72 0_2:72 Halo 2 58 0_2:72 0_2:72
+DEAL:0::Ghost 2 59 0_2:73 0_2:73 Halo 2 59 0_2:73 0_2:73
+
+DEAL:1::dim = 2
+DEAL:1::Ghost 2 2 0_2:02 0_2:02 Halo 2 2 0_2:02 0_2:02
+DEAL:1::Ghost 2 3 0_2:03 0_2:03 Halo 2 3 0_2:03 0_2:03
+DEAL:1::Ghost 2 6 0_2:12 0_2:12 Halo 2 6 0_2:12 0_2:12
+DEAL:1::Ghost 2 7 0_2:13 0_2:13 Halo 2 7 0_2:13 0_2:13
+DEAL:1::dim = 3
+DEAL:1::Ghost 2 4 0_2:04 0_2:04 Halo 2 4 0_2:04 0_2:04
+DEAL:1::Ghost 2 5 0_2:05 0_2:05 Halo 2 5 0_2:05 0_2:05
+DEAL:1::Ghost 2 6 0_2:06 0_2:06 Halo 2 6 0_2:06 0_2:06
+DEAL:1::Ghost 2 7 0_2:07 0_2:07 Halo 2 7 0_2:07 0_2:07
+DEAL:1::Ghost 2 12 0_2:14 0_2:14 Halo 2 12 0_2:14 0_2:14
+DEAL:1::Ghost 2 13 0_2:15 0_2:15 Halo 2 13 0_2:15 0_2:15
+DEAL:1::Ghost 2 14 0_2:16 0_2:16 Halo 2 14 0_2:16 0_2:16
+DEAL:1::Ghost 2 15 0_2:17 0_2:17 Halo 2 15 0_2:17 0_2:17
+DEAL:1::Ghost 2 20 0_2:24 0_2:24 Halo 2 20 0_2:24 0_2:24
+DEAL:1::Ghost 2 21 0_2:25 0_2:25 Halo 2 21 0_2:25 0_2:25
+DEAL:1::Ghost 2 22 0_2:26 0_2:26 Halo 2 22 0_2:26 0_2:26
+DEAL:1::Ghost 2 23 0_2:27 0_2:27 Halo 2 23 0_2:27 0_2:27
+DEAL:1::Ghost 2 28 0_2:34 0_2:34 Halo 2 28 0_2:34 0_2:34
+DEAL:1::Ghost 2 29 0_2:35 0_2:35 Halo 2 29 0_2:35 0_2:35
+DEAL:1::Ghost 2 30 0_2:36 0_2:36 Halo 2 30 0_2:36 0_2:36
+DEAL:1::Ghost 2 31 0_2:37 0_2:37 Halo 2 31 0_2:37 0_2:37
+
diff --git a/tests/deal.II/grid_tools_halo_layer_ghost_cells.mpirun=3.output b/tests/deal.II/grid_tools_halo_layer_ghost_cells.mpirun=3.output
new file mode 100644 (file)
index 0000000..97b0bef
--- /dev/null
@@ -0,0 +1,99 @@
+
+DEAL:0::dim = 2
+DEAL:0::Ghost 2 4 0_2:10 0_2:10 Halo 2 4 0_2:10 0_2:10
+DEAL:0::Ghost 2 6 0_2:12 0_2:12 Halo 2 6 0_2:12 0_2:12
+DEAL:0::Ghost 2 8 0_2:20 0_2:20 Halo 2 8 0_2:20 0_2:20
+DEAL:0::Ghost 2 9 0_2:21 0_2:21 Halo 2 9 0_2:21 0_2:21
+DEAL:0::Ghost 2 12 0_2:30 0_2:30 Halo 2 12 0_2:30 0_2:30
+DEAL:0::dim = 3
+DEAL:0::Ghost 2 24 0_2:30 0_2:30 Halo 2 24 0_2:30 0_2:30
+DEAL:0::Ghost 2 25 0_2:31 0_2:31 Halo 2 25 0_2:31 0_2:31
+DEAL:0::Ghost 2 26 0_2:32 0_2:32 Halo 2 26 0_2:32 0_2:32
+DEAL:0::Ghost 2 28 0_2:34 0_2:34 Halo 2 28 0_2:34 0_2:34
+DEAL:0::Ghost 2 29 0_2:35 0_2:35 Halo 2 29 0_2:35 0_2:35
+DEAL:0::Ghost 2 30 0_2:36 0_2:36 Halo 2 30 0_2:36 0_2:36
+DEAL:0::Ghost 2 32 0_2:40 0_2:40 Halo 2 32 0_2:40 0_2:40
+DEAL:0::Ghost 2 33 0_2:41 0_2:41 Halo 2 33 0_2:41 0_2:41
+DEAL:0::Ghost 2 34 0_2:42 0_2:42 Halo 2 34 0_2:42 0_2:42
+DEAL:0::Ghost 2 35 0_2:43 0_2:43 Halo 2 35 0_2:43 0_2:43
+DEAL:0::Ghost 2 40 0_2:50 0_2:50 Halo 2 40 0_2:50 0_2:50
+DEAL:0::Ghost 2 41 0_2:51 0_2:51 Halo 2 41 0_2:51 0_2:51
+DEAL:0::Ghost 2 42 0_2:52 0_2:52 Halo 2 42 0_2:52 0_2:52
+DEAL:0::Ghost 2 43 0_2:53 0_2:53 Halo 2 43 0_2:53 0_2:53
+DEAL:0::Ghost 2 48 0_2:60 0_2:60 Halo 2 48 0_2:60 0_2:60
+DEAL:0::Ghost 2 49 0_2:61 0_2:61 Halo 2 49 0_2:61 0_2:61
+DEAL:0::Ghost 2 50 0_2:62 0_2:62 Halo 2 50 0_2:62 0_2:62
+DEAL:0::Ghost 2 51 0_2:63 0_2:63 Halo 2 51 0_2:63 0_2:63
+DEAL:0::Ghost 2 56 0_2:70 0_2:70 Halo 2 56 0_2:70 0_2:70
+DEAL:0::Ghost 2 57 0_2:71 0_2:71 Halo 2 57 0_2:71 0_2:71
+DEAL:0::Ghost 2 58 0_2:72 0_2:72 Halo 2 58 0_2:72 0_2:72
+
+DEAL:1::dim = 2
+DEAL:1::Ghost 2 1 0_2:01 0_2:01 Halo 2 1 0_2:01 0_2:01
+DEAL:1::Ghost 2 2 0_2:02 0_2:02 Halo 2 2 0_2:02 0_2:02
+DEAL:1::Ghost 2 3 0_2:03 0_2:03 Halo 2 3 0_2:03 0_2:03
+DEAL:1::Ghost 2 12 0_2:30 0_2:30 Halo 2 12 0_2:30 0_2:30
+DEAL:1::Ghost 2 13 0_2:31 0_2:31 Halo 2 13 0_2:31 0_2:31
+DEAL:1::Ghost 2 14 0_2:32 0_2:32 Halo 2 14 0_2:32 0_2:32
+DEAL:1::dim = 3
+DEAL:1::Ghost 2 3 0_2:03 0_2:03 Halo 2 3 0_2:03 0_2:03
+DEAL:1::Ghost 2 4 0_2:04 0_2:04 Halo 2 4 0_2:04 0_2:04
+DEAL:1::Ghost 2 5 0_2:05 0_2:05 Halo 2 5 0_2:05 0_2:05
+DEAL:1::Ghost 2 6 0_2:06 0_2:06 Halo 2 6 0_2:06 0_2:06
+DEAL:1::Ghost 2 7 0_2:07 0_2:07 Halo 2 7 0_2:07 0_2:07
+DEAL:1::Ghost 2 10 0_2:12 0_2:12 Halo 2 10 0_2:12 0_2:12
+DEAL:1::Ghost 2 11 0_2:13 0_2:13 Halo 2 11 0_2:13 0_2:13
+DEAL:1::Ghost 2 12 0_2:14 0_2:14 Halo 2 12 0_2:14 0_2:14
+DEAL:1::Ghost 2 14 0_2:16 0_2:16 Halo 2 14 0_2:16 0_2:16
+DEAL:1::Ghost 2 15 0_2:17 0_2:17 Halo 2 15 0_2:17 0_2:17
+DEAL:1::Ghost 2 17 0_2:21 0_2:21 Halo 2 17 0_2:21 0_2:21
+DEAL:1::Ghost 2 19 0_2:23 0_2:23 Halo 2 19 0_2:23 0_2:23
+DEAL:1::Ghost 2 20 0_2:24 0_2:24 Halo 2 20 0_2:24 0_2:24
+DEAL:1::Ghost 2 21 0_2:25 0_2:25 Halo 2 21 0_2:25 0_2:25
+DEAL:1::Ghost 2 23 0_2:27 0_2:27 Halo 2 23 0_2:27 0_2:27
+DEAL:1::Ghost 2 40 0_2:50 0_2:50 Halo 2 40 0_2:50 0_2:50
+DEAL:1::Ghost 2 42 0_2:52 0_2:52 Halo 2 42 0_2:52 0_2:52
+DEAL:1::Ghost 2 43 0_2:53 0_2:53 Halo 2 43 0_2:53 0_2:53
+DEAL:1::Ghost 2 44 0_2:54 0_2:54 Halo 2 44 0_2:54 0_2:54
+DEAL:1::Ghost 2 46 0_2:56 0_2:56 Halo 2 46 0_2:56 0_2:56
+DEAL:1::Ghost 2 48 0_2:60 0_2:60 Halo 2 48 0_2:60 0_2:60
+DEAL:1::Ghost 2 49 0_2:61 0_2:61 Halo 2 49 0_2:61 0_2:61
+DEAL:1::Ghost 2 51 0_2:63 0_2:63 Halo 2 51 0_2:63 0_2:63
+DEAL:1::Ghost 2 52 0_2:64 0_2:64 Halo 2 52 0_2:64 0_2:64
+DEAL:1::Ghost 2 53 0_2:65 0_2:65 Halo 2 53 0_2:65 0_2:65
+DEAL:1::Ghost 2 56 0_2:70 0_2:70 Halo 2 56 0_2:70 0_2:70
+DEAL:1::Ghost 2 57 0_2:71 0_2:71 Halo 2 57 0_2:71 0_2:71
+DEAL:1::Ghost 2 58 0_2:72 0_2:72 Halo 2 58 0_2:72 0_2:72
+DEAL:1::Ghost 2 59 0_2:73 0_2:73 Halo 2 59 0_2:73 0_2:73
+DEAL:1::Ghost 2 60 0_2:74 0_2:74 Halo 2 60 0_2:74 0_2:74
+
+
+DEAL:2::dim = 2
+DEAL:2::Ghost 2 3 0_2:03 0_2:03 Halo 2 3 0_2:03 0_2:03
+DEAL:2::Ghost 2 6 0_2:12 0_2:12 Halo 2 6 0_2:12 0_2:12
+DEAL:2::Ghost 2 7 0_2:13 0_2:13 Halo 2 7 0_2:13 0_2:13
+DEAL:2::Ghost 2 9 0_2:21 0_2:21 Halo 2 9 0_2:21 0_2:21
+DEAL:2::Ghost 2 11 0_2:23 0_2:23 Halo 2 11 0_2:23 0_2:23
+DEAL:2::dim = 3
+DEAL:2::Ghost 2 5 0_2:05 0_2:05 Halo 2 5 0_2:05 0_2:05
+DEAL:2::Ghost 2 6 0_2:06 0_2:06 Halo 2 6 0_2:06 0_2:06
+DEAL:2::Ghost 2 7 0_2:07 0_2:07 Halo 2 7 0_2:07 0_2:07
+DEAL:2::Ghost 2 12 0_2:14 0_2:14 Halo 2 12 0_2:14 0_2:14
+DEAL:2::Ghost 2 13 0_2:15 0_2:15 Halo 2 13 0_2:15 0_2:15
+DEAL:2::Ghost 2 14 0_2:16 0_2:16 Halo 2 14 0_2:16 0_2:16
+DEAL:2::Ghost 2 15 0_2:17 0_2:17 Halo 2 15 0_2:17 0_2:17
+DEAL:2::Ghost 2 20 0_2:24 0_2:24 Halo 2 20 0_2:24 0_2:24
+DEAL:2::Ghost 2 21 0_2:25 0_2:25 Halo 2 21 0_2:25 0_2:25
+DEAL:2::Ghost 2 22 0_2:26 0_2:26 Halo 2 22 0_2:26 0_2:26
+DEAL:2::Ghost 2 23 0_2:27 0_2:27 Halo 2 23 0_2:27 0_2:27
+DEAL:2::Ghost 2 28 0_2:34 0_2:34 Halo 2 28 0_2:34 0_2:34
+DEAL:2::Ghost 2 29 0_2:35 0_2:35 Halo 2 29 0_2:35 0_2:35
+DEAL:2::Ghost 2 30 0_2:36 0_2:36 Halo 2 30 0_2:36 0_2:36
+DEAL:2::Ghost 2 31 0_2:37 0_2:37 Halo 2 31 0_2:37 0_2:37
+DEAL:2::Ghost 2 33 0_2:41 0_2:41 Halo 2 33 0_2:41 0_2:41
+DEAL:2::Ghost 2 34 0_2:42 0_2:42 Halo 2 34 0_2:42 0_2:42
+DEAL:2::Ghost 2 35 0_2:43 0_2:43 Halo 2 35 0_2:43 0_2:43
+DEAL:2::Ghost 2 37 0_2:45 0_2:45 Halo 2 37 0_2:45 0_2:45
+DEAL:2::Ghost 2 38 0_2:46 0_2:46 Halo 2 38 0_2:46 0_2:46
+DEAL:2::Ghost 2 39 0_2:47 0_2:47 Halo 2 39 0_2:47 0_2:47
+

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.