16 #include <deal.II/base/mpi.h> 17 #include <deal.II/base/utilities.h> 19 #include <deal.II/distributed/shared_tria.h> 20 #include <deal.II/distributed/tria.h> 22 #include <deal.II/grid/filtered_iterator.h> 23 #include <deal.II/grid/grid_tools.h> 24 #include <deal.II/grid/tria.h> 25 #include <deal.II/grid/tria_accessor.h> 26 #include <deal.II/grid/tria_iterator.h> 28 #include <deal.II/lac/sparsity_tools.h> 31 DEAL_II_NAMESPACE_OPEN
33 #ifdef DEAL_II_WITH_MPI 38 template <
int dim,
int spacedim>
40 MPI_Comm mpi_communicator,
41 const typename ::Triangulation<dim, spacedim>::MeshSmoothing
43 const bool allow_artificial_cells,
49 , allow_artificial_cells(allow_artificial_cells)
51 const auto partition_settings =
55 (void)partition_settings;
63 "Settings must contain exactly one type of the active cell partitioning scheme."))
66 allow_artificial_cells,
68 "construct_multigrid_hierarchy requires allow_artificial_cells to be set to true."))
73 template <
int dim,
int spacedim>
80 const unsigned int max_active_cells =
85 "A parallel::shared::Triangulation needs to be refined in the same" 86 "way on all processors, but the participating processors don't " 87 "agree on the number of active cells."))
94 # ifdef DEAL_II_TRILINOS_WITH_ZOLTAN 96 # elif defined DEAL_II_WITH_METIS 104 # ifndef DEAL_II_TRILINOS_WITH_ZOLTAN 107 "Choosing 'partition_zoltan' requires the library " 108 "to be compiled with support for Zoltan! " 109 "Instead, you might use 'partition_auto' to select " 110 "a partitioning algorithm that is supported " 111 "by your current configuration."));
119 # ifndef DEAL_II_WITH_METIS 122 "Choosing 'partition_metis' requires the library " 123 "to be compiled with support for METIS! " 124 "Instead, you might use 'partition_auto' to select " 125 "a partitioning algorithm that is supported " 126 "by your current configuration."));
150 ::GridTools::partition_multigrid_levels(*
this);
156 spacedim>::active_cell_iterator
166 active_cell_iterator &)>
171 spacedim>::active_cell_iterator>
172 active_halo_layer_vector =
173 ::GridTools::compute_active_cell_halo_layer(*
this,
175 std::set<typename parallel::shared::Triangulation<dim, spacedim>::
176 active_cell_iterator>
177 active_halo_layer(active_halo_layer_vector.begin(),
178 active_halo_layer_vector.end());
180 for (
unsigned int index = 0; cell != endc; cell++, index++)
185 if (cell->is_locally_owned() ==
false &&
186 active_halo_layer.find(cell) == active_halo_layer.end())
191 if (
settings & construct_multigrid_hierarchy)
199 for (
unsigned int lvl = 0; lvl < this->
n_levels(); ++lvl)
206 spacedim>::cell_iterator>
207 level_halo_layer_vector =
208 ::GridTools::compute_cell_halo_layer_on_level(
209 *
this, predicate, lvl);
210 std::set<
typename parallel::shared::
211 Triangulation<dim, spacedim>::cell_iterator>
212 level_halo_layer(level_halo_layer_vector.begin(),
213 level_halo_layer_vector.end());
216 cell_iterator cell = this->
begin(lvl),
217 endc = this->
end(lvl);
218 for (
unsigned int index = 0; cell != endc; cell++, index++)
223 cell->level_subdomain_id();
232 if (!cell->has_children() &&
233 cell->subdomain_id() !=
238 if (cell->has_children())
240 bool keep_cell =
false;
241 for (
unsigned int c = 0;
242 c < GeometryInfo<dim>::max_children_per_cell;
244 if (cell->child(c)->level_subdomain_id() ==
256 if (!cell->is_locally_owned_on_level() &&
257 level_halo_layer.find(cell) != level_halo_layer.end())
261 cell->set_level_subdomain_id(
270 for (
unsigned int index = 0; cell != endc; cell++, index++)
277 unsigned int n_my_cells = 0;
279 spacedim>::active_cell_iterator
282 for (; cell != endc; ++cell)
283 if (cell->is_locally_owned())
286 const unsigned int total_cells =
289 ExcMessage(
"Not all cells are assigned to a processor."))
294 if (
settings & construct_multigrid_hierarchy)
296 unsigned int n_my_cells = 0;
297 typename parallel::shared::Triangulation<dim, spacedim>::cell_iterator
298 cell = this->
begin(),
300 for (; cell != endc; ++cell)
301 if (cell->is_locally_owned_on_level())
304 const unsigned int total_cells =
307 ExcMessage(
"Not all cells are assigned to a processor."))
314 template <
int dim,
int spacedim>
323 template <
int dim,
int spacedim>
324 const std::vector<types::subdomain_id> &
332 template <
int dim,
int spacedim>
333 const std::vector<types::subdomain_id> &
335 const unsigned int level)
const 342 template <
int dim,
int spacedim>
353 template <
int dim,
int spacedim>
366 const typename ::Triangulation<dim, spacedim>::DistortedCellList
379 template <
int dim,
int spacedim>
382 const ::Triangulation<dim, spacedim> &other_tria)
386 const ::parallel::distributed::Triangulation<dim, spacedim> *
>(
387 &other_tria) ==
nullptr),
389 "Cannot use this function on parallel::distributed::Triangulation."));
399 template <
int dim,
int spacedim>
417 template <
int dim,
int spacedim>
427 template <
int dim,
int spacedim>
428 const std::vector<unsigned int> &
437 template <
int dim,
int spacedim>
438 const std::vector<unsigned int> &
440 const unsigned int)
const 453 #include "shared_tria.inst" 455 DEAL_II_NAMESPACE_CLOSE
std::vector< std::vector< types::subdomain_id > > true_level_subdomain_ids_of_cells
cell_iterator end() const
types::subdomain_id my_subdomain
std::vector< Point< spacedim > > vertices
const std::vector< types::subdomain_id > & get_true_subdomain_ids_of_cells() const
void fill_level_ghost_owners()
virtual void execute_coarsening_and_refinement() override
#define AssertThrow(cond, exc)
cell_iterator begin(const unsigned int level=0) const
unsigned int n_active_cells() const
unsigned int n_levels() const
virtual void execute_coarsening_and_refinement()
virtual MPI_Comm get_communicator() const
virtual void copy_triangulation(const ::Triangulation< dim, spacedim > &old_tria) override
static::ExceptionBase & ExcMessage(std::string arg1)
virtual void update_number_cache() override
unsigned int n_cells() const
T sum(const T &t, const MPI_Comm &mpi_communicator)
virtual void create_triangulation(const std::vector< Point< spacedim >> &vertices, const std::vector< CellData< dim >> &cells, const SubCellData &subcelldata)
#define Assert(cond, exc)
virtual void update_number_cache()
const bool allow_artificial_cells
bool with_artificial_cells() const
Triangulation(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)
const types::subdomain_id artificial_subdomain_id
std::vector< types::subdomain_id > true_subdomain_ids_of_cells
virtual void create_triangulation(const std::vector< Point< spacedim >> &vertices, const std::vector< CellData< dim >> &cells, const SubCellData &subcelldata) override
static::ExceptionBase & ExcNotImplemented()
active_cell_iterator begin_active(const unsigned int level=0) const
virtual void copy_triangulation(const ::Triangulation< dim, spacedim > &other_tria) override
T max(const T &t, const MPI_Comm &mpi_communicator)
const std::vector< types::subdomain_id > & get_true_level_subdomain_ids_of_cells(const unsigned int level) const
types::subdomain_id n_subdomains
static::ExceptionBase & ExcInternalError()