Reference documentation for deal.II version 9.4.0
\(\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
24#include <deal.II/base/mpi.h>
28#include <deal.II/base/types.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
50namespace 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
62namespace PETScWrappers
63{
64 namespace MPI
65 {
66 class Vector;
67 }
68} // namespace PETScWrappers
69# endif
70
71# ifdef DEAL_II_WITH_TRILINOS
72namespace TrilinosWrappers
73{
74 namespace MPI
75 {
76 class Vector;
77 }
78} // namespace TrilinosWrappers
79# endif
80
81# ifdef DEAL_II_WITH_CUDA
82namespace LinearAlgebra
83{
84 namespace CUDAWrappers
85 {
86 template <typename>
87 class Vector;
88 }
89} // namespace LinearAlgebra
90# endif
91#endif
92
93namespace LinearAlgebra
94{
132 template <typename Number>
134 {
135 public:
141 using value_type = Number;
143 using const_pointer = const value_type *;
145 using const_iterator = const value_type *;
150
159
164
170
175 explicit ReadWriteVector(const IndexSet &locally_stored_indices);
176
180 ~ReadWriteVector() override = default;
181
190 virtual void
191 reinit(const size_type size, const bool omit_zeroing_entries = false);
192
201 template <typename Number2>
202 void
204 const bool omit_zeroing_entries = false);
205
215 virtual void
216 reinit(const IndexSet &locally_stored_indices,
217 const bool omit_zeroing_entries = false);
218
219
220#ifdef DEAL_II_WITH_TRILINOS
232 void
234#endif
235
249 template <typename Functor>
250 void
251 apply(const Functor &func);
252
265 void
267
273
277 template <typename Number2>
280
286 operator=(const Number s);
287
297 void
298 import(const ::Vector<Number> &vec,
299 VectorOperation::values operation,
300 const std::shared_ptr<const Utilities::MPI::CommunicationPatternBase>
301 &communication_pattern = {});
302
312 void
314 VectorOperation::values operation,
315 const std::shared_ptr<const Utilities::MPI::CommunicationPatternBase>
316 &communication_pattern = {});
317
326 template <typename MemorySpace>
327 void
329 VectorOperation::values operation,
330 const std::shared_ptr<const Utilities::MPI::CommunicationPatternBase>
331 &communication_pattern = {});
332
333#ifdef DEAL_II_WITH_PETSC
342 void
343 import(const PETScWrappers::MPI::Vector &petsc_vec,
344 VectorOperation::values operation,
345 const std::shared_ptr<const Utilities::MPI::CommunicationPatternBase>
346 &communication_pattern = {});
347#endif
348
349#ifdef DEAL_II_WITH_TRILINOS
360 void
361 import(const TrilinosWrappers::MPI::Vector &trilinos_vec,
362 VectorOperation::values operation,
363 const std::shared_ptr<const Utilities::MPI::CommunicationPatternBase>
364 &communication_pattern = {});
365
366# ifdef DEAL_II_TRILINOS_WITH_TPETRA
375 void
376 import(const TpetraWrappers::Vector<Number> &tpetra_vec,
377 VectorOperation::values operation,
378 const std::shared_ptr<const Utilities::MPI::CommunicationPatternBase>
379 &communication_pattern = {});
380# endif
381
390 void
391 import(const EpetraWrappers::Vector &epetra_vec,
392 VectorOperation::values operation,
393 const std::shared_ptr<const Utilities::MPI::CommunicationPatternBase>
394 &communication_pattern = {});
395#endif
396
397#ifdef DEAL_II_WITH_CUDA
404 void
405 import(const CUDAWrappers::Vector<Number> &cuda_vec,
406 VectorOperation::values operation,
407 const std::shared_ptr<const Utilities::MPI::CommunicationPatternBase>
408 &communication_pattern = {});
409#endif
410
420 size() const;
421
431 n_elements() const;
432
439
443 const IndexSet &
445
453
459 begin() const;
460
467
473 end() const;
475
476
481
487 Number
489
495 Number &
497
505 Number
507
515 Number &
517
533 template <typename Number2>
534 void
535 extract_subvector_to(const std::vector<size_type> &indices,
536 std::vector<Number2> & values) const;
537
565 template <typename ForwardIterator, typename OutputIterator>
566 void
567 extract_subvector_to(ForwardIterator indices_begin,
568 const ForwardIterator indices_end,
569 OutputIterator values_begin) const;
570
581 Number
582 local_element(const size_type local_index) const;
583
594 Number &
595 local_element(const size_type local_index);
597
598
603
608 template <typename Number2>
609 void
610 add(const std::vector<size_type> &indices,
611 const std::vector<Number2> & values);
612
617 template <typename Number2>
618 void
619 add(const std::vector<size_type> & indices,
621
627 template <typename Number2>
628 void
630 const size_type *indices,
631 const Number2 * values);
632
636 void
637 print(std::ostream & out,
638 const unsigned int precision = 3,
639 const bool scientific = true) const;
640
644 std::size_t
647
648 protected:
649#ifdef DEAL_II_WITH_TRILINOS
650# ifdef DEAL_II_TRILINOS_WITH_TPETRA
656 void
657 import(
658 const Tpetra::Vector<Number, int, types::global_dof_index> &tpetra_vector,
659 const IndexSet & locally_owned_elements,
660 VectorOperation::values operation,
661 const MPI_Comm & mpi_comm,
662 const std::shared_ptr<const Utilities::MPI::CommunicationPatternBase>
663 &communication_pattern);
664# endif
665
671 void
672 import(const Epetra_MultiVector &multivector,
673 const IndexSet & locally_owned_elements,
674 VectorOperation::values operation,
675 const MPI_Comm & mpi_comm,
676 const std::shared_ptr<const Utilities::MPI::CommunicationPatternBase>
677 &communication_pattern);
678#endif
679
683 unsigned int
685 {
686 // the following will throw an exception if the global_index is not
687 // in the remaining_elements
688 return static_cast<unsigned int>(
690 }
691
695 void
696 resize_val(const size_type new_allocated_size);
697
698#ifdef DEAL_II_WITH_TRILINOS
699# ifdef DEAL_II_TRILINOS_WITH_TPETRA
705 create_tpetra_comm_pattern(const IndexSet &source_index_set,
706 const MPI_Comm &mpi_comm);
707# endif
708
714 create_epetra_comm_pattern(const IndexSet &source_index_set,
715 const MPI_Comm &mpi_comm);
716#endif
717
722
727
732 std::shared_ptr<Utilities::MPI::CommunicationPatternBase> comm_pattern;
733
737 std::unique_ptr<Number[], decltype(std::free) *> values;
738
743 mutable std::shared_ptr<::parallel::internal::TBBPartitioner>
745
746 // Make all other ReadWriteVector types friends.
747 template <typename Number2>
748 friend class ReadWriteVector;
749
750 private:
756 template <typename Functor>
758 {
759 public:
764
768 virtual void
770
771 private:
776
780 const Functor &functor;
781 };
782 };
783
787 /*---------------------------- Inline functions ---------------------------*/
788
789#ifndef DOXYGEN
790
791 template <typename Number>
793 : Subscriptor()
794 , values(nullptr, free)
795 {
796 // virtual functions called in constructors and destructors never use the
797 // override in a derived class
798 // for clarity be explicit on which function is called
800 }
801
802
803
804 template <typename Number>
806 const ReadWriteVector<Number> &v)
807 : Subscriptor()
808 , values(nullptr, free)
809 {
810 this->operator=(v);
811 }
812
813
814
815 template <typename Number>
817 : Subscriptor()
818 , values(nullptr, free)
819 {
820 // virtual functions called in constructors and destructors never use the
821 // override in a derived class
822 // for clarity be explicit on which function is called
824 }
825
826
827
828 template <typename Number>
830 const IndexSet &locally_stored_indices)
831 : Subscriptor()
832 , values(nullptr, free)
833 {
834 // virtual functions called in constructors and destructors never use the
835 // override in a derived class
836 // for clarity be explicit on which function is called
837 ReadWriteVector<Number>::reinit(locally_stored_indices);
838 }
839
840
841
842 template <typename Number>
845 {
846 return stored_elements.size();
847 }
848
849
850
851 template <typename Number>
854 {
856 }
857
858
859
860 template <typename Number>
863 {
865 }
866
867
868
869 template <typename Number>
870 inline const IndexSet &
872 {
873 return stored_elements;
874 }
875
876
877
878 template <typename Number>
881 {
882 return values.get();
883 }
884
885
886
887 template <typename Number>
890 {
891 return values.get();
892 }
893
894
895
896 template <typename Number>
899 {
900 return values.get() + this->n_elements();
901 }
902
903
904
905 template <typename Number>
908 {
909 return values.get() + this->n_elements();
910 }
911
912
913
914 template <typename Number>
915 inline Number
917 {
919 }
920
921
922
923 template <typename Number>
924 inline Number &
926 {
928 }
929
930
931
932 template <typename Number>
933 inline Number
935 {
936 return operator()(global_index);
937 }
938
939
940
941 template <typename Number>
942 inline Number &
944 {
945 return operator()(global_index);
946 }
947
948
949
950 template <typename Number>
951 template <typename Number2>
952 inline void
954 const std::vector<size_type> &indices,
955 std::vector<Number2> & extracted_values) const
956 {
957 for (size_type i = 0; i < indices.size(); ++i)
958 extracted_values[i] = operator()(indices[i]);
959 }
960
961
962
963 template <typename Number>
964 template <typename ForwardIterator, typename OutputIterator>
965 inline void
967 ForwardIterator indices_begin,
968 const ForwardIterator indices_end,
969 OutputIterator values_begin) const
970 {
971 while (indices_begin != indices_end)
972 {
973 *values_begin = operator()(*indices_begin);
974 indices_begin++;
975 values_begin++;
976 }
977 }
978
979
980
981 template <typename Number>
982 inline Number
984 {
985 AssertIndexRange(local_index, this->n_elements());
986
987 return values[local_index];
988 }
989
990
991
992 template <typename Number>
993 inline Number &
995 {
996 AssertIndexRange(local_index, this->n_elements());
997
998 return values[local_index];
999 }
1000
1001
1002
1003 template <typename Number>
1004 template <typename Number2>
1005 inline void
1006 ReadWriteVector<Number>::add(const std::vector<size_type> &indices,
1007 const std::vector<Number2> & values)
1008 {
1009 AssertDimension(indices.size(), values.size());
1010 add(indices.size(), indices.data(), values.data());
1011 }
1012
1013
1014
1015 template <typename Number>
1016 template <typename Number2>
1017 inline void
1018 ReadWriteVector<Number>::add(const std::vector<size_type> & indices,
1019 const ReadWriteVector<Number2> &values)
1020 {
1021 const size_type size = indices.size();
1022 for (size_type i = 0; i < size; ++i)
1023 {
1024 Assert(
1026 ExcMessage(
1027 "The given value is not finite but either infinite or Not A Number (NaN)"));
1028 this->operator()(indices[i]) += values[indices[i]];
1029 }
1030 }
1031
1032
1033
1034 template <typename Number>
1035 template <typename Number2>
1036 inline void
1038 const size_type *indices,
1039 const Number2 * values_to_add)
1040 {
1041 for (size_type i = 0; i < n_indices; ++i)
1042 {
1043 Assert(
1045 ExcMessage(
1046 "The given value is not finite but either infinite or Not A Number (NaN)"));
1047 this->operator()(indices[i]) += values_to_add[i];
1048 }
1049 }
1050
1051
1052
1053 template <typename Number>
1054 template <typename Functor>
1056 ReadWriteVector<Number> &parent,
1057 const Functor & functor)
1058 : parent(parent)
1059 , functor(functor)
1060 {}
1061
1062
1063
1064 template <typename Number>
1065 template <typename Functor>
1066 void
1068 const size_type begin,
1069 const size_type end)
1070 {
1071 for (size_type i = begin; i < end; ++i)
1072 functor(parent.values[i]);
1073 }
1074
1075#endif // ifndef DOXYGEN
1076
1077} // end of namespace LinearAlgebra
1078
1079
1080
1088template <typename Number>
1089inline void
1092{
1093 u.swap(v);
1094}
1095
1096
1098
1099#endif
size_type size() const
Definition: index_set.h:1636
size_type index_within_set(const size_type global_index) const
Definition: index_set.h:1923
size_type n_elements() const
Definition: index_set.h:1834
Definition: vector.h:109
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)
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::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)
Number & operator[](const size_type global_index)
void print(std::ostream &out, const unsigned int precision=3, const bool scientific=true) const
FunctorTemplate(ReadWriteVector< Number > &parent, const Functor &functor)
void add(const std::vector< size_type > &indices, const ReadWriteVector< Number2 > &values)
Number & local_element(const size_type local_index)
virtual void operator()(const size_type begin, const size_type end)
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)
ReadWriteVector< Number > & operator=(const Number s)
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)
types::global_dof_index size_type
const IndexSet & get_stored_elements() const
ReadWriteVector< Number > & operator=(const ReadWriteVector< Number2 > &in_vector)
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
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:539
unsigned int global_dof_index
Definition: types.h:76