Reference documentation for deal.II version GIT relicensing-218-g1f873f3929 2024-03-28 15:00:02+00:00
\(\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\}}\)
Loading...
Searching...
No Matches
partitioner.h
Go to the documentation of this file.
1// ------------------------------------------------------------------------
2//
3// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 2012 - 2023 by the deal.II authors
5//
6// This file is part of the deal.II library.
7//
8// Part of the source code is dual licensed under Apache-2.0 WITH
9// LLVM-exception OR LGPL-2.1-or-later. Detailed license information
10// governing the source code and code contributions can be found in
11// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II.
12//
13// ------------------------------------------------------------------------
14
15#ifndef dealii_partitioner_h
16#define dealii_partitioner_h
17
18#include <deal.II/base/config.h>
19
25#include <deal.II/base/types.h>
26
28
29#include <limits>
30#include <memory>
31
33
34namespace Utilities
35{
36 namespace MPI
37 {
197 {
198 public:
203
208 Partitioner(const unsigned int size);
209
225 const types::global_dof_index ghost_size,
226 const MPI_Comm communicator);
227
235 Partitioner(const IndexSet &locally_owned_indices,
236 const IndexSet &ghost_indices_in,
237 const MPI_Comm communicator_in);
238
246 Partitioner(const IndexSet &locally_owned_indices,
247 const MPI_Comm communicator_in);
248
262 virtual void
263 reinit(const IndexSet &locally_owned_indices,
264 const IndexSet &ghost_indices,
265 const MPI_Comm communicator) override;
266
270 void
271 set_owned_indices(const IndexSet &locally_owned_indices);
272
283 void
285 const IndexSet &larger_ghost_index_set = IndexSet());
286
291 size() const;
292
299 unsigned int
301
308 const IndexSet &
310
316 std::pair<types::global_dof_index, types::global_dof_index>
317 local_range() const;
318
323 bool
324 in_local_range(const types::global_dof_index global_index) const;
325
336 unsigned int
337 global_to_local(const types::global_dof_index global_index) const;
338
348 local_to_global(const unsigned int local_index) const;
349
355 bool
356 is_ghost_entry(const types::global_dof_index global_index) const;
357
361 const IndexSet &
363
368 unsigned int
370
381 const std::vector<std::pair<unsigned int, unsigned int>> &
383
389 const std::vector<std::pair<unsigned int, unsigned int>> &
391
400 const std::vector<std::pair<unsigned int, unsigned int>> &
402
407 unsigned int
409
419 const std::vector<std::pair<unsigned int, unsigned int>> &
421
432 bool
433 is_compatible(const Partitioner &part) const;
434
451 bool
453
458 unsigned int
460
465 unsigned int
467
471 virtual MPI_Comm
472 get_mpi_communicator() const override;
473
479 bool
481
482#ifdef DEAL_II_WITH_MPI
518 template <typename Number, typename MemorySpaceType = MemorySpace::Host>
519 void
521 const unsigned int communication_channel,
522 const ArrayView<const Number, MemorySpaceType> &locally_owned_array,
523 const ArrayView<Number, MemorySpaceType> &temporary_storage,
524 const ArrayView<Number, MemorySpaceType> &ghost_array,
525 std::vector<MPI_Request> &requests) const;
526
544 template <typename Number, typename MemorySpaceType = MemorySpace::Host>
545 void
547 const ArrayView<Number, MemorySpaceType> &ghost_array,
548 std::vector<MPI_Request> &requests) const;
549
587 template <typename Number, typename MemorySpaceType = MemorySpace::Host>
588 void
590 const VectorOperation::values vector_operation,
591 const unsigned int communication_channel,
592 const ArrayView<Number, MemorySpaceType> &ghost_array,
593 const ArrayView<Number, MemorySpaceType> &temporary_storage,
594 std::vector<MPI_Request> &requests) const;
595
630 template <typename Number, typename MemorySpaceType = MemorySpace::Host>
631 void
633 const VectorOperation::values vector_operation,
634 const ArrayView<const Number, MemorySpaceType> &temporary_storage,
635 const ArrayView<Number, MemorySpaceType> &locally_owned_storage,
636 const ArrayView<Number, MemorySpaceType> &ghost_array,
637 std::vector<MPI_Request> &requests) const;
638#endif
639
643 std::size_t
645
651 unsigned int,
652 << "Global index " << arg1
653 << " neither owned nor ghost on proc " << arg2 << '.');
654
659 unsigned int,
660 unsigned int,
661 unsigned int,
662 << "The size of the ghost index array (" << arg1
663 << ") must either equal the number of ghost in the "
664 << "partitioner (" << arg2
665 << ") or be equal in size to a more comprehensive index"
666 << "set which contains " << arg3
667 << " elements for this partitioner.");
668
669 private:
674 void
676
681
686
691 std::pair<types::global_dof_index, types::global_dof_index>
693
699
705
710 std::vector<std::pair<unsigned int, unsigned int>> ghost_targets_data;
711
718 std::vector<std::pair<unsigned int, unsigned int>> import_indices_data;
719
726 // The variable is mutable to enable lazy initialization in
727 // export_to_ghosted_array_start().
728 mutable std::vector<
729 Kokkos::View<unsigned int *, MemorySpace::Default::kokkos_space>>
731
737
742 std::vector<std::pair<unsigned int, unsigned int>> import_targets_data;
743
748 std::vector<unsigned int> import_indices_chunks_by_rank_data;
749
755
761
767 std::vector<std::pair<unsigned int, unsigned int>>
769
773 unsigned int my_pid;
774
778 unsigned int n_procs;
779
784
789 };
790
791
792
793 /*--------------------- Inline functions --------------------------------*/
794
795#ifndef DOXYGEN
796
798 Partitioner::size() const
799 {
800 return global_size;
801 }
802
803
804
805 inline const IndexSet &
807 {
809 }
810
811
812
813 inline std::pair<types::global_dof_index, types::global_dof_index>
815 {
816 return local_range_data;
817 }
818
819
820
821 inline unsigned int
823 {
825 local_range_data.second - local_range_data.first;
826 Assert(size <= std::numeric_limits<unsigned int>::max(),
828 return static_cast<unsigned int>(size);
829 }
830
831
832
833 inline bool
835 const types::global_dof_index global_index) const
836 {
837 return (local_range_data.first <= global_index &&
838 global_index < local_range_data.second);
839 }
840
841
842
843 inline bool
845 const types::global_dof_index global_index) const
846 {
847 // if the index is in the global range, it is trivially not a ghost
848 if (in_local_range(global_index) == true)
849 return false;
850 else
851 return ghost_indices().is_element(global_index);
852 }
853
854
855
856 inline unsigned int
858 const types::global_dof_index global_index) const
859 {
860 Assert(in_local_range(global_index) || is_ghost_entry(global_index),
861 ExcIndexNotPresent(global_index, my_pid));
862 if (in_local_range(global_index))
863 return static_cast<unsigned int>(global_index - local_range_data.first);
864 else
865 {
866 // avoid checking the ghost index set via a binary search twice by
867 // querying the index within set, which returns invalid_dof_index
868 // for non-existent entries
869 const types::global_dof_index index_within_ghosts =
871 if (index_within_ghosts == numbers::invalid_dof_index)
873 else
874 return locally_owned_size() +
875 static_cast<unsigned int>(index_within_ghosts);
876 }
877 }
878
879
880
882 Partitioner::local_to_global(const unsigned int local_index) const
883 {
884 AssertIndexRange(local_index,
886 if (local_index < locally_owned_size())
887 return local_range_data.first + types::global_dof_index(local_index);
888 else
889 return ghost_indices_data.nth_index_in_set(local_index -
891 }
892
893
894
895 inline const IndexSet &
897 {
898 return ghost_indices_data;
899 }
900
901
902
903 inline unsigned int
905 {
907 }
908
909
910
911 inline const std::vector<std::pair<unsigned int, unsigned int>> &
913 {
915 }
916
917
918
919 inline const std::vector<std::pair<unsigned int, unsigned int>> &
921 {
922 return ghost_targets_data;
923 }
924
925
926 inline const std::vector<std::pair<unsigned int, unsigned int>> &
928 {
929 return import_indices_data;
930 }
931
932
933
934 inline unsigned int
936 {
938 }
939
940
941
942 inline const std::vector<std::pair<unsigned int, unsigned int>> &
944 {
945 return import_targets_data;
946 }
947
948
949
950 inline unsigned int
952 {
953 // return the id from the variable stored in this class instead of
954 // Utilities::MPI::this_mpi_process() in order to make this query also
955 // work when MPI is not initialized.
956 return my_pid;
957 }
958
959
960
961 inline unsigned int
963 {
964 // return the number of MPI processes from the variable stored in this
965 // class instead of Utilities::MPI::n_mpi_processes() in order to make
966 // this query also work when MPI is not initialized.
967 return n_procs;
968 }
969
970
971
972 inline MPI_Comm
974 {
975 return communicator;
976 }
977
978
979
980 inline bool
982 {
983 return have_ghost_indices;
984 }
985
986#endif // ifndef DOXYGEN
987
988 } // end of namespace MPI
989
990} // end of namespace Utilities
991
992
994
995#endif
size_type index_within_set(const size_type global_index) const
Definition index_set.h:1990
bool is_element(const size_type index) const
Definition index_set.h:1883
size_type nth_index_in_set(const size_type local_index) const
Definition index_set.h:1971
unsigned int n_ghost_indices_data
Partitioner(const IndexSet &locally_owned_indices, const IndexSet &ghost_indices_in, const MPI_Comm communicator_in)
std::size_t memory_consumption() const
unsigned int global_to_local(const types::global_dof_index global_index) const
bool in_local_range(const types::global_dof_index global_index) const
std::vector< Kokkos::View< unsigned int *, MemorySpace::Default::kokkos_space > > import_indices_plain_dev
std::vector< std::pair< unsigned int, unsigned int > > import_indices_data
bool is_globally_compatible(const Partitioner &part) const
void set_owned_indices(const IndexSet &locally_owned_indices)
const IndexSet & locally_owned_range() const
unsigned int n_ghost_indices_in_larger_set
unsigned int n_ghost_indices() const
std::vector< std::pair< unsigned int, unsigned int > > ghost_targets_data
const std::vector< std::pair< unsigned int, unsigned int > > & ghost_indices_within_larger_ghost_set() const
const IndexSet & ghost_indices() const
const std::vector< std::pair< unsigned int, unsigned int > > & ghost_targets() const
Partitioner(const IndexSet &locally_owned_indices, const MPI_Comm communicator_in)
void export_to_ghosted_array_start(const unsigned int communication_channel, const ArrayView< const Number, MemorySpaceType > &locally_owned_array, const ArrayView< Number, MemorySpaceType > &temporary_storage, const ArrayView< Number, MemorySpaceType > &ghost_array, std::vector< MPI_Request > &requests) const
unsigned int this_mpi_process() const
const std::vector< std::pair< unsigned int, unsigned int > > & import_targets() const
void import_from_ghosted_array_finish(const VectorOperation::values vector_operation, const ArrayView< const Number, MemorySpaceType > &temporary_storage, const ArrayView< Number, MemorySpaceType > &locally_owned_storage, const ArrayView< Number, MemorySpaceType > &ghost_array, std::vector< MPI_Request > &requests) const
unsigned int n_mpi_processes() const
void initialize_import_indices_plain_dev() const
Partitioner(const unsigned int size)
types::global_dof_index global_size
unsigned int locally_owned_size() const
types::global_dof_index local_to_global(const unsigned int local_index) const
std::pair< types::global_dof_index, types::global_dof_index > local_range_data
std::vector< unsigned int > import_indices_chunks_by_rank_data
virtual MPI_Comm get_mpi_communicator() const override
unsigned int n_import_indices_data
void set_ghost_indices(const IndexSet &ghost_indices, const IndexSet &larger_ghost_index_set=IndexSet())
void export_to_ghosted_array_finish(const ArrayView< Number, MemorySpaceType > &ghost_array, std::vector< MPI_Request > &requests) const
Partitioner(const types::global_dof_index local_size, const types::global_dof_index ghost_size, const MPI_Comm communicator)
std::pair< types::global_dof_index, types::global_dof_index > local_range() const
unsigned int n_import_indices() const
const std::vector< std::pair< unsigned int, unsigned int > > & import_indices() const
std::vector< unsigned int > ghost_indices_subset_chunks_by_rank_data
void import_from_ghosted_array_start(const VectorOperation::values vector_operation, const unsigned int communication_channel, const ArrayView< Number, MemorySpaceType > &ghost_array, const ArrayView< Number, MemorySpaceType > &temporary_storage, std::vector< MPI_Request > &requests) const
std::vector< std::pair< unsigned int, unsigned int > > import_targets_data
types::global_dof_index size() const
virtual void reinit(const IndexSet &locally_owned_indices, const IndexSet &ghost_indices, const MPI_Comm communicator) override
bool is_ghost_entry(const types::global_dof_index global_index) const
bool is_compatible(const Partitioner &part) const
std::vector< std::pair< unsigned int, unsigned int > > ghost_indices_subset_data
bool ghost_indices_initialized() const
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:502
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:503
static ::ExceptionBase & ExcNotImplemented()
#define Assert(cond, exc)
static ::ExceptionBase & ExcIndexNotPresent(types::global_dof_index arg1, unsigned int arg2)
#define DeclException2(Exception2, type1, type2, outsequence)
Definition exceptions.h:539
static ::ExceptionBase & ExcGhostIndexArrayHasWrongSize(unsigned int arg1, unsigned int arg2, unsigned int arg3)
#define AssertIndexRange(index, range)
#define DeclException3(Exception3, type1, type2, type3, outsequence)
Definition exceptions.h:562
TrilinosWrappers::types::int_type global_index(const Epetra_BlockMap &map, const ::types::global_dof_index i)
static const unsigned int invalid_unsigned_int
Definition types.h:220
const types::global_dof_index invalid_dof_index
Definition types.h:252
unsigned int global_dof_index
Definition types.h:81