<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)
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.
*
*/
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.
*/
MPI_Comm mpi_communicator;
+ /**
+ * store the Settings.
+ */
+ Settings settings;
+
/**
* The subdomain id to be
* used for the current
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
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),
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,