Reference documentation for deal.II version GIT 6a72d26406 2023-06-07 13:05: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\}}\)
partitioner.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2011 - 2022 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 
257  virtual void
258  reinit(const IndexSet &vector_space_vector_index_set,
259  const IndexSet &read_write_vector_index_set,
260  const MPI_Comm communicator) override;
261 
265  void
266  set_owned_indices(const IndexSet &locally_owned_indices);
267 
278  void
280  const IndexSet &larger_ghost_index_set = IndexSet());
281 
286  size() const;
287 
298  unsigned int
299  local_size() const;
300 
307  unsigned int
309 
316  const IndexSet &
318 
324  std::pair<types::global_dof_index, types::global_dof_index>
325  local_range() const;
326 
331  bool
333 
344  unsigned int
346 
356  local_to_global(const unsigned int local_index) const;
357 
363  bool
365 
369  const IndexSet &
370  ghost_indices() const;
371 
376  unsigned int
378 
389  const std::vector<std::pair<unsigned int, unsigned int>> &
391 
397  const std::vector<std::pair<unsigned int, unsigned int>> &
398  ghost_targets() const;
399 
408  const std::vector<std::pair<unsigned int, unsigned int>> &
409  import_indices() const;
410 
415  unsigned int
417 
427  const std::vector<std::pair<unsigned int, unsigned int>> &
428  import_targets() const;
429 
440  bool
441  is_compatible(const Partitioner &part) const;
442 
459  bool
461 
466  unsigned int
468 
473  unsigned int
475 
479  virtual MPI_Comm
480  get_mpi_communicator() const override;
481 
487  bool
489 
490 #ifdef DEAL_II_WITH_MPI
526  template <typename Number, typename MemorySpaceType = MemorySpace::Host>
527  void
529  const unsigned int communication_channel,
530  const ArrayView<const Number, MemorySpaceType> &locally_owned_array,
531  const ArrayView<Number, MemorySpaceType> & temporary_storage,
532  const ArrayView<Number, MemorySpaceType> & ghost_array,
533  std::vector<MPI_Request> & requests) const;
534 
552  template <typename Number, typename MemorySpaceType = MemorySpace::Host>
553  void
555  const ArrayView<Number, MemorySpaceType> &ghost_array,
556  std::vector<MPI_Request> & requests) const;
557 
595  template <typename Number, typename MemorySpaceType = MemorySpace::Host>
596  void
598  const VectorOperation::values vector_operation,
599  const unsigned int communication_channel,
600  const ArrayView<Number, MemorySpaceType> &ghost_array,
601  const ArrayView<Number, MemorySpaceType> &temporary_storage,
602  std::vector<MPI_Request> & requests) const;
603 
638  template <typename Number, typename MemorySpaceType = MemorySpace::Host>
639  void
641  const VectorOperation::values vector_operation,
642  const ArrayView<const Number, MemorySpaceType> &temporary_storage,
643  const ArrayView<Number, MemorySpaceType> & locally_owned_storage,
644  const ArrayView<Number, MemorySpaceType> & ghost_array,
645  std::vector<MPI_Request> & requests) const;
646 #endif
647 
651  std::size_t
653 
659  unsigned int,
660  << "Global index " << arg1
661  << " neither owned nor ghost on proc " << arg2 << '.');
662 
667  unsigned int,
668  unsigned int,
669  unsigned int,
670  << "The size of the ghost index array (" << arg1
671  << ") must either equal the number of ghost in the "
672  << "partitioner (" << arg2
673  << ") or be equal in size to a more comprehensive index"
674  << "set which contains " << arg3
675  << " elements for this partitioner.");
676 
677  private:
682  void
684 
689 
694 
699  std::pair<types::global_dof_index, types::global_dof_index>
701 
707 
712  unsigned int n_ghost_indices_data;
713 
718  std::vector<std::pair<unsigned int, unsigned int>> ghost_targets_data;
719 
726  std::vector<std::pair<unsigned int, unsigned int>> import_indices_data;
727 
734  // The variable is mutable to enable lazy initialization in
735  // export_to_ghosted_array_start().
736  mutable std::vector<
737  Kokkos::View<unsigned int *, MemorySpace::Default::kokkos_space>>
739 
744  unsigned int n_import_indices_data;
745 
750  std::vector<std::pair<unsigned int, unsigned int>> import_targets_data;
751 
756  std::vector<unsigned int> import_indices_chunks_by_rank_data;
757 
763 
768  std::vector<unsigned int> ghost_indices_subset_chunks_by_rank_data;
769 
775  std::vector<std::pair<unsigned int, unsigned int>>
777 
781  unsigned int my_pid;
782 
786  unsigned int n_procs;
787 
791  MPI_Comm communicator;
792 
797  };
798 
799 
800 
801  /*--------------------- Inline functions --------------------------------*/
802 
803 #ifndef DOXYGEN
804 
806  Partitioner::size() const
807  {
808  return global_size;
809  }
810 
811 
812 
813  inline const IndexSet &
815  {
817  }
818 
819 
820 
821  inline std::pair<types::global_dof_index, types::global_dof_index>
823  {
824  return local_range_data;
825  }
826 
827 
828 
829  inline unsigned int
831  {
832  return locally_owned_size();
833  }
834 
835 
836 
837  inline unsigned int
839  {
841  local_range_data.second - local_range_data.first;
844  return static_cast<unsigned int>(size);
845  }
846 
847 
848 
849  inline bool
852  {
853  return (local_range_data.first <= global_index &&
854  global_index < local_range_data.second);
855  }
856 
857 
858 
859  inline bool
862  {
863  // if the index is in the global range, it is trivially not a ghost
864  if (in_local_range(global_index) == true)
865  return false;
866  else
868  }
869 
870 
871 
872  inline unsigned int
875  {
879  return static_cast<unsigned int>(global_index - local_range_data.first);
880  else
881  {
882  // avoid checking the ghost index set via a binary search twice by
883  // querying the index within set, which returns invalid_dof_index
884  // for non-existent entries
885  const types::global_dof_index index_within_ghosts =
887  if (index_within_ghosts == numbers::invalid_dof_index)
889  else
890  return locally_owned_size() +
891  static_cast<unsigned int>(index_within_ghosts);
892  }
893  }
894 
895 
896 
898  Partitioner::local_to_global(const unsigned int local_index) const
899  {
900  AssertIndexRange(local_index,
902  if (local_index < locally_owned_size())
903  return local_range_data.first + types::global_dof_index(local_index);
904  else
905  return ghost_indices_data.nth_index_in_set(local_index -
907  }
908 
909 
910 
911  inline const IndexSet &
913  {
914  return ghost_indices_data;
915  }
916 
917 
918 
919  inline unsigned int
921  {
922  return n_ghost_indices_data;
923  }
924 
925 
926 
927  inline const std::vector<std::pair<unsigned int, unsigned int>> &
929  {
931  }
932 
933 
934 
935  inline const std::vector<std::pair<unsigned int, unsigned int>> &
937  {
938  return ghost_targets_data;
939  }
940 
941 
942  inline const std::vector<std::pair<unsigned int, unsigned int>> &
944  {
945  return import_indices_data;
946  }
947 
948 
949 
950  inline unsigned int
952  {
953  return n_import_indices_data;
954  }
955 
956 
957 
958  inline const std::vector<std::pair<unsigned int, unsigned int>> &
960  {
961  return import_targets_data;
962  }
963 
964 
965 
966  inline unsigned int
968  {
969  // return the id from the variable stored in this class instead of
970  // Utilities::MPI::this_mpi_process() in order to make this query also
971  // work when MPI is not initialized.
972  return my_pid;
973  }
974 
975 
976 
977  inline unsigned int
979  {
980  // return the number of MPI processes from the variable stored in this
981  // class instead of Utilities::MPI::n_mpi_processes() in order to make
982  // this query also work when MPI is not initialized.
983  return n_procs;
984  }
985 
986 
987 
988  inline MPI_Comm
990  {
991  return communicator;
992  }
993 
994 
995 
996  inline bool
998  {
999  return have_ghost_indices;
1000  }
1001 
1002 #endif // ifndef DOXYGEN
1003 
1004  } // end of namespace MPI
1005 
1006 } // end of namespace Utilities
1007 
1008 
1010 
1011 #endif
size_type index_within_set(const size_type global_index) const
Definition: index_set.h:1883
bool is_element(const size_type index) const
Definition: index_set.h:1776
size_type nth_index_in_set(const size_type local_index) const
Definition: index_set.h:1864
const std::vector< std::pair< unsigned int, unsigned int > > & ghost_targets() const
unsigned int n_ghost_indices_data
Definition: partitioner.h:712
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:738
std::vector< std::pair< unsigned int, unsigned int > > import_indices_data
Definition: partitioner.h:726
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:762
unsigned int n_ghost_indices() const
std::vector< std::pair< unsigned int, unsigned int > > ghost_targets_data
Definition: partitioner.h:718
unsigned int local_size() const
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:688
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:700
std::vector< unsigned int > import_indices_chunks_by_rank_data
Definition: partitioner.h:756
virtual MPI_Comm get_mpi_communicator() const override
unsigned int n_import_indices_data
Definition: partitioner.h:744
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:768
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:750
types::global_dof_index size() const
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
virtual void reinit(const IndexSet &vector_space_vector_index_set, const IndexSet &read_write_vector_index_set, const MPI_Comm communicator) override
std::vector< std::pair< unsigned int, unsigned int > > ghost_indices_subset_data
Definition: partitioner.h:776
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_DEPRECATED
Definition: config.h:174
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:474
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:475
static ::ExceptionBase & ExcGhostIndexArrayHasWrongSize(unsigned int arg1, unsigned int arg2, unsigned int arg3)
#define Assert(cond, exc)
Definition: exceptions.h:1614
static ::ExceptionBase & ExcNotImplemented()
#define DeclException2(Exception2, type1, type2, outsequence)
Definition: exceptions.h:533
#define AssertIndexRange(index, range)
Definition: exceptions.h:1855
#define DeclException3(Exception3, type1, type2, type3, outsequence)
Definition: exceptions.h:556
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