Loading [MathJax]/extensions/TeX/newcommand.js
 deal.II version GIT relicensing-3087-ga35b476a78 2025-04-19 20:30:01+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\}}
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages Concepts
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
26#include <deal.II/base/types.h>
27
29
30#include <Kokkos_Core.hpp>
31
32#include <limits>
33#include <memory>
34#include <vector>
35
36
38
39namespace Utilities
40{
41 namespace MPI
42 {
202 {
203 public:
208
213 Partitioner(const unsigned int size);
214
230 const types::global_dof_index ghost_size,
231 const MPI_Comm communicator);
232
240 Partitioner(const IndexSet &locally_owned_indices,
241 const IndexSet &ghost_indices_in,
242 const MPI_Comm communicator_in);
243
251 Partitioner(const IndexSet &locally_owned_indices,
252 const MPI_Comm communicator_in);
253
267 virtual void
268 reinit(const IndexSet &locally_owned_indices,
269 const IndexSet &ghost_indices,
270 const MPI_Comm communicator) override;
271
275 void
276 set_owned_indices(const IndexSet &locally_owned_indices);
277
288 void
290 const IndexSet &larger_ghost_index_set = IndexSet());
291
296 size() const;
297
304 unsigned int
306
313 const IndexSet &
315
321 std::pair<types::global_dof_index, types::global_dof_index>
322 local_range() const;
323
328 bool
329 in_local_range(const types::global_dof_index global_index) const;
330
341 unsigned int
342 global_to_local(const types::global_dof_index global_index) const;
343
353 local_to_global(const unsigned int local_index) const;
354
360 bool
361 is_ghost_entry(const types::global_dof_index global_index) const;
362
366 const IndexSet &
368
373 unsigned int
375
386 const std::vector<std::pair<unsigned int, unsigned int>> &
388
394 const std::vector<std::pair<unsigned int, unsigned int>> &
396
405 const std::vector<std::pair<unsigned int, unsigned int>> &
407
412 unsigned int
414
424 const std::vector<std::pair<unsigned int, unsigned int>> &
426
437 bool
438 is_compatible(const Partitioner &part) const;
439
456 bool
458
463 unsigned int
465
470 unsigned int
472
476 virtual MPI_Comm
477 get_mpi_communicator() const override;
478
484 bool
486
487#ifdef DEAL_II_WITH_MPI
523 template <typename Number, typename MemorySpaceType = MemorySpace::Host>
524 void
526 const unsigned int communication_channel,
527 const ArrayView<const Number, MemorySpaceType> &locally_owned_array,
528 const ArrayView<Number, MemorySpaceType> &temporary_storage,
529 const ArrayView<Number, MemorySpaceType> &ghost_array,
530 std::vector<MPI_Request> &requests) const;
531
549 template <typename Number, typename MemorySpaceType = MemorySpace::Host>
550 void
552 const ArrayView<Number, MemorySpaceType> &ghost_array,
553 std::vector<MPI_Request> &requests) const;
554
592 template <typename Number, typename MemorySpaceType = MemorySpace::Host>
593 void
595 const VectorOperation::values vector_operation,
596 const unsigned int communication_channel,
597 const ArrayView<Number, MemorySpaceType> &ghost_array,
598 const ArrayView<Number, MemorySpaceType> &temporary_storage,
599 std::vector<MPI_Request> &requests) const;
600
635 template <typename Number, typename MemorySpaceType = MemorySpace::Host>
636 void
638 const VectorOperation::values vector_operation,
639 const ArrayView<const Number, MemorySpaceType> &temporary_storage,
640 const ArrayView<Number, MemorySpaceType> &locally_owned_storage,
641 const ArrayView<Number, MemorySpaceType> &ghost_array,
642 std::vector<MPI_Request> &requests) const;
643#endif
644
648 std::size_t
650
656 unsigned int,
657 << "Global index " << arg1
658 << " neither owned nor ghost on proc " << arg2 << '.');
659
664 unsigned int,
665 unsigned int,
666 unsigned int,
667 << "The size of the ghost index array (" << arg1
668 << ") must either equal the number of ghost in the "
669 << "partitioner (" << arg2
670 << ") or be equal in size to a more comprehensive index"
671 << "set which contains " << arg3
672 << " elements for this partitioner.");
673
674 private:
679 void
681
686
691
696 std::pair<types::global_dof_index, types::global_dof_index>
698
704
710
715 std::vector<std::pair<unsigned int, unsigned int>> ghost_targets_data;
716
723 std::vector<std::pair<unsigned int, unsigned int>> import_indices_data;
724
731 // The variable is mutable to enable lazy initialization in
732 // export_to_ghosted_array_start().
733 mutable std::vector<
734 Kokkos::View<unsigned int *, MemorySpace::Default::kokkos_space>>
736
742
747 std::vector<std::pair<unsigned int, unsigned int>> import_targets_data;
748
753 std::vector<unsigned int> import_indices_chunks_by_rank_data;
754
760
766
772 std::vector<std::pair<unsigned int, unsigned int>>
774
778 unsigned int my_pid;
779
783 unsigned int n_procs;
784
789
794 };
795
796
797
798 /*--------------------- Inline functions --------------------------------*/
799
800#ifndef DOXYGEN
801
803 Partitioner::size() const
804 {
805 return global_size;
806 }
807
808
809
810 inline const IndexSet &
812 {
814 }
815
816
817
818 inline std::pair<types::global_dof_index, types::global_dof_index>
820 {
821 return local_range_data;
822 }
823
824
825
826 inline unsigned int
828 {
830 local_range_data.second - local_range_data.first;
831 Assert(size <= std::numeric_limits<unsigned int>::max(),
833 return static_cast<unsigned int>(size);
834 }
835
836
837
838 inline bool
840 const types::global_dof_index global_index) const
841 {
842 return (local_range_data.first <= global_index &&
843 global_index < local_range_data.second);
844 }
845
846
847
848 inline bool
850 const types::global_dof_index global_index) const
851 {
852 // if the index is in the global range, it is trivially not a ghost
853 if (in_local_range(global_index) == true)
854 return false;
855 else
856 return ghost_indices().is_element(global_index);
857 }
858
859
860
861 inline unsigned int
863 const types::global_dof_index global_index) const
864 {
865 Assert(in_local_range(global_index) || is_ghost_entry(global_index),
866 ExcIndexNotPresent(global_index, my_pid));
867 if (in_local_range(global_index))
868 return static_cast<unsigned int>(global_index - local_range_data.first);
869 else
870 {
871 // avoid checking the ghost index set via a binary search twice by
872 // querying the index within set, which returns invalid_dof_index
873 // for non-existent entries
874 const types::global_dof_index index_within_ghosts =
876 if (index_within_ghosts == numbers::invalid_dof_index)
878 else
879 return locally_owned_size() +
880 static_cast<unsigned int>(index_within_ghosts);
881 }
882 }
883
884
885
887 Partitioner::local_to_global(const unsigned int local_index) const
888 {
889 AssertIndexRange(local_index,
891 if (local_index < locally_owned_size())
892 return local_range_data.first + types::global_dof_index(local_index);
893 else
894 return ghost_indices_data.nth_index_in_set(local_index -
896 }
897
898
899
900 inline const IndexSet &
902 {
903 return ghost_indices_data;
904 }
905
906
907
908 inline unsigned int
910 {
912 }
913
914
915
916 inline const std::vector<std::pair<unsigned int, unsigned int>> &
918 {
920 }
921
922
923
924 inline const std::vector<std::pair<unsigned int, unsigned int>> &
926 {
927 return ghost_targets_data;
928 }
929
930
931 inline const std::vector<std::pair<unsigned int, unsigned int>> &
933 {
934 return import_indices_data;
935 }
936
937
938
939 inline unsigned int
941 {
943 }
944
945
946
947 inline const std::vector<std::pair<unsigned int, unsigned int>> &
949 {
950 return import_targets_data;
951 }
952
953
954
955 inline unsigned int
957 {
958 // return the id from the variable stored in this class instead of
959 // Utilities::MPI::this_mpi_process() in order to make this query also
960 // work when MPI is not initialized.
961 return my_pid;
962 }
963
964
965
966 inline unsigned int
968 {
969 // return the number of MPI processes from the variable stored in this
970 // class instead of Utilities::MPI::n_mpi_processes() in order to make
971 // this query also work when MPI is not initialized.
972 return n_procs;
973 }
974
975
976
977 inline MPI_Comm
979 {
980 return communicator;
981 }
982
983
984
985 inline bool
987 {
988 return have_ghost_indices;
989 }
990
991#endif // ifndef DOXYGEN
992
993 } // end of namespace MPI
994
995} // end of namespace Utilities
996
997
999
1000#endif
size_type index_within_set(const size_type global_index) const
Definition index_set.h:2013
bool is_element(const size_type index) const
Definition index_set.h:1905
size_type nth_index_in_set(const size_type local_index) const
Definition index_set.h:1994
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:35
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:36
static ::ExceptionBase & ExcNotImplemented()
#define Assert(cond, exc)
static ::ExceptionBase & ExcIndexNotPresent(types::global_dof_index arg1, unsigned int arg2)
#define DeclException2(Exception2, type1, type2, outsequence)
static ::ExceptionBase & ExcGhostIndexArrayHasWrongSize(unsigned int arg1, unsigned int arg2, unsigned int arg3)
#define AssertIndexRange(index, range)
#define DeclException3(Exception3, type1, type2, type3, outsequence)
TrilinosWrappers::types::int_type global_index(const Epetra_BlockMap &map, const ::types::global_dof_index i)
constexpr types::global_dof_index invalid_dof_index
Definition types.h:269
constexpr unsigned int invalid_unsigned_int
Definition types.h:238
unsigned int global_dof_index
Definition types.h:94