From: heister Date: Wed, 26 Oct 2011 16:25:12 +0000 (+0000) Subject: New: parallel::distributed::Triangulation::mesh_reconstruction_after_repartition... X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=2358ce244a3b001571f89744a55b466752864039;p=dealii-svn.git New: parallel::distributed::Triangulation::mesh_reconstruction_after_repartitioning setting which is necessary for save()/load() to be deterministic. git-svn-id: https://svn.dealii.org/trunk@24678 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/doc/news/changes.h b/deal.II/doc/news/changes.h index 240bd820f8..ff36636034 100644 --- a/deal.II/doc/news/changes.h +++ b/deal.II/doc/news/changes.h @@ -48,6 +48,13 @@ enabled due to a missing include file in file

Specific improvements

    + +
  1. New: parallel::distributed::Triangulation::mesh_reconstruction_after_repartitioning +setting which is necessary for save()/load() to be deterministic. Otherwise the matrix assembly +is done in a different order depending on the order of old refinements. +
    +(Timo Heister, 2011/10/26) +
  2. New: TriaAccessor<>::minimum_vertex_distance().
    (Timo Heister, 2011/10/25) diff --git a/deal.II/include/deal.II/distributed/tria.h b/deal.II/include/deal.II/distributed/tria.h index 62cba33b66..cd3e0e6de8 100644 --- a/deal.II/include/deal.II/distributed/tria.h +++ b/deal.II/include/deal.II/distributed/tria.h @@ -190,6 +190,30 @@ namespace parallel typedef typename dealii::Triangulation::cell_iterator cell_iterator; typedef typename dealii::Triangulation::raw_cell_iterator raw_cell_iterator; + /** + * Generic settings for distributed + * Triangulations. If + * mesh_reconstruction_after_repartitioning + * is set, the deal.II mesh will be + * reconstructed from the coarse mesh + * every time a repartioning in p4est + * happens. This can be a bit more + * expensive, but guarantees the same + * memory layout and therefore cell + * ordering in the deal.II mesh. As + * assembly is done in the deal.II + * cell ordering, this flag is + * required to get reproducable + * behaviour after snapshot/resume. + */ + enum Settings + { + default_setting = 0x0, + mesh_reconstruction_after_repartitioning = 0x1, + }; + + + /** * Constructor. * @@ -238,7 +262,8 @@ namespace parallel */ Triangulation (MPI_Comm mpi_communicator, const typename dealii::Triangulation::MeshSmoothing - smooth_grid = (dealii::Triangulation::none)); + smooth_grid = (dealii::Triangulation::none), + const Settings settings = default_setting); /** * Destructor. @@ -536,6 +561,11 @@ namespace parallel */ MPI_Comm mpi_communicator; + /** + * store the Settings. + */ + Settings settings; + /** * The subdomain id to be * used for the current diff --git a/deal.II/source/distributed/tria.cc b/deal.II/source/distributed/tria.cc index 12578ddc19..2638a1dda7 100644 --- a/deal.II/source/distributed/tria.cc +++ b/deal.II/source/distributed/tria.cc @@ -1796,7 +1796,8 @@ namespace parallel template Triangulation:: Triangulation (MPI_Comm mpi_communicator, - const typename dealii::Triangulation::MeshSmoothing smooth_grid) + const typename dealii::Triangulation::MeshSmoothing smooth_grid, + const Settings settings_) : // do not check // for distorted @@ -1806,6 +1807,7 @@ namespace parallel false), mpi_communicator (Utilities::MPI:: duplicate_communicator(mpi_communicator)), + settings(settings_), my_subdomain (Utilities::MPI::this_mpi_process (this->mpi_communicator)), triangulation_has_content (false), connectivity (0), @@ -2392,6 +2394,45 @@ namespace parallel this->smooth_grid = dealii::Triangulation::none; bool mesh_changed = false; + + // remove all deal.II refinements. Note that we could skip this and + // start from our current state, because the algorithm later coarsens as + // necessary. This has the advantage of being faster when large parts + // of the local partition changes (likely) and gives a deterministic + // ordering of the cells (useful for snapshot/resume). + // TODO: is there a more efficient way to do this? + if (settings & mesh_reconstruction_after_repartitioning) + while (this->begin_active()->level() > 0) + { + for (typename Triangulation::active_cell_iterator + cell = this->begin_active(); + cell != this->end(); + ++cell) + { + cell->set_coarsen_flag(); + } + + this->prepare_coarsening_and_refinement(); + const bool saved_refinement_in_progress = refinement_in_progress; + refinement_in_progress = true; + + try + { + this->execute_coarsening_and_refinement(); + } + catch (const typename Triangulation::DistortedCellList &) + { + // the underlying + // triangulation should not + // be checking for + // distorted cells + AssertThrow (false, ExcInternalError()); + } + + refinement_in_progress = saved_refinement_in_progress; + } + + // query p4est for the ghost cells typename dealii::internal::p4est::types::ghost * ghostlayer; ghostlayer = dealii::internal::p4est::functions::ghost_new (parallel_forest,