From: Reza Rastak Date: Wed, 7 Aug 2019 16:24:23 +0000 (-0600) Subject: Cell id is the key X-Git-Tag: v9.2.0-rc1~1278^2~1 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=b087ed798a29e83f584b66e7415db579e67bb250;p=dealii.git Cell id is the key --- diff --git a/doc/news/changes/minor/20190807RezaRastakPeterMunch b/doc/news/changes/minor/20190807RezaRastakPeterMunch new file mode 100644 index 0000000000..f0e66adda5 --- /dev/null +++ b/doc/news/changes/minor/20190807RezaRastakPeterMunch @@ -0,0 +1,4 @@ +Fixed: Bug in CellDataStorage that produced error when applying adaptive +refinement in multi-material grids. +
+(Reza Rastak, Peter Munch, 2019/08/07) diff --git a/include/deal.II/base/quadrature_point_data.h b/include/deal.II/base/quadrature_point_data.h index 8656eac401..758a0a698f 100644 --- a/include/deal.II/base/quadrature_point_data.h +++ b/include/deal.II/base/quadrature_point_data.h @@ -160,9 +160,31 @@ public: private: /** - * A map to store a vector of data on a cell. + * Depending on @p CellIteratorType we need to define a triangulation type */ - std::map>> map; + using triangulation_type = + Triangulation; + + /** + * The key is used as identified for a cell data structure. + * We need to use CellId as part of the key because it remains unique during + * adaptive refinement. We need to use a pointer to the triangulation in the + * key because we want to allow the data storage class to keep information on + * multiple meshes simultaneously. + */ + using key_type = std::pair; + + /** + * A map to store a vector of data on each cell. + */ + std::map>> map; + + /** + * Computes the key from a cell iterator + */ + static key_type + get_key(const CellIteratorType &cell); /** * @addtogroup Exceptions @@ -522,6 +544,14 @@ namespace parallel // CellDataStorage //-------------------------------------------------------------------- +template +typename CellDataStorage::key_type +CellDataStorage::get_key( + const CellIteratorType &cell) +{ + return {&cell->get_triangulation(), cell->id()}; +} + template template void @@ -532,13 +562,14 @@ CellDataStorage::initialize( static_assert(std::is_base_of::value, "User's T class should be derived from user's DataType class"); - if (map.find(cell) == map.end()) + const auto key = get_key(cell); + if (map.find(key) == map.end()) { - map[cell] = std::vector>(n_q_points); + map[key] = std::vector>(n_q_points); // we need to initialize one-by-one as the std::vector<>(q, T()) // will end with a single same T object stored in each element of the // vector: - auto it = map.find(cell); + const auto it = map.find(key); for (unsigned int q = 0; q < n_q_points; q++) it->second[q] = std::make_shared(); } @@ -565,7 +596,8 @@ template bool CellDataStorage::erase(const CellIteratorType &cell) { - const auto it = map.find(cell); + const auto key = get_key(cell); + const auto it = map.find(key); if (it == map.end()) return false; for (unsigned int i = 0; i < it->second.size(); i++) @@ -576,7 +608,7 @@ CellDataStorage::erase(const CellIteratorType &cell) "Can not erase the cell data multiple objects reference its data.")); } - return (map.erase(cell) == 1); + return (map.erase(key) == 1); } @@ -614,7 +646,7 @@ CellDataStorage::get_data( static_assert(std::is_base_of::value, "User's T class should be derived from user's DataType class"); - auto it = map.find(cell); + const auto it = map.find(get_key(cell)); Assert(it != map.end(), ExcMessage("Could not find data for the cell")); // It would be nice to have a specialized version of this function for @@ -642,7 +674,7 @@ CellDataStorage::get_data( static_assert(std::is_base_of::value, "User's T class should be derived from user's DataType class"); - auto it = map.find(cell); + const auto it = map.find(get_key(cell)); Assert(it != map.end(), ExcMessage("Could not find QP data for the cell")); // Cast base class to the desired class. This has to be done irrespectively of diff --git a/tests/base/cell_data_storage_03.cc b/tests/base/cell_data_storage_03.cc new file mode 100644 index 0000000000..27c38493e7 --- /dev/null +++ b/tests/base/cell_data_storage_03.cc @@ -0,0 +1,128 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2019 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE.md at +// the top level directory of deal.II. +// +// --------------------------------------------------------------------- + +// It checks whether the cell_data_storage, specifically the methods +// initialize() and get_data(), works with two different materials +// types storing two different data structures. +// It makes sure the results are not affected during local refinement and +// and coarsening. + +#include + +#include + +#include "../tests.h" + +// history structs for each material +struct MaterialBase +{ + virtual ~MaterialBase() + {} +}; +struct Mat0 : MaterialBase +{ + double x{5.}; +}; +struct Mat1 : MaterialBase +{ + double y{8.}; +}; +constexpr unsigned int n_data_points_per_cell = 4; + + +void setup_history( + Triangulation<2> &tria, + CellDataStorage::active_cell_iterator, MaterialBase> + &storage) +{ + deallog << "initializing history" << std::endl; + for (auto cell : tria.active_cell_iterators()) + { + deallog << "initializing cell " << cell->id() << std::endl; + if (cell->material_id() == 0) + storage.template initialize(cell, n_data_points_per_cell); + else if (cell->material_id() == 1) + storage.template initialize(cell, n_data_points_per_cell); + } +} + +void read_history( + Triangulation<2> &tria, + CellDataStorage::active_cell_iterator, MaterialBase> + &storage) +{ + deallog << "reading history" << std::endl; + for (auto cell : tria.active_cell_iterators()) + { + if (cell->material_id() == 0) + { + auto data_vec = storage.template get_data(cell); + deallog << "Cell with material id = 0 contains the data " + << data_vec[0]->x << std::endl; + } + else if (cell->material_id() == 1) + { + auto data_vec = storage.template get_data(cell); + deallog << "Cell with material id = 1 contains the data " + << data_vec[0]->y << std::endl; + } + } +} + +int +main() +{ + initlog(); + + // create a mesh with two cells + Triangulation<2> tria; + GridGenerator::subdivided_hyper_rectangle(tria, + /*repetitions*/ {2, 1}, + /*point1*/ {0., 0.}, + /*point2*/ {2., 1.}); + + // assign material ids + auto cell0 = tria.begin_active(); + cell0->set_material_id(0); + auto cell1 = cell0; + ++cell1; + cell1->set_material_id(1); + + // refine one cell + cell0->set_refine_flag(); + tria.execute_coarsening_and_refinement(); + deallog << "first refinement done" << std::endl; + + // create a history structure and populate it + CellDataStorage::active_cell_iterator, MaterialBase> + storage; + setup_history(tria, storage); + read_history(tria, storage); + + // coarsen the cell and refine the other cell + for (auto cell : tria.cell_iterators_on_level(1)) + cell->set_coarsen_flag(); + cell1->set_refine_flag(); + tria.execute_coarsening_and_refinement(); + deallog << "second refinement done" << std::endl; + + // initialize history for newly created cells + setup_history(tria, storage); + read_history(tria, storage); + + deallog << "OK" << std::endl; + + return 0; +} diff --git a/tests/base/cell_data_storage_03.output b/tests/base/cell_data_storage_03.output new file mode 100644 index 0000000000..36297befb6 --- /dev/null +++ b/tests/base/cell_data_storage_03.output @@ -0,0 +1,28 @@ + +DEAL::first refinement done +DEAL::initializing history +DEAL::initializing cell 1_0: +DEAL::initializing cell 0_1:0 +DEAL::initializing cell 0_1:1 +DEAL::initializing cell 0_1:2 +DEAL::initializing cell 0_1:3 +DEAL::reading history +DEAL::Cell with material id = 1 contains the data 8.00000 +DEAL::Cell with material id = 0 contains the data 5.00000 +DEAL::Cell with material id = 0 contains the data 5.00000 +DEAL::Cell with material id = 0 contains the data 5.00000 +DEAL::Cell with material id = 0 contains the data 5.00000 +DEAL::second refinement done +DEAL::initializing history +DEAL::initializing cell 0_0: +DEAL::initializing cell 1_1:0 +DEAL::initializing cell 1_1:1 +DEAL::initializing cell 1_1:2 +DEAL::initializing cell 1_1:3 +DEAL::reading history +DEAL::Cell with material id = 0 contains the data 5.00000 +DEAL::Cell with material id = 1 contains the data 8.00000 +DEAL::Cell with material id = 1 contains the data 8.00000 +DEAL::Cell with material id = 1 contains the data 8.00000 +DEAL::Cell with material id = 1 contains the data 8.00000 +DEAL::OK