16#include <deal.II/base/mpi.templates.h>
35#ifdef DEAL_II_WITH_MPI
40 template <
int dim,
int spacedim>
54 const auto partition_settings =
58 (void)partition_settings;
64 ExcMessage(
"Settings must contain exactly one type of the active "
65 "cell partitioning scheme."));
69 ExcMessage(
"construct_multigrid_hierarchy requires "
70 "allow_artificial_cells to be set to true."));
75 template <
int dim,
int spacedim>
85 template <
int dim,
int spacedim>
92 const unsigned int max_active_cells =
97 "A parallel::shared::Triangulation needs to be refined in the same "
98 "way on all processors, but the participating processors don't "
99 "agree on the number of active cells."));
106# ifdef DEAL_II_TRILINOS_WITH_ZOLTAN
108# elif defined DEAL_II_WITH_METIS
116# ifndef DEAL_II_TRILINOS_WITH_ZOLTAN
119 "Choosing 'partition_zoltan' requires the library "
120 "to be compiled with support for Zoltan! "
121 "Instead, you might use 'partition_auto' to select "
122 "a partitioning algorithm that is supported "
123 "by your current configuration."));
131# ifndef DEAL_II_WITH_METIS
134 "Choosing 'partition_metis' requires the library "
135 "to be compiled with support for METIS! "
136 "Instead, you might use 'partition_auto' to select "
137 "a partitioning algorithm that is supported "
138 "by your current configuration."));
184 active_halo_layer_vector =
187 std::set<typename parallel::shared::Triangulation<dim, spacedim>::
188 active_cell_iterator>
189 active_halo_layer(active_halo_layer_vector.begin(),
190 active_halo_layer_vector.end());
192 for (
unsigned int index = 0; cell != endc; cell++, index++)
197 if (cell->is_locally_owned() ==
false &&
198 active_halo_layer.find(cell) == active_halo_layer.end())
211 for (
unsigned int lvl = 0; lvl < this->
n_levels(); ++lvl)
219 level_halo_layer_vector =
221 *
this, predicate, lvl);
222 std::set<
typename parallel::shared::
223 Triangulation<dim, spacedim>::cell_iterator>
224 level_halo_layer(level_halo_layer_vector.begin(),
225 level_halo_layer_vector.end());
228 cell_iterator cell = this->
begin(lvl),
229 endc = this->
end(lvl);
230 for (
unsigned int index = 0; cell != endc; cell++, index++)
235 cell->level_subdomain_id();
244 if (cell->is_active() &&
245 cell->subdomain_id() !=
250 if (cell->has_children())
252 bool keep_cell =
false;
253 for (
unsigned int c = 0;
254 c < GeometryInfo<dim>::max_children_per_cell;
256 if (cell->child(c)->level_subdomain_id() ==
268 if (!cell->is_locally_owned_on_level() &&
269 level_halo_layer.find(cell) != level_halo_layer.
end())
273 cell->set_level_subdomain_id(
282 for (
unsigned int index = 0; cell != endc; cell++,
index++)
283 true_subdomain_ids_of_cells[
index] = cell->subdomain_id();
289 const unsigned int n_my_cells = std::count_if(
293 [](
const auto &i) {
return (i.is_locally_owned()); });
295 const unsigned int total_cells =
298 ExcMessage(
"Not all cells are assigned to a processor."));
303 if (settings & construct_multigrid_hierarchy)
305 const unsigned int n_my_cells =
306 std::count_if(this->
begin(), this->
end(), [](
const auto &i) {
307 return (i.is_locally_owned_on_level());
311 const unsigned int total_cells =
314 ExcMessage(
"Not all cells are assigned to a processor."));
321 template <
int dim,
int spacedim>
330 template <
int dim,
int spacedim>
332 const std::vector<types::subdomain_id>
340 template <
int dim,
int spacedim>
342 const std::vector<types::subdomain_id>
344 const unsigned int level)
const
356 template <
int dim,
int spacedim>
372 using int_type = std::underlying_type_t<
374 static_assert(
sizeof(int_type) ==
sizeof(std::uint8_t),
375 "Internal type mismatch.");
377 std::vector<int_type> refinement_configurations(this->
n_active_cells() *
381 if (cell->is_locally_owned())
383 refinement_configurations[cell->active_cell_index() * 2 + 0] =
384 static_cast<int_type
>(cell->refine_flag_set());
385 refinement_configurations[cell->active_cell_index() * 2 + 1] =
386 static_cast<int_type
>(cell->coarsen_flag_set() ? 1 : 0);
390 this->get_communicator(),
391 refinement_configurations);
395 cell->clear_refine_flag();
396 cell->clear_coarsen_flag();
399 (refinement_configurations[cell->active_cell_index() * 2 + 0] >
403 refinement_configurations[cell->active_cell_index() * 2 +
407 "Refinement/coarsening flags of cells are not consistent in parallel!"));
409 if (refinement_configurations[cell->active_cell_index() * 2 + 0] !=
412 refinement_configurations[cell->active_cell_index() * 2 + 0]));
414 if (refinement_configurations[cell->active_cell_index() * 2 + 1] >
416 cell->set_coarsen_flag();
422 this->update_number_cache();
427 template <
int dim,
int spacedim>
440 const typename ::Triangulation<dim, spacedim>::DistortedCellList
453 template <
int dim,
int spacedim>
459 (void)construction_data;
466 template <
int dim,
int spacedim>
473 const ::parallel::DistributedTriangulationBase<dim, spacedim>
474 *
>(&other_tria) ==
nullptr),
476 "Cannot use this function on parallel::distributed::Triangulation."));
492 template <
int dim,
int spacedim>
502 template <
int dim,
int spacedim>
512 template <
int dim,
int spacedim>
514 const std::vector<unsigned int>
518 return true_subdomain_ids_of_cells;
523 template <
int dim,
int spacedim>
525 const std::vector<unsigned int>
527 const unsigned int)
const
530 return true_level_subdomain_ids_of_cells;
546 template <
int dim,
int spacedim>
558 const std::vector<types::subdomain_id> &true_subdomain_ids =
559 shared_tria->get_true_subdomain_ids_of_cells();
561 saved_subdomain_ids.resize(shared_tria->n_active_cells());
562 for (const auto &cell : shared_tria->active_cell_iterators())
564 const unsigned int index = cell->active_cell_index();
565 saved_subdomain_ids[index] = cell->subdomain_id();
566 cell->set_subdomain_id(true_subdomain_ids[index]);
573 template <
int dim,
int spacedim>
580 for (
const auto &cell :
shared_tria->active_cell_iterators())
582 const unsigned int index = cell->active_cell_index();
593#include "shared_tria.inst"
cell_iterator begin(const unsigned int level=0) const
unsigned int n_active_cells() const
cell_iterator end() const
unsigned int n_cells() const
active_cell_iterator begin_active(const unsigned int level=0) const
cell_iterator begin(const unsigned int level=0) const
std::vector< Point< spacedim > > vertices
unsigned int n_active_cells() const
unsigned int n_levels() const
cell_iterator end() const
unsigned int n_cells() const
MeshSmoothing smooth_grid
active_cell_iterator begin_active(const unsigned int level=0) const
std::vector< unsigned int > saved_subdomain_ids
const SmartPointer< const ::parallel::shared::Triangulation< dim, spacedim > > shared_tria
TemporarilyRestoreSubdomainIds(const Triangulation< dim, spacedim > &tria)
types::subdomain_id n_subdomains
const MPI_Comm mpi_communicator
virtual void update_number_cache()
TriangulationBase(const MPI_Comm mpi_communicator, const typename ::Triangulation< dim, spacedim >::MeshSmoothing smooth_grid=(::Triangulation< dim, spacedim >::none), const bool check_for_distorted_cells=false)
virtual MPI_Comm get_communicator() const override
types::subdomain_id my_subdomain
virtual void copy_triangulation(const ::Triangulation< dim, spacedim > &old_tria) override
virtual void execute_coarsening_and_refinement() override
std::vector< types::subdomain_id > true_subdomain_ids_of_cells
typename ::Triangulation< dim, spacedim >::cell_iterator cell_iterator
virtual void copy_triangulation(const ::Triangulation< dim, spacedim > &other_tria) override
const std::vector< types::subdomain_id > & get_true_subdomain_ids_of_cells() const
std::vector< std::vector< types::subdomain_id > > true_level_subdomain_ids_of_cells
@ partition_custom_signal
@ construct_multigrid_hierarchy
bool with_artificial_cells() const
const std::vector< types::subdomain_id > & get_true_level_subdomain_ids_of_cells(const unsigned int level) const
Triangulation(const MPI_Comm mpi_communicator, const typename ::Triangulation< dim, spacedim >::MeshSmoothing=(::Triangulation< dim, spacedim >::none), const bool allow_artificial_cells=false, const Settings settings=partition_auto)
typename ::Triangulation< dim, spacedim >::active_cell_iterator active_cell_iterator
virtual bool is_multilevel_hierarchy_constructed() const override
const bool allow_artificial_cells
virtual void create_triangulation(const std::vector< Point< spacedim > > &vertices, const std::vector< CellData< dim > > &cells, const SubCellData &subcelldata) override
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_CXX20_REQUIRES(condition)
#define DEAL_II_NAMESPACE_CLOSE
IteratorRange< active_cell_iterator > active_cell_iterators() const
#define Assert(cond, exc)
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)
#define DEAL_II_ASSERT_UNREACHABLE()
#define DEAL_II_NOT_IMPLEMENTED()
T sum(const T &t, const MPI_Comm mpi_communicator)
T max(const T &t, const MPI_Comm mpi_communicator)
const types::subdomain_id artificial_subdomain_id