Reference documentation for deal.II version 9.3.3
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
shared_tria.cc
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2015 - 2021 by the deal.II authors
4//
5// This file is part of the deal.II library.
6//
7// The deal.II library is free software; you can use it, redistribute
8// it, and/or modify it under the terms of the GNU Lesser General
9// Public License as published by the Free Software Foundation; either
10// version 2.1 of the License, or (at your option) any later version.
11// The full text of the license can be found in the file LICENSE.md at
12// the top level directory of deal.II.
13//
14// ---------------------------------------------------------------------
15
16#include <deal.II/base/mpi.h>
18
21
24#include <deal.II/grid/tria.h>
27
29
30
32
33#ifdef DEAL_II_WITH_MPI
34namespace parallel
35{
36 namespace shared
37 {
38 template <int dim, int spacedim>
40 const MPI_Comm &mpi_communicator,
41 const typename ::Triangulation<dim, spacedim>::MeshSmoothing
42 smooth_grid,
43 const bool allow_artificial_cells,
44 const Settings settings)
45 : ::parallel::TriangulationBase<dim, spacedim>(mpi_communicator,
46 smooth_grid,
47 false)
48 , settings(settings)
49 , allow_artificial_cells(allow_artificial_cells)
50 {
51 const auto partition_settings =
55 (void)partition_settings;
56 Assert(partition_settings == partition_auto ||
57 partition_settings == partition_metis ||
58 partition_settings == partition_zoltan ||
59 partition_settings == partition_zorder ||
60 partition_settings == partition_custom_signal,
61 ExcMessage("Settings must contain exactly one type of the active "
62 "cell partitioning scheme."));
63
66 ExcMessage("construct_multigrid_hierarchy requires "
67 "allow_artificial_cells to be set to true."));
68 }
69
70
71
72 template <int dim, int spacedim>
73 bool
75 {
76 return (settings & construct_multigrid_hierarchy);
77 }
78
79
80
81 template <int dim, int spacedim>
82 void
84 {
85# ifdef DEBUG
86 // Check that all meshes are the same (or at least have the same
87 // total number of active cells):
88 const unsigned int max_active_cells =
89 Utilities::MPI::max(this->n_active_cells(), this->get_communicator());
90 Assert(
91 max_active_cells == this->n_active_cells(),
93 "A parallel::shared::Triangulation needs to be refined in the same "
94 "way on all processors, but the participating processors don't "
95 "agree on the number of active cells."));
96# endif
97
98 auto partition_settings = (partition_zoltan | partition_metis |
99 partition_zorder | partition_custom_signal) &
100 settings;
101 if (partition_settings == partition_auto)
102# ifdef DEAL_II_TRILINOS_WITH_ZOLTAN
103 partition_settings = partition_zoltan;
104# elif defined DEAL_II_WITH_METIS
105 partition_settings = partition_metis;
106# else
107 partition_settings = partition_zorder;
108# endif
109
110 if (partition_settings == partition_zoltan)
111 {
112# ifndef DEAL_II_TRILINOS_WITH_ZOLTAN
113 AssertThrow(false,
115 "Choosing 'partition_zoltan' requires the library "
116 "to be compiled with support for Zoltan! "
117 "Instead, you might use 'partition_auto' to select "
118 "a partitioning algorithm that is supported "
119 "by your current configuration."));
120# else
122 this->n_subdomains, *this, SparsityTools::Partitioner::zoltan);
123# endif
124 }
125 else if (partition_settings == partition_metis)
126 {
127# ifndef DEAL_II_WITH_METIS
128 AssertThrow(false,
130 "Choosing 'partition_metis' requires the library "
131 "to be compiled with support for METIS! "
132 "Instead, you might use 'partition_auto' to select "
133 "a partitioning algorithm that is supported "
134 "by your current configuration."));
135# else
136 GridTools::partition_triangulation(this->n_subdomains,
137 *this,
139# endif
140 }
141 else if (partition_settings == partition_zorder)
142 {
143 GridTools::partition_triangulation_zorder(this->n_subdomains, *this);
144 }
145 else if (partition_settings == partition_custom_signal)
146 {
147 // User partitions mesh manually
148 }
149 else
150 {
152 }
153
154 // do not partition multigrid levels if user is
155 // defining a custom partition
156 if ((settings & construct_multigrid_hierarchy) &&
157 !(settings & partition_custom_signal))
159
160 true_subdomain_ids_of_cells.resize(this->n_active_cells());
161
162 // loop over all cells and mark artificial:
164 spacedim>::active_cell_iterator
165 cell = this->begin_active(),
166 endc = this->end();
167
168 if (allow_artificial_cells)
169 {
170 // get active halo layer of (ghost) cells
171 // parallel::shared::Triangulation<dim>::
172 std::function<bool(
175 predicate = IteratorFilters::SubdomainEqualTo(this->my_subdomain);
176
177 const std::vector<typename parallel::shared::Triangulation<
178 dim,
179 spacedim>::active_cell_iterator>
180 active_halo_layer_vector =
182 predicate);
183 std::set<typename parallel::shared::Triangulation<dim, spacedim>::
184 active_cell_iterator>
185 active_halo_layer(active_halo_layer_vector.begin(),
186 active_halo_layer_vector.end());
187
188 for (unsigned int index = 0; cell != endc; cell++, index++)
189 {
190 // store original/true subdomain ids:
191 true_subdomain_ids_of_cells[index] = cell->subdomain_id();
192
193 if (cell->is_locally_owned() == false &&
194 active_halo_layer.find(cell) == active_halo_layer.end())
195 cell->set_subdomain_id(numbers::artificial_subdomain_id);
196 }
197
198 // loop over all cells in multigrid hierarchy and mark artificial:
199 if (settings & construct_multigrid_hierarchy)
200 {
201 true_level_subdomain_ids_of_cells.resize(this->n_levels());
202
203 std::function<bool(
205 cell_iterator &)>
207 for (unsigned int lvl = 0; lvl < this->n_levels(); ++lvl)
208 {
209 true_level_subdomain_ids_of_cells[lvl].resize(
210 this->n_cells(lvl));
211
212 const std::vector<typename parallel::shared::Triangulation<
213 dim,
214 spacedim>::cell_iterator>
215 level_halo_layer_vector =
217 *this, predicate, lvl);
220 level_halo_layer(level_halo_layer_vector.begin(),
221 level_halo_layer_vector.end());
222
224 cell_iterator cell = this->begin(lvl),
225 endc = this->end(lvl);
226 for (unsigned int index = 0; cell != endc; cell++, index++)
227 {
228 // Store true level subdomain IDs before setting
229 // artificial
230 true_level_subdomain_ids_of_cells[lvl][index] =
231 cell->level_subdomain_id();
232
233 // for active cells, we must have knowledge of level
234 // subdomain ids of all neighbors to our subdomain, not
235 // just neighbors on the same level. if the cells
236 // subdomain id was not set to artitficial above, we will
237 // also keep its level subdomain id since it is either
238 // owned by this processor or in the ghost layer of the
239 // active mesh.
240 if (cell->is_active() &&
241 cell->subdomain_id() !=
243 continue;
244
245 // we must have knowledge of our parent in the hierarchy
246 if (cell->has_children())
247 {
248 bool keep_cell = false;
249 for (unsigned int c = 0;
250 c < GeometryInfo<dim>::max_children_per_cell;
251 ++c)
252 if (cell->child(c)->level_subdomain_id() ==
253 this->my_subdomain)
254 {
255 keep_cell = true;
256 break;
257 }
258 if (keep_cell)
259 continue;
260 }
261
262 // we must have knowledge of our neighbors on the same
263 // level
264 if (!cell->is_locally_owned_on_level() &&
265 level_halo_layer.find(cell) != level_halo_layer.end())
266 continue;
267
268 // mark all other cells to artificial
269 cell->set_level_subdomain_id(
271 }
272 }
273 }
274 }
275 else
276 {
277 // just store true subdomain ids
278 for (unsigned int index = 0; cell != endc; cell++, index++)
279 true_subdomain_ids_of_cells[index] = cell->subdomain_id();
280 }
281
282# ifdef DEBUG
283 {
284 // Assert that each cell is owned by a processor
285 unsigned int n_my_cells = 0;
287 spacedim>::active_cell_iterator
288 cell = this->begin_active(),
289 endc = this->end();
290 for (; cell != endc; ++cell)
291 if (cell->is_locally_owned())
292 n_my_cells += 1;
293
294 const unsigned int total_cells =
295 Utilities::MPI::sum(n_my_cells, this->get_communicator());
296 Assert(total_cells == this->n_active_cells(),
297 ExcMessage("Not all cells are assigned to a processor."));
298 }
299
300 // If running with multigrid, assert that each level
301 // cell is owned by a processor
302 if (settings & construct_multigrid_hierarchy)
303 {
304 unsigned int n_my_cells = 0;
306 cell = this->begin(),
307 endc = this->end();
308 for (; cell != endc; ++cell)
309 if (cell->is_locally_owned_on_level())
310 n_my_cells += 1;
311
312 const unsigned int total_cells =
313 Utilities::MPI::sum(n_my_cells, this->get_communicator());
314 Assert(total_cells == this->n_cells(),
315 ExcMessage("Not all cells are assigned to a processor."));
316 }
317# endif
318 }
319
320
321
322 template <int dim, int spacedim>
323 bool
325 {
326 return allow_artificial_cells;
327 }
328
329
330
331 template <int dim, int spacedim>
332 const std::vector<types::subdomain_id> &
334 {
335 return true_subdomain_ids_of_cells;
336 }
337
338
339
340 template <int dim, int spacedim>
341 const std::vector<types::subdomain_id> &
343 const unsigned int level) const
344 {
345 Assert(level < true_level_subdomain_ids_of_cells.size(),
347 Assert(true_level_subdomain_ids_of_cells[level].size() ==
348 this->n_cells(level),
350 return true_level_subdomain_ids_of_cells[level];
351 }
352
353
354
355 template <int dim, int spacedim>
356 void
358 {
360 partition();
361 this->update_number_cache();
362 }
363
364
365
366 template <int dim, int spacedim>
367 void
369 const std::vector<Point<spacedim>> &vertices,
370 const std::vector<CellData<dim>> & cells,
371 const SubCellData & subcelldata)
372 {
373 try
374 {
376 vertices, cells, subcelldata);
377 }
378 catch (
379 const typename ::Triangulation<dim, spacedim>::DistortedCellList
380 &)
381 {
382 // the underlying triangulation should not be checking for distorted
383 // cells
384 Assert(false, ExcInternalError());
385 }
386 partition();
387 this->update_number_cache();
388 }
389
390
391
392 template <int dim, int spacedim>
393 void
396 &construction_data)
397 {
398 (void)construction_data;
399
400 Assert(false, ExcInternalError());
401 }
402
403
404
405 template <int dim, int spacedim>
406 void
408 const ::Triangulation<dim, spacedim> &other_tria)
409 {
410 Assert(
411 (dynamic_cast<
412 const ::parallel::distributed::Triangulation<dim, spacedim> *>(
413 &other_tria) == nullptr),
415 "Cannot use this function on parallel::distributed::Triangulation."));
416
418 other_tria);
419 partition();
420 this->update_number_cache();
421 }
422 } // namespace shared
423} // namespace parallel
424
425#else
426
427namespace parallel
428{
429 namespace shared
430 {
431 template <int dim, int spacedim>
432 bool
434 {
435 Assert(false, ExcNotImplemented());
436 return true;
437 }
438
439 template <int dim, int spacedim>
440 bool
442 {
443 return false;
444 }
445
446 template <int dim, int spacedim>
447 const std::vector<unsigned int> &
449 {
450 Assert(false, ExcNotImplemented());
451 return true_subdomain_ids_of_cells;
452 }
453
454
455
456 template <int dim, int spacedim>
457 const std::vector<unsigned int> &
459 const unsigned int) const
460 {
461 Assert(false, ExcNotImplemented());
462 return true_level_subdomain_ids_of_cells;
463 }
464 } // namespace shared
465} // namespace parallel
466
467
468#endif
469
470
471
472namespace internal
473{
474 namespace parallel
475 {
476 namespace shared
477 {
478 template <int dim, int spacedim>
481 : shared_tria(
482 dynamic_cast<
483 const ::parallel::shared::Triangulation<dim, spacedim> *>(
484 &tria))
485 {
486 if (shared_tria && shared_tria->with_artificial_cells())
487 {
488 // Save the current set of subdomain IDs, and set subdomain IDs
489 // to the "true" owner of each cell.
490 const std::vector<types::subdomain_id> &true_subdomain_ids =
491 shared_tria->get_true_subdomain_ids_of_cells();
492
493 saved_subdomain_ids.resize(shared_tria->n_active_cells());
494 for (const auto &cell : shared_tria->active_cell_iterators())
495 {
496 const unsigned int index = cell->active_cell_index();
497 saved_subdomain_ids[index] = cell->subdomain_id();
498 cell->set_subdomain_id(true_subdomain_ids[index]);
499 }
500 }
501 }
502
503
504
505 template <int dim, int spacedim>
508 {
509 if (shared_tria && shared_tria->with_artificial_cells())
510 {
511 // Undo the subdomain modification.
512 for (const auto &cell : shared_tria->active_cell_iterators())
513 {
514 const unsigned int index = cell->active_cell_index();
515 cell->set_subdomain_id(saved_subdomain_ids[index]);
516 }
517 }
518 }
519 } // namespace shared
520 } // namespace parallel
521} // namespace internal
522
523
524/*-------------- Explicit Instantiations -------------------------------*/
525#include "shared_tria.inst"
526
virtual void create_triangulation(const std::vector< Point< spacedim > > &vertices, const std::vector< CellData< dim > > &cells, const SubCellData &subcelldata)
cell_iterator end() const
virtual void execute_coarsening_and_refinement()
const SmartPointer< const ::parallel::shared::Triangulation< dim, spacedim > > shared_tria
Definition: shared_tria.h:535
TemporarilyRestoreSubdomainIds(const Triangulation< dim, spacedim > &tria)
Definition: shared_tria.cc:480
types::subdomain_id my_subdomain
Definition: tria_base.h:316
virtual void copy_triangulation(const ::Triangulation< dim, spacedim > &old_tria) override
Definition: tria_base.cc:65
virtual void execute_coarsening_and_refinement() override
Definition: shared_tria.cc:357
typename ::Triangulation< dim, spacedim >::active_cell_iterator active_cell_iterator
Definition: shared_tria.h:110
virtual void copy_triangulation(const ::Triangulation< dim, spacedim > &other_tria) override
Definition: shared_tria.cc:407
const std::vector< types::subdomain_id > & get_true_subdomain_ids_of_cells() const
Definition: shared_tria.cc:333
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)
Definition: shared_tria.cc:39
typename ::Triangulation< dim, spacedim >::cell_iterator cell_iterator
Definition: shared_tria.h:112
const std::vector< types::subdomain_id > & get_true_level_subdomain_ids_of_cells(const unsigned int level) const
Definition: shared_tria.cc:342
virtual bool is_multilevel_hierarchy_constructed() const override
Definition: shared_tria.cc:74
virtual void create_triangulation(const std::vector< Point< spacedim > > &vertices, const std::vector< CellData< dim > > &cells, const SubCellData &subcelldata) override
Definition: shared_tria.cc:368
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:402
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:403
Point< 3 > vertices[4]
unsigned int level
Definition: grid_out.cc:4590
__global__ void set(Number *val, const Number s, const size_type N)
static ::ExceptionBase & ExcNotImplemented()
#define Assert(cond, exc)
Definition: exceptions.h:1465
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)
Definition: exceptions.h:1575
void partition_triangulation_zorder(const unsigned int n_partitions, Triangulation< dim, spacedim > &triangulation, const bool group_siblings=true)
Definition: grid_tools.cc:3989
std::vector< typename MeshType::active_cell_iterator > compute_active_cell_halo_layer(const MeshType &mesh, const std::function< bool(const typename MeshType::active_cell_iterator &)> &predicate)
void partition_multigrid_levels(Triangulation< dim, spacedim > &triangulation)
Definition: grid_tools.cc:4094
std::vector< typename MeshType::cell_iterator > compute_cell_halo_layer_on_level(const MeshType &mesh, const std::function< bool(const typename MeshType::cell_iterator &)> &predicate, const unsigned int level)
void partition_triangulation(const unsigned int n_partitions, Triangulation< dim, spacedim > &triangulation, const SparsityTools::Partitioner partitioner=SparsityTools::Partitioner::metis)
Definition: grid_tools.cc:3762
void partition(const SparsityPattern &sparsity_pattern, const unsigned int n_partitions, std::vector< unsigned int > &partition_indices, const Partitioner partitioner=Partitioner::metis)
VectorType::value_type * end(VectorType &V)
VectorType::value_type * begin(VectorType &V)
T sum(const T &t, const MPI_Comm &mpi_communicator)
T max(const T &t, const MPI_Comm &mpi_communicator)
unsigned int n_active_cells(const internal::TriangulationImplementation::NumberCache< 1 > &c)
Definition: tria.cc:12594
unsigned int n_cells(const internal::TriangulationImplementation::NumberCache< 1 > &c)
Definition: tria.cc:12587
const types::subdomain_id artificial_subdomain_id
Definition: types.h:293