Reference documentation for deal.II version GIT c14369f203 2023-10-01 07:40: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\}}\)
partitioner.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2011 - 2023 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 #ifndef dealii_partitioner_h
17 #define dealii_partitioner_h
18 
19 #include <deal.II/base/config.h>
20 
23 #include <deal.II/base/index_set.h>
25 #include <deal.II/base/mpi_stub.h>
26 #include <deal.II/base/types.h>
27 
29 
30 #include <limits>
31 #include <memory>
32 
34 
35 namespace Utilities
36 {
37  namespace MPI
38  {
198  {
199  public:
204 
209  Partitioner(const unsigned int size);
210 
226  const types::global_dof_index ghost_size,
227  const MPI_Comm communicator);
228 
236  Partitioner(const IndexSet &locally_owned_indices,
237  const IndexSet &ghost_indices_in,
238  const MPI_Comm communicator_in);
239 
247  Partitioner(const IndexSet &locally_owned_indices,
248  const MPI_Comm communicator_in);
249 
263  virtual void
264  reinit(const IndexSet &locally_owned_indices,
265  const IndexSet &ghost_indices,
266  const MPI_Comm communicator) override;
267 
271  void
272  set_owned_indices(const IndexSet &locally_owned_indices);
273 
284  void
286  const IndexSet &larger_ghost_index_set = IndexSet());
287 
292  size() const;
293 
300  unsigned int
302 
309  const IndexSet &
311 
317  std::pair<types::global_dof_index, types::global_dof_index>
318  local_range() const;
319 
324  bool
326 
337  unsigned int
339 
349  local_to_global(const unsigned int local_index) const;
350 
356  bool
358 
362  const IndexSet &
363  ghost_indices() const;
364 
369  unsigned int
371 
382  const std::vector<std::pair<unsigned int, unsigned int>> &
384 
390  const std::vector<std::pair<unsigned int, unsigned int>> &
391  ghost_targets() const;
392 
401  const std::vector<std::pair<unsigned int, unsigned int>> &
402  import_indices() const;
403 
408  unsigned int
410 
420  const std::vector<std::pair<unsigned int, unsigned int>> &
421  import_targets() const;
422 
433  bool
434  is_compatible(const Partitioner &part) const;
435 
452  bool
454 
459  unsigned int
461 
466  unsigned int
468 
472  virtual MPI_Comm
473  get_mpi_communicator() const override;
474 
480  bool
482 
483 #ifdef DEAL_II_WITH_MPI
519  template <typename Number, typename MemorySpaceType = MemorySpace::Host>
520  void
522  const unsigned int communication_channel,
523  const ArrayView<const Number, MemorySpaceType> &locally_owned_array,
524  const ArrayView<Number, MemorySpaceType> &temporary_storage,
525  const ArrayView<Number, MemorySpaceType> &ghost_array,
526  std::vector<MPI_Request> &requests) const;
527 
545  template <typename Number, typename MemorySpaceType = MemorySpace::Host>
546  void
548  const ArrayView<Number, MemorySpaceType> &ghost_array,
549  std::vector<MPI_Request> &requests) const;
550 
588  template <typename Number, typename MemorySpaceType = MemorySpace::Host>
589  void
591  const VectorOperation::values vector_operation,
592  const unsigned int communication_channel,
593  const ArrayView<Number, MemorySpaceType> &ghost_array,
594  const ArrayView<Number, MemorySpaceType> &temporary_storage,
595  std::vector<MPI_Request> &requests) const;
596 
631  template <typename Number, typename MemorySpaceType = MemorySpace::Host>
632  void
634  const VectorOperation::values vector_operation,
635  const ArrayView<const Number, MemorySpaceType> &temporary_storage,
636  const ArrayView<Number, MemorySpaceType> &locally_owned_storage,
637  const ArrayView<Number, MemorySpaceType> &ghost_array,
638  std::vector<MPI_Request> &requests) const;
639 #endif
640 
644  std::size_t
646 
652  unsigned int,
653  << "Global index " << arg1
654  << " neither owned nor ghost on proc " << arg2 << '.');
655 
660  unsigned int,
661  unsigned int,
662  unsigned int,
663  << "The size of the ghost index array (" << arg1
664  << ") must either equal the number of ghost in the "
665  << "partitioner (" << arg2
666  << ") or be equal in size to a more comprehensive index"
667  << "set which contains " << arg3
668  << " elements for this partitioner.");
669 
670  private:
675  void
677 
682 
687 
692  std::pair<types::global_dof_index, types::global_dof_index>
694 
700 
705  unsigned int n_ghost_indices_data;
706 
711  std::vector<std::pair<unsigned int, unsigned int>> ghost_targets_data;
712 
719  std::vector<std::pair<unsigned int, unsigned int>> import_indices_data;
720 
727  // The variable is mutable to enable lazy initialization in
728  // export_to_ghosted_array_start().
729  mutable std::vector<
730  Kokkos::View<unsigned int *, MemorySpace::Default::kokkos_space>>
732 
737  unsigned int n_import_indices_data;
738 
743  std::vector<std::pair<unsigned int, unsigned int>> import_targets_data;
744 
749  std::vector<unsigned int> import_indices_chunks_by_rank_data;
750 
756 
761  std::vector<unsigned int> ghost_indices_subset_chunks_by_rank_data;
762 
768  std::vector<std::pair<unsigned int, unsigned int>>
770 
774  unsigned int my_pid;
775 
779  unsigned int n_procs;
780 
784  MPI_Comm communicator;
785 
790  };
791 
792 
793 
794  /*--------------------- Inline functions --------------------------------*/
795 
796 #ifndef DOXYGEN
797 
799  Partitioner::size() const
800  {
801  return global_size;
802  }
803 
804 
805 
806  inline const IndexSet &
808  {
810  }
811 
812 
813 
814  inline std::pair<types::global_dof_index, types::global_dof_index>
816  {
817  return local_range_data;
818  }
819 
820 
821 
822  inline unsigned int
824  {
826  local_range_data.second - local_range_data.first;
829  return static_cast<unsigned int>(size);
830  }
831 
832 
833 
834  inline bool
837  {
838  return (local_range_data.first <= global_index &&
839  global_index < local_range_data.second);
840  }
841 
842 
843 
844  inline bool
847  {
848  // if the index is in the global range, it is trivially not a ghost
849  if (in_local_range(global_index) == true)
850  return false;
851  else
853  }
854 
855 
856 
857  inline unsigned int
860  {
864  return static_cast<unsigned int>(global_index - local_range_data.first);
865  else
866  {
867  // avoid checking the ghost index set via a binary search twice by
868  // querying the index within set, which returns invalid_dof_index
869  // for non-existent entries
870  const types::global_dof_index index_within_ghosts =
872  if (index_within_ghosts == numbers::invalid_dof_index)
874  else
875  return locally_owned_size() +
876  static_cast<unsigned int>(index_within_ghosts);
877  }
878  }
879 
880 
881 
883  Partitioner::local_to_global(const unsigned int local_index) const
884  {
885  AssertIndexRange(local_index,
887  if (local_index < locally_owned_size())
888  return local_range_data.first + types::global_dof_index(local_index);
889  else
890  return ghost_indices_data.nth_index_in_set(local_index -
892  }
893 
894 
895 
896  inline const IndexSet &
898  {
899  return ghost_indices_data;
900  }
901 
902 
903 
904  inline unsigned int
906  {
907  return n_ghost_indices_data;
908  }
909 
910 
911 
912  inline const std::vector<std::pair<unsigned int, unsigned int>> &
914  {
916  }
917 
918 
919 
920  inline const std::vector<std::pair<unsigned int, unsigned int>> &
922  {
923  return ghost_targets_data;
924  }
925 
926 
927  inline const std::vector<std::pair<unsigned int, unsigned int>> &
929  {
930  return import_indices_data;
931  }
932 
933 
934 
935  inline unsigned int
937  {
938  return n_import_indices_data;
939  }
940 
941 
942 
943  inline const std::vector<std::pair<unsigned int, unsigned int>> &
945  {
946  return import_targets_data;
947  }
948 
949 
950 
951  inline unsigned int
953  {
954  // return the id from the variable stored in this class instead of
955  // Utilities::MPI::this_mpi_process() in order to make this query also
956  // work when MPI is not initialized.
957  return my_pid;
958  }
959 
960 
961 
962  inline unsigned int
964  {
965  // return the number of MPI processes from the variable stored in this
966  // class instead of Utilities::MPI::n_mpi_processes() in order to make
967  // this query also work when MPI is not initialized.
968  return n_procs;
969  }
970 
971 
972 
973  inline MPI_Comm
975  {
976  return communicator;
977  }
978 
979 
980 
981  inline bool
983  {
984  return have_ghost_indices;
985  }
986 
987 #endif // ifndef DOXYGEN
988 
989  } // end of namespace MPI
990 
991 } // end of namespace Utilities
992 
993 
995 
996 #endif
size_type index_within_set(const size_type global_index) const
Definition: index_set.h:1921
bool is_element(const size_type index) const
Definition: index_set.h:1814
size_type nth_index_in_set(const size_type local_index) const
Definition: index_set.h:1902
const std::vector< std::pair< unsigned int, unsigned int > > & ghost_targets() const
unsigned int n_ghost_indices_data
Definition: partitioner.h:705
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
Definition: partitioner.h:731
std::vector< std::pair< unsigned int, unsigned int > > import_indices_data
Definition: partitioner.h:719
bool is_globally_compatible(const Partitioner &part) const
void set_owned_indices(const IndexSet &locally_owned_indices)
const std::vector< std::pair< unsigned int, unsigned int > > & import_indices() const
unsigned int n_ghost_indices_in_larger_set
Definition: partitioner.h:755
unsigned int n_ghost_indices() const
std::vector< std::pair< unsigned int, unsigned int > > ghost_targets_data
Definition: partitioner.h:711
Partitioner(const IndexSet &locally_owned_indices, const MPI_Comm communicator_in)
const std::vector< std::pair< unsigned int, unsigned int > > & ghost_indices_within_larger_ghost_set() const
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
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
Definition: partitioner.h:681
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
Definition: partitioner.h:693
std::vector< unsigned int > import_indices_chunks_by_rank_data
Definition: partitioner.h:749
virtual MPI_Comm get_mpi_communicator() const override
unsigned int n_import_indices_data
Definition: partitioner.h:737
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)
unsigned int n_import_indices() const
std::vector< unsigned int > ghost_indices_subset_chunks_by_rank_data
Definition: partitioner.h:761
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
Definition: partitioner.h:743
types::global_dof_index size() const
virtual void reinit(const IndexSet &locally_owned_indices, const IndexSet &ghost_indices, const MPI_Comm communicator) override
const IndexSet & locally_owned_range() const
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
Definition: partitioner.h:769
bool ghost_indices_initialized() const
std::pair< types::global_dof_index, types::global_dof_index > local_range() const
const std::vector< std::pair< unsigned int, unsigned int > > & import_targets() const
const IndexSet & ghost_indices() const
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:477
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:478
static ::ExceptionBase & ExcGhostIndexArrayHasWrongSize(unsigned int arg1, unsigned int arg2, unsigned int arg3)
#define Assert(cond, exc)
Definition: exceptions.h:1616
static ::ExceptionBase & ExcNotImplemented()
#define DeclException2(Exception2, type1, type2, outsequence)
Definition: exceptions.h:535
#define AssertIndexRange(index, range)
Definition: exceptions.h:1857
#define DeclException3(Exception3, type1, type2, type3, outsequence)
Definition: exceptions.h:558
static ::ExceptionBase & ExcIndexNotPresent(types::global_dof_index arg1, unsigned int arg2)
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:213
const types::global_dof_index invalid_dof_index
Definition: types.h:233
unsigned int global_dof_index
Definition: types.h:82