Reference documentation for deal.II version GIT relicensing-214-g6e74dec06b 2024-03-27 18:10: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\}}\)
Loading...
Searching...
No Matches
read_write_vector.h
Go to the documentation of this file.
1// ------------------------------------------------------------------------
2//
3// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 2015 - 2024 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_read_write_vector_h
16#define dealii_read_write_vector_h
17
18#include <deal.II/base/config.h>
19
28#include <deal.II/base/types.h>
29
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
50template <typename>
51class Vector;
52
53namespace LinearAlgebra
54{
55 template <typename>
56 class Vector;
57 namespace distributed
58 {
59 template <typename, typename>
60 class Vector;
61 } // namespace distributed
62} // namespace LinearAlgebra
63
64# ifdef DEAL_II_WITH_PETSC
65namespace PETScWrappers
66{
67 namespace MPI
68 {
69 class Vector;
70 }
71} // namespace PETScWrappers
72# endif
73
74# ifdef DEAL_II_WITH_TRILINOS
75namespace TrilinosWrappers
76{
77 namespace MPI
78 {
79 class Vector;
80 }
81} // namespace TrilinosWrappers
82# endif
83
84# ifdef DEAL_II_WITH_CUDA
85namespace LinearAlgebra
86{
87 namespace CUDAWrappers
88 {
89 template <typename>
90 class Vector;
91 }
92} // namespace LinearAlgebra
93# endif
94#endif
95
96namespace LinearAlgebra
97{
135 template <typename Number>
136 class ReadWriteVector : public Subscriptor, public ReadVector<Number>
137 {
138 public:
144 using value_type = Number;
146 using const_pointer = const value_type *;
148 using const_iterator = const value_type *;
153
162
167
173
178 explicit ReadWriteVector(const IndexSet &locally_stored_indices);
179
183 ~ReadWriteVector() override = default;
184
193 virtual void
194 reinit(const size_type size, const bool omit_zeroing_entries = false);
195
204 template <typename Number2>
205 void
207 const bool omit_zeroing_entries = false);
208
218 virtual void
219 reinit(const IndexSet &locally_stored_indices,
220 const bool omit_zeroing_entries = false);
221
222
223#ifdef DEAL_II_WITH_TRILINOS
235 void
237#endif
238
252 template <typename Functor>
253 void
254 apply(const Functor &func);
255
268 void
270
276
280 template <typename Number2>
283
289 operator=(const Number s);
290
300 void
302 const ::Vector<Number> &vec,
303 VectorOperation::values operation,
304 const std::shared_ptr<const Utilities::MPI::CommunicationPatternBase>
305 &communication_pattern = {});
306
311 void
312 import(const ::Vector<Number> &V,
313 VectorOperation::values operation,
314 const std::shared_ptr<const Utilities::MPI::CommunicationPatternBase>
315 &communication_pattern = {})
316 {
317 import_elements(V, operation, communication_pattern);
318 }
319
328 template <typename MemorySpace>
329 void
332 VectorOperation::values operation,
333 const std::shared_ptr<const Utilities::MPI::CommunicationPatternBase>
334 &communication_pattern = {});
335
339 template <typename MemorySpace>
342 VectorOperation::values operation,
343 const std::shared_ptr<const Utilities::MPI::CommunicationPatternBase>
344 &communication_pattern = {})
345 {
346 import_elements(V, operation, communication_pattern);
347 }
348
349
350#ifdef DEAL_II_WITH_PETSC
359 void
361 const PETScWrappers::MPI::Vector &petsc_vec,
362 VectorOperation::values operation,
363 const std::shared_ptr<const Utilities::MPI::CommunicationPatternBase>
364 &communication_pattern = {});
365
370 void
372 VectorOperation::values operation,
373 const std::shared_ptr<const Utilities::MPI::CommunicationPatternBase>
374 &communication_pattern = {})
375 {
376 import_elements(V, operation, communication_pattern);
377 }
378#endif
379
380#ifdef DEAL_II_WITH_TRILINOS
391 void
393 const TrilinosWrappers::MPI::Vector &trilinos_vec,
394 VectorOperation::values operation,
395 const std::shared_ptr<const Utilities::MPI::CommunicationPatternBase>
396 &communication_pattern = {});
397
402 void
404 VectorOperation::values operation,
405 const std::shared_ptr<const Utilities::MPI::CommunicationPatternBase>
406 &communication_pattern = {})
407 {
408 import_elements(V, operation, communication_pattern);
409 }
410
411# ifdef DEAL_II_TRILINOS_WITH_TPETRA
420 template <typename Dummy = Number>
421 std::enable_if_t<std::is_same_v<Dummy, Number> &&
424 const TpetraWrappers::Vector<Number> &tpetra_vec,
425 VectorOperation::values operation,
426 const std::shared_ptr<const Utilities::MPI::CommunicationPatternBase>
427 &communication_pattern = {});
428
432 template <typename Dummy = Number>
433 DEAL_II_DEPRECATED std::enable_if_t<std::is_same_v<Dummy, Number> &&
436 VectorOperation::values operation,
437 const std::shared_ptr<const Utilities::MPI::CommunicationPatternBase>
438 &communication_pattern = {})
439 {
440 import_elements(V, operation, communication_pattern);
441 }
442# endif
443
452 void
454 const EpetraWrappers::Vector &epetra_vec,
455 VectorOperation::values operation,
456 const std::shared_ptr<const Utilities::MPI::CommunicationPatternBase>
457 &communication_pattern = {});
458
463 void
464 import(const EpetraWrappers::Vector &V,
465 VectorOperation::values operation,
466 const std::shared_ptr<const Utilities::MPI::CommunicationPatternBase>
467 &communication_pattern = {})
468 {
469 import_elements(V, operation, communication_pattern);
470 }
471#endif
472
473#ifdef DEAL_II_WITH_CUDA
480 void
482 const CUDAWrappers::Vector<Number> &cuda_vec,
483 VectorOperation::values operation,
484 const std::shared_ptr<const Utilities::MPI::CommunicationPatternBase>
485 &communication_pattern = {});
486
487
492 void
494 VectorOperation::values operation,
495 const std::shared_ptr<const Utilities::MPI::CommunicationPatternBase>
496 &communication_pattern = {})
497 {
498 import_elements(V, operation, communication_pattern);
499 }
500#endif
501
511 size() const override;
512
519
523 const IndexSet &
525
533
539 begin() const;
540
547
553 end() const;
567 Number
568 operator()(const size_type global_index) const;
569
575 Number &
576 operator()(const size_type global_index);
577
585 Number
586 operator[](const size_type global_index) const;
587
595 Number &
596 operator[](const size_type global_index);
597
613 template <typename Number2>
614 void
615 extract_subvector_to(const std::vector<size_type> &indices,
616 std::vector<Number2> &values) const;
617
621 virtual void
624 ArrayView<Number> &entries) const override;
625
653 template <typename ForwardIterator, typename OutputIterator>
654 void
655 extract_subvector_to(ForwardIterator indices_begin,
656 const ForwardIterator indices_end,
657 OutputIterator values_begin) const;
658
669 Number
670 local_element(const size_type local_index) const;
671
682 Number &
683 local_element(const size_type local_index);
696 template <typename Number2>
697 void
698 add(const std::vector<size_type> &indices,
699 const std::vector<Number2> &values);
700
705 template <typename Number2>
706 void
707 add(const std::vector<size_type> &indices,
709
715 template <typename Number2>
716 void
717 add(const size_type n_elements,
718 const size_type *indices,
719 const Number2 *values);
720
724 void
725 print(std::ostream &out,
726 const unsigned int precision = 3,
727 const bool scientific = true) const;
728
732 std::size_t
736 protected:
737#ifdef DEAL_II_WITH_TRILINOS
738# ifdef DEAL_II_TRILINOS_WITH_TPETRA
744 template <typename Dummy = Number>
745 std::enable_if_t<std::is_same_v<Dummy, Number> &&
748 const Tpetra::Vector<Number, int, types::signed_global_dof_index>
749 &tpetra_vector,
750 const IndexSet &locally_owned_elements,
751 VectorOperation::values operation,
752 const MPI_Comm mpi_comm,
753 const std::shared_ptr<const Utilities::MPI::CommunicationPatternBase>
754 &communication_pattern);
755# endif
756
762 void
764 const Epetra_MultiVector &multivector,
765 const IndexSet &locally_owned_elements,
766 VectorOperation::values operation,
767 const MPI_Comm mpi_comm,
768 const std::shared_ptr<const Utilities::MPI::CommunicationPatternBase>
769 &communication_pattern);
770#endif
771
775 unsigned int
776 global_to_local(const types::global_dof_index global_index) const;
777
778#ifdef DEAL_II_WITH_TRILINOS
779# ifdef DEAL_II_TRILINOS_WITH_TPETRA
785 create_tpetra_comm_pattern(const IndexSet &source_index_set,
786 const MPI_Comm mpi_comm);
787# endif
788
794 create_epetra_comm_pattern(const IndexSet &source_index_set,
795 const MPI_Comm mpi_comm);
796#endif
797
802
807
812 std::shared_ptr<Utilities::MPI::CommunicationPatternBase> comm_pattern;
813
818
823 mutable std::shared_ptr<::parallel::internal::TBBPartitioner>
825
826 // Make all other ReadWriteVector types friends.
827 template <typename Number2>
828 friend class ReadWriteVector;
829
830 private:
836 template <typename Functor>
838 {
839 public:
844
848 virtual void
850
851 private:
856
860 const Functor &functor;
861 };
862 };
863
867 /*---------------------------- Inline functions ---------------------------*/
868
869#ifndef DOXYGEN
870
871
872 template <typename Number>
874 {
875 // virtual functions called in constructors and destructors never use the
876 // override in a derived class
877 // for clarity be explicit on which function is called
879 }
880
881
882
883 template <typename Number>
885 const ReadWriteVector<Number> &v)
886 : Subscriptor()
887 {
888 this->operator=(v);
889 }
890
891
892
893 template <typename Number>
895 {
896 // virtual functions called in constructors and destructors never use the
897 // override in a derived class
898 // for clarity be explicit on which function is called
900 }
901
902
903
904 template <typename Number>
906 const IndexSet &locally_stored_indices)
907 {
908 // virtual functions called in constructors and destructors never use the
909 // override in a derived class
910 // for clarity be explicit on which function is called
911 ReadWriteVector<Number>::reinit(locally_stored_indices);
912 }
913
914
915
916 template <typename Number>
919 {
920 return stored_elements.size();
921 }
922
923
924
925 template <typename Number>
928 {
930 }
931
932
933
934 template <typename Number>
935 inline const IndexSet &
937 {
938 return stored_elements;
939 }
940
941
942
943 template <typename Number>
946 {
947 return values.begin();
948 }
949
950
951
952 template <typename Number>
955 {
956 return values.begin();
957 }
958
959
960
961 template <typename Number>
964 {
965 return values.end();
966 }
967
968
969
970 template <typename Number>
973 {
974 return values.end();
975 }
976
977
978
979 template <typename Number>
980 inline Number
981 ReadWriteVector<Number>::operator()(const size_type global_index) const
982 {
983 return values[global_to_local(global_index)];
984 }
985
986
987
988 template <typename Number>
989 inline Number &
991 {
992 return values[global_to_local(global_index)];
993 }
994
995
996
997 template <typename Number>
998 inline Number
999 ReadWriteVector<Number>::operator[](const size_type global_index) const
1000 {
1001 return operator()(global_index);
1002 }
1003
1004
1005
1006 template <typename Number>
1007 inline Number &
1009 {
1010 return operator()(global_index);
1011 }
1012
1013
1014
1015 template <typename Number>
1016 template <typename Number2>
1017 inline void
1019 const std::vector<size_type> &indices,
1020 std::vector<Number2> &extracted_values) const
1021 {
1022 for (size_type i = 0; i < indices.size(); ++i)
1023 extracted_values[i] = operator()(indices[i]);
1024 }
1025
1026
1027
1028 template <typename Number>
1029 void
1032 ArrayView<Number> &entries) const
1033 {
1034 AssertDimension(indices.size(), entries.size());
1035 for (unsigned int i = 0; i < indices.size(); ++i)
1036 entries[i] = (*this)[indices[i]];
1037 }
1038
1039
1040
1041 template <typename Number>
1042 template <typename ForwardIterator, typename OutputIterator>
1043 inline void
1045 ForwardIterator indices_begin,
1046 const ForwardIterator indices_end,
1047 OutputIterator values_begin) const
1048 {
1049 while (indices_begin != indices_end)
1050 {
1051 *values_begin = operator()(*indices_begin);
1052 ++indices_begin;
1053 ++values_begin;
1054 }
1055 }
1056
1057
1058
1059 template <typename Number>
1060 inline Number
1061 ReadWriteVector<Number>::local_element(const size_type local_index) const
1062 {
1063 AssertIndexRange(local_index, this->locally_owned_size());
1064
1065 return values[local_index];
1066 }
1067
1068
1069
1070 template <typename Number>
1071 inline Number &
1073 {
1074 AssertIndexRange(local_index, this->locally_owned_size());
1075
1076 return values[local_index];
1077 }
1078
1079
1080
1081 template <typename Number>
1082 template <typename Number2>
1083 inline void
1084 ReadWriteVector<Number>::add(const std::vector<size_type> &indices,
1085 const std::vector<Number2> &values)
1086 {
1087 AssertDimension(indices.size(), values.size());
1088 add(indices.size(), indices.data(), values.data());
1089 }
1090
1091
1092
1093 template <typename Number>
1094 template <typename Number2>
1095 inline void
1096 ReadWriteVector<Number>::add(const std::vector<size_type> &indices,
1097 const ReadWriteVector<Number2> &values)
1098 {
1099 const size_type size = indices.size();
1100 for (size_type i = 0; i < size; ++i)
1101 {
1102 Assert(
1104 ExcMessage(
1105 "The given value is not finite but either infinite or Not A Number (NaN)"));
1106 this->operator()(indices[i]) += values[indices[i]];
1107 }
1108 }
1109
1110
1111
1112 template <typename Number>
1113 template <typename Number2>
1114 inline void
1116 const size_type *indices,
1117 const Number2 *values_to_add)
1118 {
1119 for (size_type i = 0; i < n_indices; ++i)
1120 {
1121 Assert(
1123 ExcMessage(
1124 "The given value is not finite but either infinite or Not A Number (NaN)"));
1125 this->operator()(indices[i]) += values_to_add[i];
1126 }
1127 }
1128
1129
1130
1131 template <typename Number>
1132 inline unsigned int
1134 const types::global_dof_index global_index) const
1135 {
1136 // the following will throw an exception if the global_index is not
1137 // in the remaining_elements
1138 return static_cast<unsigned int>(
1139 stored_elements.index_within_set(global_index));
1140 }
1141
1142
1143
1144 template <typename Number>
1145 template <typename Functor>
1147 ReadWriteVector<Number> &parent,
1148 const Functor &functor)
1149 : parent(parent)
1150 , functor(functor)
1151 {}
1152
1153
1154
1155 template <typename Number>
1156 template <typename Functor>
1157 void
1159 const size_type begin,
1160 const size_type end)
1161 {
1162 for (size_type i = begin; i < end; ++i)
1163 functor(parent.values[i]);
1164 }
1165
1166#endif // ifndef DOXYGEN
1167
1168} // end of namespace LinearAlgebra
1169
1170
1171
1179template <typename Number>
1180inline void
1186
1187
1189
1190#endif
iterator end()
iterator begin()
pointer data()
size_type size() const
std::size_t size() const
Definition array_view.h:684
size_type size() const
Definition index_set.h:1765
size_type index_within_set(const size_type global_index) const
Definition index_set.h:1990
size_type n_elements() const
Definition index_set.h:1923
FunctorTemplate(ReadWriteVector< Number > &parent, const Functor &functor)
virtual void operator()(const size_type begin, const size_type end)
void import_elements(const ::Vector< Number > &vec, VectorOperation::values operation, const std::shared_ptr< const Utilities::MPI::CommunicationPatternBase > &communication_pattern={})
void reinit(const ReadWriteVector< Number2 > &in_vector, const bool omit_zeroing_entries=false)
Number & operator()(const size_type global_index)
ReadWriteVector(const ReadWriteVector< Number > &in_vector)
ReadWriteVector(const size_type size)
std::enable_if_t< std::is_same_v< Dummy, Number > &&::is_tpetra_type< Number >::value > import_elements(const TpetraWrappers::Vector< Number > &tpetra_vec, VectorOperation::values operation, const std::shared_ptr< const Utilities::MPI::CommunicationPatternBase > &communication_pattern={})
void import_elements(const TrilinosWrappers::MPI::Vector &trilinos_vec, VectorOperation::values operation, const std::shared_ptr< const Utilities::MPI::CommunicationPatternBase > &communication_pattern={})
size_type size() const override
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
std::enable_if_t< std::is_same_v< Dummy, Number > &&::is_tpetra_type< Number >::value > import_elements(const Tpetra::Vector< Number, int, types::signed_global_dof_index > &tpetra_vector, const IndexSet &locally_owned_elements, VectorOperation::values operation, const MPI_Comm mpi_comm, const std::shared_ptr< const Utilities::MPI::CommunicationPatternBase > &communication_pattern)
void import_elements(const CUDAWrappers::Vector< Number > &cuda_vec, VectorOperation::values operation, const std::shared_ptr< const Utilities::MPI::CommunicationPatternBase > &communication_pattern={})
Number & operator[](const size_type global_index)
void print(std::ostream &out, const unsigned int precision=3, const bool scientific=true) const
void add(const std::vector< size_type > &indices, const ReadWriteVector< Number2 > &values)
Number & local_element(const size_type local_index)
void import_elements(const EpetraWrappers::Vector &epetra_vec, VectorOperation::values operation, const std::shared_ptr< const Utilities::MPI::CommunicationPatternBase > &communication_pattern={})
AlignedVector< Number > values
unsigned int global_to_local(const types::global_dof_index global_index) const
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)
virtual void extract_subvector_to(const ArrayView< const types::global_dof_index > &indices, ArrayView< Number > &entries) const override
void import_elements(const PETScWrappers::MPI::Vector &petsc_vec, VectorOperation::values operation, const std::shared_ptr< const Utilities::MPI::CommunicationPatternBase > &communication_pattern={})
ReadWriteVector< Number > & operator=(const Number s)
TpetraWrappers::CommunicationPattern create_tpetra_comm_pattern(const IndexSet &source_index_set, const MPI_Comm mpi_comm)
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)
size_type locally_owned_size() const
const_iterator end() const
ReadWriteVector< Number > & operator=(const ReadWriteVector< Number > &in_vector)
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)
EpetraWrappers::CommunicationPattern create_epetra_comm_pattern(const IndexSet &source_index_set, const MPI_Comm mpi_comm)
types::global_dof_index size_type
const IndexSet & get_stored_elements() const
ReadWriteVector< Number > & operator=(const ReadWriteVector< Number2 > &in_vector)
void import_elements(const distributed::Vector< Number, MemorySpace > &vec, VectorOperation::values operation, const std::shared_ptr< const Utilities::MPI::CommunicationPatternBase > &communication_pattern={})
void swap(ReadWriteVector< Number > &v)
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
void import_elements(const Epetra_MultiVector &multivector, const IndexSet &locally_owned_elements, VectorOperation::values operation, const MPI_Comm mpi_comm, const std::shared_ptr< const Utilities::MPI::CommunicationPatternBase > &communication_pattern)
void swap(LinearAlgebra::ReadWriteVector< Number > &u, LinearAlgebra::ReadWriteVector< Number > &v)
#define DEAL_II_DEPRECATED
Definition config.h:207
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:502
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:503
#define Assert(cond, exc)
#define AssertDimension(dim1, dim2)
#define AssertIndexRange(index, range)
static ::ExceptionBase & ExcMessage(std::string arg1)
bool is_finite(const double x)
Definition numbers.h:538
unsigned int global_dof_index
Definition types.h:81