From: Jean-Paul Pelteret Date: Thu, 20 Aug 2015 15:24:43 +0000 (+0200) Subject: Implementation of general cell halo layer function inside GridTools. X-Git-Tag: v8.4.0-rc2~545^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=refs%2Fpull%2F1396%2Fhead;p=dealii.git Implementation of general cell halo layer function inside GridTools. 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. --- diff --git a/doc/news/changes.h b/doc/news/changes.h index de656a9d3e..2e6bdc88bd 100644 --- a/doc/news/changes.h +++ b/doc/news/changes.h @@ -191,6 +191,14 @@ inconvenience this causes.
    +
  1. 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. +
    + (Jean-Paul Pelteret, Denis Davydov, Wolfgang Bangerth, 2015/08/21) +
  2. +
  3. Documentation: How to set up a testsuite in a user project is now properly documented.
    diff --git a/include/deal.II/grid/filtered_iterator.h b/include/deal.II/grid/filtered_iterator.h index 7b4c804056..39ddaae813 100644 --- a/include/deal.II/grid/filtered_iterator.h +++ b/include/deal.II/grid/filtered_iterator.h @@ -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 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 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 + bool operator () (const Iterator &i) const; + + protected: + /** + * Stored value to compare the material id with. + */ + const std::set 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 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 + bool operator () (const Iterator &i) const; + + protected: + /** + * Stored value to compare the material id with. + */ + const std::set 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 material_ids, + const bool only_locally_owned) + : + material_ids (material_ids), + only_locally_owned (only_locally_owned) + {} + + + + template + 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 active_fe_indices, + const bool only_locally_owned) + : + active_fe_indices (active_fe_indices), + only_locally_owned (only_locally_owned) + {} + + + + template + 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(); + } + } diff --git a/include/deal.II/grid/grid_tools.h b/include/deal.II/grid/grid_tools.h index daf233effd..a3155cf56d 100644 --- a/include/deal.II/grid/grid_tools.h +++ b/include/deal.II/grid/grid_tools.h @@ -556,6 +556,80 @@ namespace GridTools get_active_neighbors (const typename Container::active_cell_iterator &cell, std::vector &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 + * bool + * pred_mat_id(const typename Triangulation::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); + * @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 + std::vector + compute_active_cell_halo_layer (const Container &container, + const std_cxx11::function &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 + std::vector + compute_ghost_cell_halo_layer (const Container &container); + + /*@}*/ /** * @name Partitions and subdomains of triangulations @@ -1280,6 +1354,7 @@ namespace GridTools namespace GridTools { + template void transform (const Predicate &predicate, Triangulation &triangulation) @@ -1429,7 +1504,6 @@ namespace GridTools - // declaration of explicit specializations template <> diff --git a/source/grid/grid_tools.cc b/source/grid/grid_tools.cc index 7a04c594b9..48cc49d8db 100644 --- a/source/grid/grid_tools.cc +++ b/source/grid/grid_tools.cc @@ -26,6 +26,7 @@ #include #include #include +#include #include #include #include @@ -1501,6 +1502,95 @@ next_cell: } + namespace + { + + template + bool + contains_locally_owned_cells (const std::vector &cells) + { + for (typename std::vector::const_iterator + it = cells.begin(); it != cells.end(); ++it) + { + if ((*it)->is_locally_owned()) + return true; + } + return false; + } + + template + bool + contains_artificial_cells (const std::vector &cells) + { + for (typename std::vector::const_iterator + it = cells.begin(); it != cells.end(); ++it) + { + if ((*it)->is_artificial()) + return true; + } + return false; + } + + } + + + + template + std::vector + compute_active_cell_halo_layer (const Container &container, + const std_cxx11::function &predicate) + { + std::vector active_halo_layer; + std::vector 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::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::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 + std::vector + compute_ghost_cell_halo_layer (const Container &container) + { + const std::vector + 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(active_halo_layer) == false, + ExcMessage("Halo layer contains locally owned cells")); + Assert(contains_artificial_cells(active_halo_layer) == false, + ExcMessage("Halo layer contains artificial cells")); + + return active_halo_layer; + } + + template void diff --git a/source/grid/grid_tools.inst.in b/source/grid/grid_tools.inst.in index c1c5fe585c..7123da054c 100644 --- a/source/grid/grid_tools.inst.in +++ b/source/grid/grid_tools.inst.in @@ -40,6 +40,15 @@ for (X : TRIANGULATION_AND_DOFHANDLERS; deal_II_dimension : DIMENSIONS ; deal_II const X &, const Point &); + template + std::vector::type> + compute_active_cell_halo_layer (const X &, + const std_cxx11::function::type&)> &); + + template + std::vector::type> + compute_ghost_cell_halo_layer (const X &); + template std::list > 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 &, const Point &); + template + std::vector>::type> + compute_active_cell_halo_layer (const parallel::distributed::Triangulation &, + const std_cxx11::function>::type&)> &); + + template + std::vector>::type> + compute_ghost_cell_halo_layer (const parallel::distributed::Triangulation &); + template std::list::cell_iterator, parallel::distributed::Triangulation::cell_iterator> > get_finest_common_cells (const parallel::distributed::Triangulation &mesh_1, @@ -190,8 +209,6 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS const hp::DoFHandler &, const Point &); - - template void get_subdomain_association (const Triangulation &, std::vector &); 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 index 0000000000..911899c39c --- /dev/null +++ b/tests/deal.II/grid_tools_halo_layer_01.cc @@ -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 +#include +#include +#include +#include + +template +bool +pred_mat_id(const typename Triangulation::active_cell_iterator & cell) +{ + return cell->material_id() == 2; +} + +template +void +write_mat_id_to_file (const Triangulation & tria) +{ + int count = 0; + typename Triangulation::active_cell_iterator + cell = tria.begin_active(), + endc = tria.end(); + for (; cell != endc; ++cell, ++count) + { + deallog + << count << " " + << static_cast(cell->material_id()) + << std::endl; + } + deallog << std::endl; +} + + +template +void test () +{ + deallog << "dim = " << dim << std::endl; + + Triangulation tria; + GridGenerator::hyper_cube(tria); + tria.refine_global(2); + + typedef typename Triangulation::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 active_halo_layer + = GridTools::compute_active_cell_halo_layer(tria, pred_mat_id); // General predicate + AssertThrow(active_halo_layer.size() > 0, ExcMessage("No halo layer found.")); + for (typename std::vector::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 index 0000000000..e91b7c208b --- /dev/null +++ b/tests/deal.II/grid_tools_halo_layer_01.output @@ -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 index 0000000000..8cf6744f6d --- /dev/null +++ b/tests/deal.II/grid_tools_halo_layer_02.cc @@ -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 +#include +#include +#include +#include +#include + +template +void +write_mat_id_to_file (const Triangulation & tria) +{ + int count = 0; + typename Triangulation::active_cell_iterator + cell = tria.begin_active(), + endc = tria.end(); + for (; cell != endc; ++cell, ++count) + { + deallog + << count << " " + << static_cast(cell->material_id()) + << std::endl; + } + deallog << std::endl; +} + + +template +void test () +{ + deallog << "dim = " << dim << std::endl; + + Triangulation tria; + GridGenerator::hyper_cube(tria); + tria.refine_global(2); + + typedef typename Triangulation::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 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::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 index 0000000000..e91b7c208b --- /dev/null +++ b/tests/deal.II/grid_tools_halo_layer_02.output @@ -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 index 0000000000..c2de9aef7c --- /dev/null +++ b/tests/deal.II/grid_tools_halo_layer_03.cc @@ -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 +#include +#include +#include +#include +#include + +template +void +write_mat_id_to_file (const Triangulation & tria) +{ + int count = 0; + typename Triangulation::active_cell_iterator + cell = tria.begin_active(), + endc = tria.end(); + for (; cell != endc; ++cell, ++count) + { + deallog + << count << " " + << static_cast(cell->material_id()) + << std::endl; + } + deallog << std::endl; +} + + +template +void test () +{ + deallog << "dim = " << dim << std::endl; + + Triangulation tria; + GridGenerator::hyper_cube(tria); + tria.refine_global(2); + + typedef typename Triangulation::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 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::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 index 0000000000..b8b799e6ff --- /dev/null +++ b/tests/deal.II/grid_tools_halo_layer_03.output @@ -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 index 0000000000..9db8a250ff --- /dev/null +++ b/tests/deal.II/grid_tools_halo_layer_04.cc @@ -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 +#include +#include +#include +#include +#include +#include + +template +void +write_active_fe_index_to_file (const hp::DoFHandler & dof_handler) +{ + int count = 0; + typename hp::DoFHandler::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 +void +write_vtk (const hp::DoFHandler & dof_handler, const std::string filename) +{ + Vector active_fe_index (dof_handler.get_tria().n_active_cells()); + int count = 0; + typename hp::DoFHandler::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 + data_component_interpretation + (1, DataComponentInterpretation::component_is_scalar); + const std::vector data_names (1, "active_fe_index"); + + DataOut > data_out; + data_out.attach_dof_handler (dof_handler); + data_out.add_data_vector (active_fe_index, data_names, + DataOut >::type_cell_data, + data_component_interpretation); + data_out.build_patches (); + + std::ofstream output (filename.c_str()); + data_out.write_vtk (output); +} + + +template +void test () +{ + deallog << "dim = " << dim << std::endl; + + Triangulation tria; + GridGenerator::hyper_cube(tria); + tria.refine_global(2); + hp::DoFHandler dof_handler (tria); + + typedef typename hp::DoFHandler::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 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::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 index 0000000000..e91b7c208b --- /dev/null +++ b/tests/deal.II/grid_tools_halo_layer_04.output @@ -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 index 0000000000..8edbba36f8 --- /dev/null +++ b/tests/deal.II/grid_tools_halo_layer_05.cc @@ -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 +#include +#include +#include +#include +#include +#include + +template +void +write_active_fe_index_to_file (const hp::DoFHandler & dof_handler) +{ + int count = 0; + typename hp::DoFHandler::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 +void +write_vtk (const hp::DoFHandler & dof_handler, const std::string filename) +{ + Vector active_fe_index (dof_handler.get_tria().n_active_cells()); + int count = 0; + typename hp::DoFHandler::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 + data_component_interpretation + (1, DataComponentInterpretation::component_is_scalar); + const std::vector data_names (1, "active_fe_index"); + + DataOut > data_out; + data_out.attach_dof_handler (dof_handler); + data_out.add_data_vector (active_fe_index, data_names, + DataOut >::type_cell_data, + data_component_interpretation); + data_out.build_patches (); + + std::ofstream output (filename.c_str()); + data_out.write_vtk (output); +} + + +template +void test () +{ + deallog << "dim = " << dim << std::endl; + + Triangulation tria; + GridGenerator::hyper_cube(tria); + tria.refine_global(2); + hp::DoFHandler dof_handler (tria); + + typedef typename hp::DoFHandler::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 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::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 index 0000000000..b8b799e6ff --- /dev/null +++ b/tests/deal.II/grid_tools_halo_layer_05.output @@ -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 index 0000000000..f8fc22f479 --- /dev/null +++ b/tests/deal.II/grid_tools_halo_layer_ghost_cells.cc @@ -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 +#include +#include +#include +#include +#include + +#include + +template +void test () +{ + const MPI_Comm &mpi_communicator = MPI_COMM_WORLD; + deallog << "dim = " << dim << std::endl; + + parallel::distributed::Triangulation 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::active_cell_iterator cell_iterator; + + // Mark a small block at the corner of the hypercube + std::vector 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 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::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 index 0000000000..a50adba60c --- /dev/null +++ b/tests/deal.II/grid_tools_halo_layer_ghost_cells.mpirun=2.output @@ -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 index 0000000000..97b0bef1a4 --- /dev/null +++ b/tests/deal.II/grid_tools_halo_layer_ghost_cells.mpirun=3.output @@ -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 +