]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
New: parallel::distributed::Triangulation<dim>::mesh_reconstruction_after_repartition...
authorheister <heister@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 26 Oct 2011 16:25:12 +0000 (16:25 +0000)
committerheister <heister@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 26 Oct 2011 16:25:12 +0000 (16:25 +0000)
git-svn-id: https://svn.dealii.org/trunk@24678 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/doc/news/changes.h
deal.II/include/deal.II/distributed/tria.h
deal.II/source/distributed/tria.cc

index 240bd820f81b655a1e67521502332b1ea645286a..ff36636034dff29dbb8af4b390eb0fc3f6664989 100644 (file)
@@ -48,6 +48,13 @@ enabled due to a missing include file in file
 <h3>Specific improvements</h3>
 
 <ol>
+
+<li> New: parallel::distributed::Triangulation<dim>::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.
+<br>
+(Timo Heister, 2011/10/26)
+
 <li> New: TriaAccessor<>::minimum_vertex_distance().
 <br>
 (Timo Heister, 2011/10/25)
index 62cba33b66af4d1ba2d9d7a365871cee6a82f926..cd3e0e6de8f06cd61cb8e9649d126762799e1684 100644 (file)
@@ -190,6 +190,30 @@ namespace parallel
        typedef typename dealii::Triangulation<dim,spacedim>::cell_iterator        cell_iterator;
        typedef typename dealii::Triangulation<dim,spacedim>::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<dim,spacedim>::MeshSmoothing
-                      smooth_grid = (dealii::Triangulation<dim,spacedim>::none));
+                      smooth_grid = (dealii::Triangulation<dim,spacedim>::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
index 12578ddc1933dc071cc9bf34ea8adee20fb53c71..2638a1dda7b4da80f0594620403adba66fd344a6 100644 (file)
@@ -1796,7 +1796,8 @@ namespace parallel
     template <int dim, int spacedim>
     Triangulation<dim,spacedim>::
     Triangulation (MPI_Comm mpi_communicator,
-                  const typename dealii::Triangulation<dim,spacedim>::MeshSmoothing smooth_grid)
+                  const typename dealii::Triangulation<dim,spacedim>::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<dim,spacedim>::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<dim, spacedim>::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<dim, spacedim>::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<dim>::ghost * ghostlayer;
       ghostlayer = dealii::internal::p4est::functions<dim>::ghost_new (parallel_forest,

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.