Reference documentation for deal.II version GIT c415589cf0 2022-08-14 18:50: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\}}\)
read_write_vector.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2015 - 2021 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_read_write_vector_h
17 #define dealii_read_write_vector_h
18 
19 #include <deal.II/base/config.h>
20 
22 #include <deal.II/base/index_set.h>
24 #include <deal.II/base/mpi.h>
25 #include <deal.II/base/parallel.h>
28 #include <deal.II/base/types.h>
29 #include <deal.II/base/utilities.h>
30 
32 
33 #include <cstdlib>
34 #include <cstring>
35 #include <iomanip>
36 
37 #ifdef DEAL_II_WITH_TRILINOS
41 
42 # include <Epetra_MultiVector.h>
43 
44 #endif
45 
47 
48 // Forward declarations
49 #ifndef DOXYGEN
50 namespace LinearAlgebra
51 {
52  template <typename>
53  class Vector;
54  namespace distributed
55  {
56  template <typename, typename>
57  class Vector;
58  } // namespace distributed
59 } // namespace LinearAlgebra
60 
61 # ifdef DEAL_II_WITH_PETSC
62 namespace PETScWrappers
63 {
64  namespace MPI
65  {
66  class Vector;
67  }
68 } // namespace PETScWrappers
69 # endif
70 
71 # ifdef DEAL_II_WITH_TRILINOS
72 namespace TrilinosWrappers
73 {
74  namespace MPI
75  {
76  class Vector;
77  }
78 } // namespace TrilinosWrappers
79 # endif
80 
81 # ifdef DEAL_II_WITH_CUDA
82 namespace LinearAlgebra
83 {
84  namespace CUDAWrappers
85  {
86  template <typename>
87  class Vector;
88  }
89 } // namespace LinearAlgebra
90 # endif
91 #endif
92 
93 namespace LinearAlgebra
94 {
133  template <typename Number>
135  {
136  public:
142  using value_type = Number;
143  using pointer = value_type *;
144  using const_pointer = const value_type *;
145  using iterator = value_type *;
146  using const_iterator = const value_type *;
148  using const_reference = const value_type &;
151 
160 
165 
170  explicit ReadWriteVector(const size_type size);
171 
176  explicit ReadWriteVector(const IndexSet &locally_stored_indices);
177 
181  ~ReadWriteVector() override = default;
182 
191  virtual void
192  reinit(const size_type size, const bool omit_zeroing_entries = false);
193 
202  template <typename Number2>
203  void
205  const bool omit_zeroing_entries = false);
206 
216  virtual void
217  reinit(const IndexSet &locally_stored_indices,
218  const bool omit_zeroing_entries = false);
219 
220 
221 #ifdef DEAL_II_WITH_TRILINOS
233  void
234  reinit(const TrilinosWrappers::MPI::Vector &trilinos_vec);
235 #endif
236 
250  template <typename Functor>
251  void
252  apply(const Functor &func);
253 
266  void
268 
274 
278  template <typename Number2>
281 
287  operator=(const Number s);
288 
298  void
299  import(const ::Vector<Number> &vec,
300  VectorOperation::values operation,
301  const std::shared_ptr<const Utilities::MPI::CommunicationPatternBase>
302  &communication_pattern = {});
303 
313  void
314  import(const LinearAlgebra::Vector<Number> &vec,
315  VectorOperation::values operation,
316  const std::shared_ptr<const Utilities::MPI::CommunicationPatternBase>
317  &communication_pattern = {});
318 
327  template <typename MemorySpace>
328  void
330  VectorOperation::values operation,
331  const std::shared_ptr<const Utilities::MPI::CommunicationPatternBase>
332  &communication_pattern = {});
333 
334 #ifdef DEAL_II_WITH_PETSC
343  void
344  import(const PETScWrappers::MPI::Vector &petsc_vec,
345  VectorOperation::values operation,
346  const std::shared_ptr<const Utilities::MPI::CommunicationPatternBase>
347  &communication_pattern = {});
348 #endif
349 
350 #ifdef DEAL_II_WITH_TRILINOS
361  void
362  import(const TrilinosWrappers::MPI::Vector &trilinos_vec,
363  VectorOperation::values operation,
364  const std::shared_ptr<const Utilities::MPI::CommunicationPatternBase>
365  &communication_pattern = {});
366 
367 # ifdef DEAL_II_TRILINOS_WITH_TPETRA
376  void
377  import(const TpetraWrappers::Vector<Number> &tpetra_vec,
378  VectorOperation::values operation,
379  const std::shared_ptr<const Utilities::MPI::CommunicationPatternBase>
380  &communication_pattern = {});
381 # endif
382 
391  void
392  import(const EpetraWrappers::Vector &epetra_vec,
393  VectorOperation::values operation,
394  const std::shared_ptr<const Utilities::MPI::CommunicationPatternBase>
395  &communication_pattern = {});
396 #endif
397 
398 #ifdef DEAL_II_WITH_CUDA
405  void
406  import(const CUDAWrappers::Vector<Number> &cuda_vec,
407  VectorOperation::values operation,
408  const std::shared_ptr<const Utilities::MPI::CommunicationPatternBase>
409  &communication_pattern = {});
410 #endif
411 
420  size_type
421  size() const;
422 
431  size_type
432  n_elements() const;
433 
438  size_type
440 
444  const IndexSet &
446 
452  iterator
453  begin();
454 
460  begin() const;
461 
466  iterator
467  end();
468 
474  end() const;
488  Number
490 
496  Number &
498 
506  Number
508 
516  Number &
518 
534  template <typename Number2>
535  void
536  extract_subvector_to(const std::vector<size_type> &indices,
537  std::vector<Number2> & values) const;
538 
566  template <typename ForwardIterator, typename OutputIterator>
567  void
568  extract_subvector_to(ForwardIterator indices_begin,
569  const ForwardIterator indices_end,
570  OutputIterator values_begin) const;
571 
582  Number
583  local_element(const size_type local_index) const;
584 
595  Number &
596  local_element(const size_type local_index);
609  template <typename Number2>
610  void
611  add(const std::vector<size_type> &indices,
612  const std::vector<Number2> & values);
613 
618  template <typename Number2>
619  void
620  add(const std::vector<size_type> & indices,
622 
628  template <typename Number2>
629  void
631  const size_type *indices,
632  const Number2 * values);
633 
637  void
638  print(std::ostream & out,
639  const unsigned int precision = 3,
640  const bool scientific = true) const;
641 
645  std::size_t
649  protected:
650 #ifdef DEAL_II_WITH_TRILINOS
651 # ifdef DEAL_II_TRILINOS_WITH_TPETRA
657  void
658  import(
659  const Tpetra::Vector<Number, int, types::global_dof_index> &tpetra_vector,
660  const IndexSet & locally_owned_elements,
661  VectorOperation::values operation,
662  const MPI_Comm & mpi_comm,
663  const std::shared_ptr<const Utilities::MPI::CommunicationPatternBase>
664  &communication_pattern);
665 # endif
666 
672  void
673  import(const Epetra_MultiVector &multivector,
674  const IndexSet & locally_owned_elements,
675  VectorOperation::values operation,
676  const MPI_Comm & mpi_comm,
677  const std::shared_ptr<const Utilities::MPI::CommunicationPatternBase>
678  &communication_pattern);
679 #endif
680 
684  unsigned int
686  {
687  // the following will throw an exception if the global_index is not
688  // in the remaining_elements
689  return static_cast<unsigned int>(
691  }
692 
696  void
697  resize_val(const size_type new_allocated_size);
698 
699 #ifdef DEAL_II_WITH_TRILINOS
700 # ifdef DEAL_II_TRILINOS_WITH_TPETRA
706  create_tpetra_comm_pattern(const IndexSet &source_index_set,
707  const MPI_Comm &mpi_comm);
708 # endif
709 
715  create_epetra_comm_pattern(const IndexSet &source_index_set,
716  const MPI_Comm &mpi_comm);
717 #endif
718 
723 
728 
733  std::shared_ptr<Utilities::MPI::CommunicationPatternBase> comm_pattern;
734 
738  std::unique_ptr<Number[], decltype(std::free) *> values;
739 
744  mutable std::shared_ptr<::parallel::internal::TBBPartitioner>
746 
747  // Make all other ReadWriteVector types friends.
748  template <typename Number2>
749  friend class ReadWriteVector;
750 
751  private:
757  template <typename Functor>
759  {
760  public:
765 
769  virtual void
771 
772  private:
777 
781  const Functor &functor;
782  };
783  };
784 
788  /*---------------------------- Inline functions ---------------------------*/
789 
790 #ifndef DOXYGEN
791 
792  template <typename Number>
794  : Subscriptor()
795  , values(nullptr, free)
796  {
797  // virtual functions called in constructors and destructors never use the
798  // override in a derived class
799  // for clarity be explicit on which function is called
801  }
802 
803 
804 
805  template <typename Number>
807  const ReadWriteVector<Number> &v)
808  : Subscriptor()
809  , values(nullptr, free)
810  {
811  this->operator=(v);
812  }
813 
814 
815 
816  template <typename Number>
818  : Subscriptor()
819  , values(nullptr, free)
820  {
821  // virtual functions called in constructors and destructors never use the
822  // override in a derived class
823  // for clarity be explicit on which function is called
825  }
826 
827 
828 
829  template <typename Number>
831  const IndexSet &locally_stored_indices)
832  : Subscriptor()
833  , values(nullptr, free)
834  {
835  // virtual functions called in constructors and destructors never use the
836  // override in a derived class
837  // for clarity be explicit on which function is called
838  ReadWriteVector<Number>::reinit(locally_stored_indices);
839  }
840 
841 
842 
843  template <typename Number>
844  inline typename ReadWriteVector<Number>::size_type
846  {
847  return stored_elements.size();
848  }
849 
850 
851 
852  template <typename Number>
853  inline typename ReadWriteVector<Number>::size_type
855  {
856  return stored_elements.n_elements();
857  }
858 
859 
860 
861  template <typename Number>
862  inline typename ReadWriteVector<Number>::size_type
864  {
865  return stored_elements.n_elements();
866  }
867 
868 
869 
870  template <typename Number>
871  inline const IndexSet &
873  {
874  return stored_elements;
875  }
876 
877 
878 
879  template <typename Number>
880  inline typename ReadWriteVector<Number>::iterator
882  {
883  return values.get();
884  }
885 
886 
887 
888  template <typename Number>
891  {
892  return values.get();
893  }
894 
895 
896 
897  template <typename Number>
898  inline typename ReadWriteVector<Number>::iterator
900  {
901  return values.get() + this->n_elements();
902  }
903 
904 
905 
906  template <typename Number>
909  {
910  return values.get() + this->n_elements();
911  }
912 
913 
914 
915  template <typename Number>
916  inline Number
918  {
920  }
921 
922 
923 
924  template <typename Number>
925  inline Number &
927  {
929  }
930 
931 
932 
933  template <typename Number>
934  inline Number
936  {
937  return operator()(global_index);
938  }
939 
940 
941 
942  template <typename Number>
943  inline Number &
945  {
946  return operator()(global_index);
947  }
948 
949 
950 
951  template <typename Number>
952  template <typename Number2>
953  inline void
955  const std::vector<size_type> &indices,
956  std::vector<Number2> & extracted_values) const
957  {
958  for (size_type i = 0; i < indices.size(); ++i)
959  extracted_values[i] = operator()(indices[i]);
960  }
961 
962 
963 
964  template <typename Number>
965  template <typename ForwardIterator, typename OutputIterator>
966  inline void
968  ForwardIterator indices_begin,
969  const ForwardIterator indices_end,
970  OutputIterator values_begin) const
971  {
972  while (indices_begin != indices_end)
973  {
974  *values_begin = operator()(*indices_begin);
975  indices_begin++;
976  values_begin++;
977  }
978  }
979 
980 
981 
982  template <typename Number>
983  inline Number
984  ReadWriteVector<Number>::local_element(const size_type local_index) const
985  {
986  AssertIndexRange(local_index, this->n_elements());
987 
988  return values[local_index];
989  }
990 
991 
992 
993  template <typename Number>
994  inline Number &
996  {
997  AssertIndexRange(local_index, this->n_elements());
998 
999  return values[local_index];
1000  }
1001 
1002 
1003 
1004  template <typename Number>
1005  template <typename Number2>
1006  inline void
1007  ReadWriteVector<Number>::add(const std::vector<size_type> &indices,
1008  const std::vector<Number2> & values)
1009  {
1010  AssertDimension(indices.size(), values.size());
1011  add(indices.size(), indices.data(), values.data());
1012  }
1013 
1014 
1015 
1016  template <typename Number>
1017  template <typename Number2>
1018  inline void
1019  ReadWriteVector<Number>::add(const std::vector<size_type> & indices,
1020  const ReadWriteVector<Number2> &values)
1021  {
1022  const size_type size = indices.size();
1023  for (size_type i = 0; i < size; ++i)
1024  {
1025  Assert(
1027  ExcMessage(
1028  "The given value is not finite but either infinite or Not A Number (NaN)"));
1029  this->operator()(indices[i]) += values[indices[i]];
1030  }
1031  }
1032 
1033 
1034 
1035  template <typename Number>
1036  template <typename Number2>
1037  inline void
1038  ReadWriteVector<Number>::add(const size_type n_indices,
1039  const size_type *indices,
1040  const Number2 * values_to_add)
1041  {
1042  for (size_type i = 0; i < n_indices; ++i)
1043  {
1044  Assert(
1046  ExcMessage(
1047  "The given value is not finite but either infinite or Not A Number (NaN)"));
1048  this->operator()(indices[i]) += values_to_add[i];
1049  }
1050  }
1051 
1052 
1053 
1054  template <typename Number>
1055  template <typename Functor>
1057  ReadWriteVector<Number> &parent,
1058  const Functor & functor)
1059  : parent(parent)
1060  , functor(functor)
1061  {}
1062 
1063 
1064 
1065  template <typename Number>
1066  template <typename Functor>
1067  void
1069  const size_type begin,
1070  const size_type end)
1071  {
1072  for (size_type i = begin; i < end; ++i)
1073  functor(parent.values[i]);
1074  }
1075 
1076 #endif // ifndef DOXYGEN
1077 
1078 } // end of namespace LinearAlgebra
1079 
1080 
1081 
1089 template <typename Number>
1090 inline void
1093 {
1094  u.swap(v);
1095 }
1096 
1097 
1099 
1100 #endif
size_type size() const
Definition: index_set.h:1635
size_type index_within_set(const size_type global_index) const
Definition: index_set.h:1850
size_type n_elements() const
Definition: index_set.h:1783
FunctorTemplate(ReadWriteVector< Number > &parent, const Functor &functor)
virtual void operator()(const size_type begin, const size_type end)
ReadWriteVector< Number > & operator=(const ReadWriteVector< Number2 > &in_vector)
void reinit(const ReadWriteVector< Number2 > &in_vector, const bool omit_zeroing_entries=false)
ReadWriteVector(const ReadWriteVector< Number > &in_vector)
ReadWriteVector(const size_type size)
std::unique_ptr< Number[], decltype(std::free) * > values
void add(const std::vector< size_type > &indices, const std::vector< Number2 > &values)
Number operator()(const size_type global_index) const
void apply(const Functor &func)
std::size_t memory_consumption() const
TpetraWrappers::CommunicationPattern create_tpetra_comm_pattern(const IndexSet &source_index_set, const MPI_Comm &mpi_comm)
void print(std::ostream &out, const unsigned int precision=3, const bool scientific=true) const
Number & operator()(const size_type global_index)
void add(const std::vector< size_type > &indices, const ReadWriteVector< Number2 > &values)
ReadWriteVector< Number > & operator=(const ReadWriteVector< Number > &in_vector)
ReadWriteVector< Number > & operator=(const Number s)
unsigned int global_to_local(const types::global_dof_index global_index) const
void resize_val(const size_type new_allocated_size)
Number operator[](const size_type global_index) const
Number local_element(const size_type local_index) const
void add(const size_type n_elements, const size_type *indices, const Number2 *values)
EpetraWrappers::CommunicationPattern create_epetra_comm_pattern(const IndexSet &source_index_set, const MPI_Comm &mpi_comm)
Number & local_element(const size_type local_index)
void reinit(const TrilinosWrappers::MPI::Vector &trilinos_vec)
virtual void reinit(const IndexSet &locally_stored_indices, const bool omit_zeroing_entries=false)
void extract_subvector_to(ForwardIterator indices_begin, const ForwardIterator indices_end, OutputIterator values_begin) const
virtual void reinit(const size_type size, const bool omit_zeroing_entries=false)
const IndexSet & get_stored_elements() const
size_type locally_owned_size() const
const_iterator end() const
void extract_subvector_to(const std::vector< size_type > &indices, std::vector< Number2 > &values) const
const_iterator begin() const
~ReadWriteVector() override=default
ReadWriteVector(const IndexSet &locally_stored_indices)
types::global_dof_index size_type
void swap(ReadWriteVector< Number > &v)
size_type n_elements() const
std::shared_ptr< Utilities::MPI::CommunicationPatternBase > comm_pattern
std::shared_ptr<::parallel::internal::TBBPartitioner > thread_loop_partitioner
typename numbers::NumberTraits< Number >::real_type real_type
Number & operator[](const size_type global_index)
Definition: vector.h:110
void swap(LinearAlgebra::ReadWriteVector< Number > &u, LinearAlgebra::ReadWriteVector< Number > &v)
#define DEAL_II_DEPRECATED
Definition: config.h:164
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:442
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:443
#define Assert(cond, exc)
Definition: exceptions.h:1473
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:1667
#define AssertIndexRange(index, range)
Definition: exceptions.h:1732
static ::ExceptionBase & ExcMessage(std::string arg1)
TrilinosWrappers::types::int_type global_index(const Epetra_BlockMap &map, const ::types::global_dof_index i)
void free(T *&pointer)
Definition: cuda.h:97
bool is_finite(const double x)
Definition: numbers.h:553
unsigned int global_dof_index
Definition: types.h:76