From cf465f79e0300d43292803c2d4e00bd746cb1121 Mon Sep 17 00:00:00 2001 From: heister Date: Fri, 14 Oct 2011 15:57:21 +0000 Subject: [PATCH] fix bug in Triangulation::load() and improve documentation git-svn-id: https://svn.dealii.org/trunk@24602 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/include/deal.II/distributed/tria.h | 26 ++++- deal.II/source/distributed/tria.cc | 113 ++++++++++++--------- 2 files changed, 86 insertions(+), 53 deletions(-) diff --git a/deal.II/include/deal.II/distributed/tria.h b/deal.II/include/deal.II/distributed/tria.h index 1208278c45..62cba33b66 100644 --- a/deal.II/include/deal.II/distributed/tria.h +++ b/deal.II/include/deal.II/distributed/tria.h @@ -603,10 +603,25 @@ namespace parallel /** - * + * number of bytes that get attached to the Triangulation + * through register_data_attach() for example + * SolutionTransfer. */ unsigned int attached_data_size; + + /** + * number of functions that get attached to the Triangulation + * through register_data_attach() for example + * SolutionTransfer. + */ unsigned int n_attached_datas; + + /** + * number of functions that need to unpack their data + * after a call from load() + */ + unsigned int n_attached_deserialize; + typedef std_cxx1x::function< void(typename Triangulation::cell_iterator, CellStatus, void*) > pack_callback_t; @@ -614,12 +629,15 @@ namespace parallel typedef std::pair callback_pair_t; typedef std::list callback_list_t; + + /** + * List of callback functions registered by + * register_data_attach() that are going to be called + * for packing data. + */ callback_list_t attached_data_pack_callbacks; - - - /** * Two arrays that store which p4est * tree corresponds to which coarse diff --git a/deal.II/source/distributed/tria.cc b/deal.II/source/distributed/tria.cc index 63db8611f3..12578ddc19 100644 --- a/deal.II/source/distributed/tria.cc +++ b/deal.II/source/distributed/tria.cc @@ -1812,7 +1812,8 @@ namespace parallel parallel_forest (0), refinement_in_progress (false), attached_data_size(0), - n_attached_datas(0) + n_attached_datas(0), + n_attached_deserialize(0) { // initialize p4est. do this in a separate function since it has // to happen only once, even if we have triangulation objects @@ -1953,6 +1954,8 @@ namespace parallel Triangulation:: save(const char* filename) const { + Assert(n_attached_deserialize==0, + ExcMessage ("not all SolutionTransfer's got deserialized after the last load()")); int real_data_size = 0; if (attached_data_size>0) real_data_size = attached_data_size+sizeof(CellStatus); @@ -2009,6 +2012,7 @@ namespace parallel attached_data_size = 0; n_attached_datas = 0; + n_attached_deserialize = attached_count; parallel_forest = dealii::internal::p4est::functions::load ( filename, mpi_communicator, @@ -2840,61 +2844,72 @@ namespace parallel const CellStatus, const void*)> & unpack_callback) { - Assert(offset < attached_data_size, ExcMessage("invalid offset in notify_ready_to_unpack()")); - Assert(n_attached_datas>0, ExcMessage("notify_ready_to_unpack() called too often")); - - // Recurse over p4est and hand the caller the data back - for (typename Triangulation::cell_iterator - cell = this->begin(0); - cell != this->end(0); - ++cell) - { - //skip coarse cells, that are not ours - if (tree_exists_locally(parallel_forest, - coarse_cell_to_p4est_tree_permutation[cell->index()]) - == false) - continue; - - typename dealii::internal::p4est::types::quadrant p4est_coarse_cell; - typename dealii::internal::p4est::types::tree *tree = - init_tree(cell->index()); - - dealii::internal::p4est::init_coarse_quadrant (p4est_coarse_cell); - - // parent_cell is not correct here, - // but is only used in a refined - // cell - post_mesh_data_recursively(*tree, - cell, - cell, - p4est_coarse_cell, - offset, - unpack_callback); - } + Assert (offset < attached_data_size, ExcMessage ("invalid offset in notify_ready_to_unpack()")); + Assert (n_attached_datas > 0, ExcMessage ("notify_ready_to_unpack() called too often")); + + // Recurse over p4est and hand the caller the data back + for (typename Triangulation::cell_iterator + cell = this->begin (0); + cell != this->end (0); + ++cell) + { + //skip coarse cells, that are not ours + if (tree_exists_locally (parallel_forest, + coarse_cell_to_p4est_tree_permutation[cell->index() ]) + == false) + continue; + + typename dealii::internal::p4est::types::quadrant p4est_coarse_cell; + typename dealii::internal::p4est::types::tree *tree = + init_tree (cell->index()); + + dealii::internal::p4est::init_coarse_quadrant (p4est_coarse_cell); + + // parent_cell is not correct here, + // but is only used in a refined + // cell + post_mesh_data_recursively (*tree, + cell, + cell, + p4est_coarse_cell, + offset, + unpack_callback); + } --n_attached_datas; + if (n_attached_deserialize > 0) + { + --n_attached_deserialize; + attached_data_pack_callbacks.pop_front(); + } - if (!n_attached_datas) - { - // everybody got his data, time for cleanup! - attached_data_size=0; - attached_data_pack_callbacks.clear(); - - // and release the data - void * userptr = parallel_forest->user_pointer; - dealii::internal::p4est::functions::reset_data (parallel_forest, 0, NULL, NULL); - parallel_forest->user_pointer = userptr; - } + // important: only remove data if we are not in the deserialization + // process. There, each SolutionTransfer registers and unpacks + // before the next one does this, so n_attached_datas is only 1 here. + // This would destroy the saved data before the second SolutionTransfer + // can get it. This created a bug that is documented in + // tests/mpi/p4est_save_03 with more than one SolutionTransfer. + if (!n_attached_datas && n_attached_deserialize == 0) + { + // everybody got his data, time for cleanup! + attached_data_size = 0; + attached_data_pack_callbacks.clear(); + + // and release the data + void * userptr = parallel_forest->user_pointer; + dealii::internal::p4est::functions::reset_data (parallel_forest, 0, NULL, NULL); + parallel_forest->user_pointer = userptr; + } } - template - const std::vector & - Triangulation::get_p4est_tree_to_coarse_cell_permutation() const - { - return p4est_tree_to_coarse_cell_permutation; - } + template + const std::vector & + Triangulation::get_p4est_tree_to_coarse_cell_permutation() const + { + return p4est_tree_to_coarse_cell_permutation; + } -- 2.39.5