Reference documentation for deal.II version 9.2.0
|
#include <deal.II/distributed/fully_distributed_tria.h>
Public Types | |
using | cell_iterator = typename ::Triangulation< dim, spacedim >::cell_iterator |
using | active_cell_iterator = typename ::Triangulation< dim, spacedim >::active_cell_iterator |
using | CellStatus = typename ::Triangulation< dim, spacedim >::CellStatus |
Public Types inherited from Triangulation< dim, dim > | |
enum | MeshSmoothing |
using | cell_iterator = TriaIterator< CellAccessor< dim, spacedim > > |
using | active_cell_iterator = TriaActiveIterator< CellAccessor< dim, spacedim > > |
using | face_iterator = TriaIterator< TriaAccessor< dim - 1, dim, spacedim > > |
using | active_face_iterator = TriaActiveIterator< TriaAccessor< dim - 1, dim, spacedim > > |
using | vertex_iterator = TriaIterator<::TriaAccessor< 0, dim, spacedim > > |
using | active_vertex_iterator = TriaActiveIterator<::TriaAccessor< 0, dim, spacedim > > |
using | line_iterator = typename IteratorSelector::line_iterator |
using | active_line_iterator = typename IteratorSelector::active_line_iterator |
using | quad_iterator = typename IteratorSelector::quad_iterator |
using | active_quad_iterator = typename IteratorSelector::active_quad_iterator |
using | hex_iterator = typename IteratorSelector::hex_iterator |
using | active_hex_iterator = typename IteratorSelector::active_hex_iterator |
enum | CellStatus |
Public Member Functions | |
Triangulation (MPI_Comm mpi_communicator) | |
virtual | ~Triangulation ()=default |
void | create_triangulation (const TriangulationDescription::Description< dim, spacedim > &construction_data) override |
virtual void | create_triangulation (const std::vector< Point< spacedim >> &vertices, const std::vector<::CellData< dim >> &cells, const SubCellData &subcelldata) override |
void | copy_triangulation (const ::Triangulation< dim, spacedim > &other_tria) override |
void | set_partitioner (const std::function< void(::Triangulation< dim, spacedim > &, const unsigned int)> &partitioner, const TriangulationDescription::Settings &settings) |
virtual void | execute_coarsening_and_refinement () override |
virtual bool | prepare_coarsening_and_refinement () override |
virtual bool | has_hanging_nodes () const override |
virtual std::size_t | memory_consumption () const override |
virtual bool | is_multilevel_hierarchy_constructed () const override |
Public Member Functions inherited from parallel::DistributedTriangulationBase< dim, dim > | |
DistributedTriangulationBase (MPI_Comm mpi_communicator, const typename ::Triangulation< dim, spacedim >::MeshSmoothing smooth_grid=(::Triangulation< dim, spacedim >::none), const bool check_for_distorted_cells=false) | |
Public Member Functions inherited from parallel::TriangulationBase< dim, spacedim > | |
TriangulationBase (MPI_Comm mpi_communicator, const typename ::Triangulation< dim, spacedim >::MeshSmoothing smooth_grid=(::Triangulation< dim, spacedim >::none), const bool check_for_distorted_cells=false) | |
virtual | ~TriangulationBase () override |
virtual const MPI_Comm & | get_communicator () const |
virtual void | copy_triangulation (const ::Triangulation< dim, spacedim > &old_tria) override |
unsigned int | n_locally_owned_active_cells () const |
virtual types::global_cell_index | n_global_active_cells () const override |
virtual unsigned int | n_global_levels () const override |
types::subdomain_id | locally_owned_subdomain () const override |
const std::set< types::subdomain_id > & | ghost_owners () const |
const std::set< types::subdomain_id > & | level_ghost_owners () const |
virtual std::map< unsigned int, std::set<::types::subdomain_id > > | compute_vertices_with_ghost_neighbors () const |
virtual std::vector< types::boundary_id > | get_boundary_ids () const override |
virtual std::vector< types::manifold_id > | get_manifold_ids () const override |
Public Member Functions inherited from Triangulation< dim, dim > | |
Triangulation (const MeshSmoothing smooth_grid=none, const bool check_for_distorted_cells=false) | |
Triangulation (const Triangulation< dim, spacedim > &)=delete | |
Triangulation (Triangulation< dim, spacedim > &&tria) noexcept | |
Triangulation & | operator= (Triangulation< dim, spacedim > &&tria) noexcept |
virtual void | clear () |
virtual void | set_mesh_smoothing (const MeshSmoothing mesh_smoothing) |
virtual const MeshSmoothing & | get_mesh_smoothing () const |
void | set_manifold (const types::manifold_id number, const Manifold< dim, spacedim > &manifold_object) |
void | reset_manifold (const types::manifold_id manifold_number) |
void | reset_all_manifolds () |
void | set_all_manifold_ids (const types::manifold_id number) |
void | set_all_manifold_ids_on_boundary (const types::manifold_id number) |
void | set_all_manifold_ids_on_boundary (const types::boundary_id b_id, const types::manifold_id number) |
const Manifold< dim, spacedim > & | get_manifold (const types::manifold_id number) const |
virtual std::vector< types::boundary_id > | get_boundary_ids () const |
virtual std::vector< types::manifold_id > | get_manifold_ids () const |
virtual void | copy_triangulation (const Triangulation< dim, spacedim > &other_tria) |
virtual void | create_triangulation (const std::vector< Point< spacedim >> &vertices, const std::vector< CellData< dim >> &cells, const SubCellData &subcelldata) |
virtual void | create_triangulation_compatibility (const std::vector< Point< spacedim >> &vertices, const std::vector< CellData< dim >> &cells, const SubCellData &subcelldata) |
void | flip_all_direction_flags () |
bool | prepare_coarsening_and_refinement () |
bool | prepare_coarsening_and_refinement () |
bool | prepare_coarsening_and_refinement () |
unsigned int | n_quads () const |
unsigned int | n_quads (const unsigned int) const |
unsigned int | n_quads () const |
unsigned int | n_quads (const unsigned int) const |
unsigned int | n_quads () const |
unsigned int | n_quads (const unsigned int) const |
unsigned int | n_active_quads (const unsigned int) const |
unsigned int | n_active_quads () const |
unsigned int | n_active_quads (const unsigned int) const |
unsigned int | n_active_quads () const |
unsigned int | n_active_quads (const unsigned int) const |
unsigned int | n_active_quads () const |
unsigned int | n_hexs () const |
unsigned int | n_hexs (const unsigned int level) const |
unsigned int | n_active_hexs () const |
unsigned int | n_active_hexs (const unsigned int level) const |
unsigned int | max_adjacent_cells () const |
unsigned int | max_adjacent_cells () const |
unsigned int | max_adjacent_cells () const |
unsigned int | n_raw_lines (const unsigned int level) const |
unsigned int | n_raw_lines () const |
unsigned int | n_raw_lines (const unsigned int level) const |
unsigned int | n_raw_lines () const |
unsigned int | n_raw_lines (const unsigned int level) const |
unsigned int | n_raw_lines () const |
unsigned int | n_raw_quads (const unsigned int) const |
unsigned int | n_raw_quads (const unsigned int) const |
unsigned int | n_raw_quads (const unsigned int) const |
unsigned int | n_raw_quads (const unsigned int level) const |
unsigned int | n_raw_quads (const unsigned int level) const |
unsigned int | n_raw_quads (const unsigned int) const |
unsigned int | n_raw_quads () const |
unsigned int | n_raw_hexs (const unsigned int) const |
unsigned int | n_raw_hexs (const unsigned int) const |
unsigned int | n_raw_hexs (const unsigned int) const |
unsigned int | n_raw_hexs (const unsigned int level) const |
void | set_all_refine_flags () |
void | refine_global (const unsigned int times=1) |
void | save_refine_flags (std::ostream &out) const |
void | save_refine_flags (std::vector< bool > &v) const |
void | load_refine_flags (std::istream &in) |
void | load_refine_flags (const std::vector< bool > &v) |
void | save_coarsen_flags (std::ostream &out) const |
void | save_coarsen_flags (std::vector< bool > &v) const |
void | load_coarsen_flags (std::istream &out) |
void | load_coarsen_flags (const std::vector< bool > &v) |
bool | get_anisotropic_refinement_flag () const |
void | clear_user_flags () |
void | save_user_flags (std::ostream &out) const |
void | save_user_flags (std::vector< bool > &v) const |
void | load_user_flags (std::istream &in) |
void | load_user_flags (const std::vector< bool > &v) |
void | clear_user_flags_line () |
void | save_user_flags_line (std::ostream &out) const |
void | save_user_flags_line (std::vector< bool > &v) const |
void | load_user_flags_line (std::istream &in) |
void | load_user_flags_line (const std::vector< bool > &v) |
void | clear_user_flags_quad () |
void | save_user_flags_quad (std::ostream &out) const |
void | save_user_flags_quad (std::vector< bool > &v) const |
void | load_user_flags_quad (std::istream &in) |
void | load_user_flags_quad (const std::vector< bool > &v) |
void | clear_user_flags_hex () |
void | save_user_flags_hex (std::ostream &out) const |
void | save_user_flags_hex (std::vector< bool > &v) const |
void | load_user_flags_hex (std::istream &in) |
void | load_user_flags_hex (const std::vector< bool > &v) |
void | clear_user_data () |
void | save_user_indices (std::vector< unsigned int > &v) const |
void | load_user_indices (const std::vector< unsigned int > &v) |
void | save_user_pointers (std::vector< void * > &v) const |
void | load_user_pointers (const std::vector< void * > &v) |
void | save_user_indices_line (std::vector< unsigned int > &v) const |
void | load_user_indices_line (const std::vector< unsigned int > &v) |
void | save_user_indices_quad (std::vector< unsigned int > &v) const |
void | load_user_indices_quad (const std::vector< unsigned int > &v) |
void | save_user_indices_hex (std::vector< unsigned int > &v) const |
void | load_user_indices_hex (const std::vector< unsigned int > &v) |
void | save_user_pointers_line (std::vector< void * > &v) const |
void | load_user_pointers_line (const std::vector< void * > &v) |
void | save_user_pointers_quad (std::vector< void * > &v) const |
void | load_user_pointers_quad (const std::vector< void * > &v) |
void | save_user_pointers_hex (std::vector< void * > &v) const |
void | load_user_pointers_hex (const std::vector< void * > &v) |
cell_iterator | begin (const unsigned int level=0) const |
active_cell_iterator | begin_active (const unsigned int level=0) const |
cell_iterator | end () const |
cell_iterator | end (const unsigned int level) const |
active_cell_iterator | end_active (const unsigned int level) const |
cell_iterator | last () const |
active_cell_iterator | last_active () const |
IteratorRange< cell_iterator > | cell_iterators () const |
IteratorRange< active_cell_iterator > | active_cell_iterators () const |
IteratorRange< cell_iterator > | cell_iterators_on_level (const unsigned int level) const |
IteratorRange< active_cell_iterator > | active_cell_iterators_on_level (const unsigned int level) const |
face_iterator | begin_face () const |
active_face_iterator | begin_active_face () const |
face_iterator | end_face () const |
IteratorRange< active_face_iterator > | active_face_iterators () const |
vertex_iterator | begin_vertex () const |
active_vertex_iterator | begin_active_vertex () const |
vertex_iterator | end_vertex () const |
unsigned int | n_lines () const |
unsigned int | n_lines (const unsigned int level) const |
unsigned int | n_active_lines () const |
unsigned int | n_active_lines (const unsigned int level) const |
unsigned int | n_quads () const |
unsigned int | n_quads (const unsigned int level) const |
unsigned int | n_active_quads () const |
unsigned int | n_active_quads (const unsigned int level) const |
unsigned int | n_hexs () const |
unsigned int | n_hexs (const unsigned int level) const |
unsigned int | n_active_hexs () const |
unsigned int | n_active_hexs (const unsigned int level) const |
unsigned int | n_cells () const |
unsigned int | n_cells (const unsigned int level) const |
unsigned int | n_active_cells () const |
unsigned int | n_active_cells (const unsigned int level) const |
virtual types::global_cell_index | n_global_active_cells () const |
unsigned int | n_faces () const |
unsigned int | n_active_faces () const |
unsigned int | n_levels () const |
virtual unsigned int | n_global_levels () const |
virtual bool | has_hanging_nodes () const |
unsigned int | n_vertices () const |
const std::vector< Point< spacedim > > & | get_vertices () const |
unsigned int | n_used_vertices () const |
bool | vertex_used (const unsigned int index) const |
const std::vector< bool > & | get_used_vertices () const |
unsigned int | max_adjacent_cells () const |
virtual types::subdomain_id | locally_owned_subdomain () const |
Triangulation< dim, spacedim > & | get_triangulation () |
const Triangulation< dim, spacedim > & | get_triangulation () const |
unsigned int | n_raw_lines () const |
unsigned int | n_raw_lines (const unsigned int level) const |
unsigned int | n_raw_quads () const |
unsigned int | n_raw_quads (const unsigned int level) const |
unsigned int | n_raw_hexs (const unsigned int level) const |
unsigned int | n_raw_cells (const unsigned int level) const |
unsigned int | n_raw_faces () const |
virtual std::size_t | memory_consumption () const |
void | save (Archive &ar, const unsigned int version) const |
void | load (Archive &ar, const unsigned int version) |
virtual void | add_periodicity (const std::vector< GridTools::PeriodicFacePair< cell_iterator >> &) |
const std::map< std::pair< cell_iterator, unsigned int >, std::pair< std::pair< cell_iterator, unsigned int >, std::bitset< 3 > > > & | get_periodic_face_map () const |
void | serialize (Archive &archive, const unsigned int version) |
Public Member Functions inherited from Subscriptor | |
Subscriptor () | |
Subscriptor (const Subscriptor &) | |
Subscriptor (Subscriptor &&) noexcept | |
virtual | ~Subscriptor () |
Subscriptor & | operator= (const Subscriptor &) |
Subscriptor & | operator= (Subscriptor &&) noexcept |
void | subscribe (std::atomic< bool > *const validity, const std::string &identifier="") const |
void | unsubscribe (std::atomic< bool > *const validity, const std::string &identifier="") const |
unsigned int | n_subscriptions () const |
template<typename StreamType > | |
void | list_subscribers (StreamType &stream) const |
void | list_subscribers () const |
template<class Archive > | |
void | serialize (Archive &ar, const unsigned int version) |
Private Member Functions | |
virtual unsigned int | coarse_cell_id_to_coarse_cell_index (const types::coarse_cell_id coarse_cell_id) const override |
virtual types::coarse_cell_id | coarse_cell_index_to_coarse_cell_id (const unsigned int coarse_cell_index) const override |
Private Attributes | |
TriangulationDescription::Settings | settings |
std::function< void(::Triangulation< dim, spacedim > &, const unsigned int)> | partitioner |
std::vector< std::pair< types::coarse_cell_id, unsigned int > > | coarse_cell_id_to_coarse_cell_index_vector |
std::vector< types::coarse_cell_id > | coarse_cell_index_to_coarse_cell_id_vector |
bool | currently_processing_create_triangulation_for_internal_usage |
bool | currently_processing_prepare_coarsening_and_refinement_for_internal_usage |
Additional Inherited Members | |
Static Public Member Functions inherited from Triangulation< dim, dim > | |
static ::ExceptionBase & | ExcInvalidLevel (int arg1, int arg2) |
static ::ExceptionBase & | ExcTriangulationNotEmpty (int arg1, int arg2) |
static ::ExceptionBase & | ExcGridReadError () |
static ::ExceptionBase & | ExcFacesHaveNoLevel () |
static ::ExceptionBase & | ExcEmptyLevel (int arg1) |
static ::ExceptionBase & | ExcNonOrientableTriangulation () |
static ::ExceptionBase & | ExcBoundaryIdNotFound (types::boundary_id arg1) |
static ::ExceptionBase & | ExcInconsistentCoarseningFlags () |
Static Public Member Functions inherited from Subscriptor | |
static ::ExceptionBase & | ExcInUse (int arg1, std::string arg2, std::string arg3) |
static ::ExceptionBase & | ExcNoSubscriber (std::string arg1, std::string arg2) |
Public Attributes inherited from Triangulation< dim, dim > | |
Signals | signals |
Static Public Attributes inherited from Triangulation< dim, dim > | |
static const unsigned int | dimension |
static const unsigned int | space_dimension |
Protected Member Functions inherited from parallel::TriangulationBase< dim, spacedim > | |
virtual void | update_number_cache () |
Protected Member Functions inherited from Triangulation< dim, dim > | |
void | update_periodic_face_map () |
Static Protected Member Functions inherited from Triangulation< dim, dim > | |
static void | write_bool_vector (const unsigned int magic_number1, const std::vector< bool > &v, const unsigned int magic_number2, std::ostream &out) |
static void | read_bool_vector (const unsigned int magic_number1, std::vector< bool > &v, const unsigned int magic_number2, std::istream &in) |
Protected Attributes inherited from parallel::TriangulationBase< dim, spacedim > | |
MPI_Comm | mpi_communicator |
types::subdomain_id | my_subdomain |
types::subdomain_id | n_subdomains |
NumberCache | number_cache |
Protected Attributes inherited from Triangulation< dim, dim > | |
MeshSmoothing | smooth_grid |
A distributed triangulation with a distributed coarse grid.
The motivation for parallel::fullydistributed::Triangulation has its origins in the following observations about complex geometries and/or about given meshes created by an external mesh generator. We regard complex geometries as geometries that can be meshed only with a non-negligible number of coarse cells (>10,000):
To be able to construct a fully partitioned triangulation that distributes the coarse grid and gives flexibility regarding partitioning, the following ingredients are required:
The ingredients listed above are bundled in the struct TriangulationDescription::Description. The user has to fill this data structure - in a pre-processing step - before actually creating the triangulation. Predefined functions to create TriangulationDescription::Description can be found in the namespace TriangulationDescription::Utilities.
Once the TriangulationDescription::Description construction_data
has been constructed, the triangulation tria
can be created by calling tria.create_triangulation(construction_data);
.
Definition at line 115 of file fully_distributed_tria.h.
using parallel::fullydistributed::Triangulation< dim, spacedim >::cell_iterator = typename ::Triangulation<dim, spacedim>::cell_iterator |
Definition at line 120 of file fully_distributed_tria.h.
using parallel::fullydistributed::Triangulation< dim, spacedim >::active_cell_iterator = typename ::Triangulation<dim, spacedim>::active_cell_iterator |
Definition at line 123 of file fully_distributed_tria.h.
using parallel::fullydistributed::Triangulation< dim, spacedim >::CellStatus = typename ::Triangulation<dim, spacedim>::CellStatus |
Definition at line 126 of file fully_distributed_tria.h.
|
explicit |
Constructor.
mpi_communicator | The MPI communicator to be used for the triangulation. |
Definition at line 43 of file fully_distributed_tria.cc.
|
virtualdefault |
Destructor.
Reimplemented from Triangulation< dim, dim >.
|
overridevirtual |
Create a triangulation from a list of vertices and a list of cells, each of the latter being a list of 1<<dim
vertex indices. The triangulation must be empty upon calling this function and the cell list should be useful (connected domain, etc.). The result of calling this function is a coarse mesh.
Material data for the cells is given within the cells
array, while boundary information is given in the subcelldata
field.
The numbering of vertices within the cells
array is subject to some constraints; see the general class documentation for this.
For conditions when this function can generate a valid triangulation, see the documentation of this class, and the GridIn and GridReordering class.
If the check_for_distorted_cells
flag was specified upon creation of this object, at the very end of its operation, the current function walks over all cells and verifies that none of the cells is deformed (see the entry on distorted cells in the glossary), where we call a cell deformed if the determinant of the Jacobian of the mapping from reference cell to real cell is negative at least at one of the vertices (this computation is done using the GeometryInfo::jacobian_determinants_at_vertices function). If there are deformed cells, this function throws an exception of kind DistortedCellList. Since this happens after all data structures have been set up, you can catch and ignore this exception if you know what you do – for example, it may be that the determinant is zero (indicating that you have collapsed edges in a cell) but that this is ok because you didn't intend to integrate on this cell anyway. On the other hand, deformed cells are often a sign of a mesh that is too coarse to resolve the geometry of the domain, and in this case ignoring the exception is probably unwise.
Reimplemented from Triangulation< dim, dim >.
Definition at line 59 of file fully_distributed_tria.cc.
|
overridevirtual |
Definition at line 213 of file fully_distributed_tria.cc.
|
override |
Implementation of the same function as in the base class.
other_tria | The triangulation to be copied. It can be a serial Triangulation or a parallel::distributed::Triangulation. Both can have been refined already. |
Definition at line 232 of file fully_distributed_tria.cc.
void Triangulation< dim, spacedim >::set_partitioner | ( | const std::function< void(::Triangulation< dim, spacedim > &, const unsigned int)> & | partitioner, |
const TriangulationDescription::Settings & | settings | ||
) |
Register a partitioner, which is used within the method copy_triangulation.
partitioner | A partitioning function, which takes as input argument a reference to the triangulation to be partitioned and the number of partitions to be created. The function needs to set subdomain ids for each active cell of the given triangulation, with values between zero (inclusive) and the second argument to the function (exclusive). |
settings | See the description of the Settings enumerator. |
Definition at line 286 of file fully_distributed_tria.cc.
|
overridevirtual |
Coarsen and refine the mesh according to refinement and coarsening flags set.
Reimplemented from Triangulation< dim, dim >.
Definition at line 299 of file fully_distributed_tria.cc.
|
overridevirtual |
Override the implementation of prepare_coarsening_and_refinement from the base class.
Reimplemented from Triangulation< dim, dim >.
Definition at line 308 of file fully_distributed_tria.cc.
|
overridevirtual |
Return true if the triangulation has hanging nodes.
Definition at line 322 of file fully_distributed_tria.cc.
|
overridevirtual |
Return the local memory consumption in bytes.
Reimplemented from parallel::TriangulationBase< dim, spacedim >.
Definition at line 332 of file fully_distributed_tria.cc.
|
overridevirtual |
Return if multilevel hierarchy is supported and has been constructed.
Implements parallel::TriangulationBase< dim, spacedim >.
Definition at line 348 of file fully_distributed_tria.cc.
|
overrideprivatevirtual |
Translate the unique id of a coarse cell to its index. See the glossary entry on coarse cell IDs for more information.
coarse_cell_id | Unique id of the coarse cell. |
Reimplemented from Triangulation< dim, dim >.
Definition at line 359 of file fully_distributed_tria.cc.
|
overrideprivatevirtual |
Translate the index of coarse cell to its unique id. See the glossary entry on coarse cell IDs for more information.
coarse_cell_index | Index of the coarse cell. |
Reimplemented from Triangulation< dim, dim >.
Definition at line 378 of file fully_distributed_tria.cc.
|
private |
store the Settings.
Definition at line 248 of file fully_distributed_tria.h.
|
private |
Partitioner used in copy_triangulation().
Definition at line 255 of file fully_distributed_tria.h.
|
private |
Sorted list of pairs of coarse-cell ids and their indices.
Definition at line 261 of file fully_distributed_tria.h.
|
private |
List of the coarse-cell id for each coarse cell (stored at cell->index()).
Definition at line 268 of file fully_distributed_tria.h.
|
private |
Boolean indicating that the function create_triangulation() was called for internal usage.
Definition at line 274 of file fully_distributed_tria.h.
|
private |
Boolean indicating that the function prepare_coarsening_and_refinement() was called for internal usage.
Definition at line 281 of file fully_distributed_tria.h.