Reference documentation for deal.II version Git 3f1f337db3 2021-10-23 13:19:02 -0600
\(\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\}}\)
vector.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1999 - 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_vector_h
17 #define dealii_vector_h
18 
19 
20 #include <deal.II/base/config.h>
21 
24 #include <deal.II/base/index_set.h>
26 
28 
31 
32 #include <boost/serialization/split_member.hpp>
33 
34 #include <algorithm>
35 #include <initializer_list>
36 #include <iosfwd>
37 #include <iterator>
38 #include <vector>
39 
41 
42 
43 // Forward declarations
44 #ifndef DOXYGEN
45 # ifdef DEAL_II_WITH_PETSC
46 namespace PETScWrappers
47 {
48  class VectorBase;
49 }
50 # endif
51 
52 # ifdef DEAL_II_WITH_TRILINOS
53 namespace TrilinosWrappers
54 {
55  namespace MPI
56  {
57  class Vector;
58  }
59 } // namespace TrilinosWrappers
60 # endif
61 
62 template <typename number>
63 class LAPACKFullMatrix;
64 
65 template <typename>
66 class BlockVector;
67 
68 namespace parallel
69 {
70  namespace internal
71  {
72  class TBBPartitioner;
73  }
74 } // namespace parallel
75 #endif
76 
77 
108 template <typename Number>
109 class Vector : public Subscriptor
110 {
111 public:
112  // The assertion in vector.templates.h for whether or not a number is
113  // finite is not compatible for AD number types.
114  static_assert(
116  "The Vector class does not support auto-differentiable numbers.");
117 
122  using value_type = Number;
123  using pointer = value_type *;
124  using const_pointer = const value_type *;
125  using iterator = value_type *;
126  using const_iterator = const value_type *;
128  using const_reference = const value_type &;
130 
141 
149  Vector();
150 
160  Vector(const Vector<Number> &v);
161 
166  Vector(Vector<Number> &&v) noexcept = default;
167 
176  template <typename OtherNumber>
177  explicit Vector(const Vector<OtherNumber> &v);
178 
194  template <typename OtherNumber>
195  explicit Vector(const std::initializer_list<OtherNumber> &v);
196 
197 #ifdef DEAL_II_WITH_PETSC
198 
209  explicit Vector(const PETScWrappers::VectorBase &v);
210 #endif
211 
212 #ifdef DEAL_II_WITH_TRILINOS
213 
227  explicit Vector(const TrilinosWrappers::MPI::Vector &v);
228 #endif
229 
239  explicit Vector(const size_type n);
240 
245  template <typename InputIterator>
246  Vector(const InputIterator first, const InputIterator last);
247 
252  virtual ~Vector() override = default;
253 
267  void
268  compress(::VectorOperation::values operation =
269  ::VectorOperation::unknown) const;
270 
286  virtual void
287  reinit(const size_type N, const bool omit_zeroing_entries = false);
288 
312  void
313  grow_or_shrink(const size_type N);
314 
322  void
323  apply_givens_rotation(const std::array<Number, 3> &csr,
324  const size_type i,
325  const size_type k);
326 
334  template <typename Number2>
335  void
336  reinit(const Vector<Number2> &V, const bool omit_zeroing_entries = false);
337 
353  virtual void
354  swap(Vector<Number> &v);
355 
366  Vector<Number> &
367  operator=(const Number s);
368 
374  Vector<Number> &
375  operator=(const Vector<Number> &v);
376 
382  Vector<Number> &
383  operator=(Vector<Number> &&v) noexcept = default;
384 
390  template <typename Number2>
391  Vector<Number> &
392  operator=(const Vector<Number2> &v);
393 
397  Vector<Number> &
398  operator=(const BlockVector<Number> &v);
399 
400 #ifdef DEAL_II_WITH_PETSC
401 
412  Vector<Number> &
413  operator=(const PETScWrappers::VectorBase &v);
414 #endif
415 
416 
417 #ifdef DEAL_II_WITH_TRILINOS
418 
433  Vector<Number> &
434  operator=(const TrilinosWrappers::MPI::Vector &v);
435 #endif
436 
442  template <typename Number2>
443  bool
444  operator==(const Vector<Number2> &v) const;
445 
451  template <typename Number2>
452  bool
453  operator!=(const Vector<Number2> &v) const;
454 
456 
457 
462 
476  template <typename Number2>
477  Number
478  operator*(const Vector<Number2> &V) const;
479 
487  real_type
488  norm_sqr() const;
489 
497  Number
498  mean_value() const;
499 
507  real_type
508  l1_norm() const;
509 
518  real_type
519  l2_norm() const;
520 
529  real_type
530  lp_norm(const real_type p) const;
531 
535  real_type
536  linfty_norm() const;
537 
562  Number
563  add_and_dot(const Number a, const Vector<Number> &V, const Vector<Number> &W);
564 
566 
567 
572 
576  pointer
577  data();
578 
583  data() const;
584 
590  iterator
591  begin();
592 
597  begin() const;
598 
602  iterator
603  end();
604 
610  end() const;
611 
615  Number
616  operator()(const size_type i) const;
617 
621  Number &
622  operator()(const size_type i);
623 
629  Number
630  operator[](const size_type i) const;
631 
637  Number &
638  operator[](const size_type i);
639 
655  template <typename OtherNumber>
656  void
657  extract_subvector_to(const std::vector<size_type> &indices,
658  std::vector<OtherNumber> & values) const;
659 
687  template <typename ForwardIterator, typename OutputIterator>
688  void
689  extract_subvector_to(ForwardIterator indices_begin,
690  const ForwardIterator indices_end,
691  OutputIterator values_begin) const;
693 
694 
699 
705  Vector<Number> &
706  operator+=(const Vector<Number> &V);
707 
713  Vector<Number> &
714  operator-=(const Vector<Number> &V);
715 
720  template <typename OtherNumber>
721  void
722  add(const std::vector<size_type> & indices,
723  const std::vector<OtherNumber> &values);
724 
729  template <typename OtherNumber>
730  void
731  add(const std::vector<size_type> &indices, const Vector<OtherNumber> &values);
732 
738  template <typename OtherNumber>
739  void
740  add(const size_type n_elements,
741  const size_type * indices,
742  const OtherNumber *values);
743 
750  void
751  add(const Number s);
752 
758  void
759  add(const Number a,
760  const Vector<Number> &V,
761  const Number b,
762  const Vector<Number> &W);
763 
769  void
770  add(const Number a, const Vector<Number> &V);
771 
777  void
778  sadd(const Number s, const Vector<Number> &V);
779 
785  void
786  sadd(const Number s, const Number a, const Vector<Number> &V);
787 
793  Vector<Number> &
794  operator*=(const Number factor);
795 
801  Vector<Number> &
802  operator/=(const Number factor);
803 
811  void
812  scale(const Vector<Number> &scaling_factors);
813 
819  template <typename Number2>
820  void
821  scale(const Vector<Number2> &scaling_factors);
822 
828  void
829  equ(const Number a, const Vector<Number> &u);
830 
834  template <typename Number2>
835  void
836  equ(const Number a, const Vector<Number2> &u);
837 
842  void
843  update_ghost_values() const;
845 
846 
857  void
858  print(std::ostream & out,
859  const unsigned int precision = 3,
860  const bool scientific = true,
861  const bool across = true) const;
862 
868  void
869  block_write(std::ostream &out) const;
870 
882  void
883  block_read(std::istream &in);
884 
890  template <class Archive>
891  void
892  save(Archive &ar, const unsigned int version) const;
893 
899  template <class Archive>
900  void
901  load(Archive &ar, const unsigned int version);
902 
903 #ifdef DOXYGEN
904 
909  template <class Archive>
910  void
911  serialize(Archive &archive, const unsigned int version);
912 #else
913  // This macro defines the serialize() method that is compatible with
914  // the templated save() and load() method that have been implemented.
915  BOOST_SERIALIZATION_SPLIT_MEMBER()
916 #endif
917 
926 
932  bool
933  in_local_range(const size_type global_index) const;
934 
950  IndexSet
951  locally_owned_elements() const;
952 
956  size_type
957  size() const;
958 
966  size_type
967  locally_owned_size() const;
968 
974  bool
975  all_zero() const;
976 
986  bool
987  is_non_negative() const;
988 
993  std::size_t
994  memory_consumption() const;
995 
1001  bool
1002  has_ghost_elements() const;
1004 
1005 private:
1010 
1016  void
1017  maybe_reset_thread_partitioner();
1018 
1022  void
1023  do_reinit(const size_type new_size,
1024  const bool omit_zeroing_entries,
1025  const bool reset_partitioner);
1026 
1031  mutable std::shared_ptr<parallel::internal::TBBPartitioner>
1033 
1034  // Make all other vector types friends.
1035  template <typename Number2>
1036  friend class Vector;
1037 };
1038 
1040 /*----------------------- Inline functions ----------------------------------*/
1041 
1042 
1043 #ifndef DOXYGEN
1044 
1045 
1046 //------------------------ declarations for explicit specializations
1047 template <>
1049 Vector<int>::lp_norm(const real_type) const;
1050 
1051 
1052 //------------------------ inline functions
1053 
1054 template <typename Number>
1055 inline Vector<Number>::Vector()
1056 {
1057  // virtual functions called in constructors and destructors never use the
1058  // override in a derived class
1059  // for clarity be explicit on which function is called
1061 }
1062 
1063 
1064 
1065 template <typename Number>
1066 template <typename OtherNumber>
1067 Vector<Number>::Vector(const std::initializer_list<OtherNumber> &v)
1068  : Vector(v.begin(), v.end())
1069 {}
1070 
1071 
1072 
1073 template <typename Number>
1074 template <typename InputIterator>
1075 Vector<Number>::Vector(const InputIterator first, const InputIterator last)
1076 {
1077  // allocate memory. do not initialize it, as we will copy over to it in a
1078  // second
1079  reinit(std::distance(first, last), true);
1080  std::copy(first, last, begin());
1081 }
1082 
1083 
1084 
1085 template <typename Number>
1086 inline Vector<Number>::Vector(const size_type n)
1087 {
1088  // virtual functions called in constructors and destructors never use the
1089  // override in a derived class
1090  // for clarity be explicit on which function is called
1091  Vector<Number>::reinit(n, false);
1092 }
1093 
1094 
1095 
1096 template <typename Number>
1097 inline typename Vector<Number>::size_type
1098 Vector<Number>::size() const
1099 {
1100  return values.size();
1101 }
1102 
1103 
1104 
1105 template <typename Number>
1106 inline typename Vector<Number>::size_type
1108 {
1109  return values.size();
1110 }
1111 
1112 
1113 
1114 template <typename Number>
1115 inline bool
1117 {
1118  return true;
1119 }
1120 
1121 
1122 
1123 template <typename Number>
1124 inline typename Vector<Number>::pointer
1126 {
1127  return values.data();
1128 }
1129 
1130 
1131 
1132 template <typename Number>
1133 inline typename Vector<Number>::const_pointer
1134 Vector<Number>::data() const
1135 {
1136  return values.data();
1137 }
1138 
1139 
1140 
1141 template <typename Number>
1142 inline typename Vector<Number>::iterator
1144 {
1145  return values.begin();
1146 }
1147 
1148 
1149 
1150 template <typename Number>
1151 inline typename Vector<Number>::const_iterator
1152 Vector<Number>::begin() const
1153 {
1154  return values.begin();
1155 }
1156 
1157 
1158 
1159 template <typename Number>
1160 inline typename Vector<Number>::iterator
1162 {
1163  return values.end();
1164 }
1165 
1166 
1167 
1168 template <typename Number>
1169 inline typename Vector<Number>::const_iterator
1170 Vector<Number>::end() const
1171 {
1172  return values.end();
1173 }
1174 
1175 
1176 
1177 template <typename Number>
1178 inline Number
1180 {
1181  AssertIndexRange(i, size());
1182  return values[i];
1183 }
1184 
1185 
1186 
1187 template <typename Number>
1188 inline Number &
1190 {
1191  AssertIndexRange(i, size());
1192  return values[i];
1193 }
1194 
1195 
1196 
1197 template <typename Number>
1198 inline Number
1200 {
1201  return operator()(i);
1202 }
1203 
1204 
1205 
1206 template <typename Number>
1207 inline Number &
1209 {
1210  return operator()(i);
1211 }
1212 
1213 
1214 
1215 template <typename Number>
1216 template <typename OtherNumber>
1217 inline void
1218 Vector<Number>::extract_subvector_to(const std::vector<size_type> &indices,
1219  std::vector<OtherNumber> & values) const
1220 {
1221  for (size_type i = 0; i < indices.size(); ++i)
1222  values[i] = operator()(indices[i]);
1223 }
1224 
1225 
1226 
1227 template <typename Number>
1228 template <typename ForwardIterator, typename OutputIterator>
1229 inline void
1230 Vector<Number>::extract_subvector_to(ForwardIterator indices_begin,
1231  const ForwardIterator indices_end,
1232  OutputIterator values_begin) const
1233 {
1234  while (indices_begin != indices_end)
1235  {
1236  *values_begin = operator()(*indices_begin);
1237  indices_begin++;
1238  values_begin++;
1239  }
1240 }
1241 
1242 
1243 
1244 template <typename Number>
1245 inline Vector<Number> &
1246 Vector<Number>::operator/=(const Number factor)
1247 {
1248  AssertIsFinite(factor);
1249  Assert(factor != Number(0.), ExcZero());
1250 
1251  this->operator*=(Number(1.) / factor);
1252  return *this;
1253 }
1254 
1255 
1256 
1257 template <typename Number>
1258 template <typename OtherNumber>
1259 inline void
1260 Vector<Number>::add(const std::vector<size_type> & indices,
1261  const std::vector<OtherNumber> &values)
1262 {
1263  Assert(indices.size() == values.size(),
1264  ExcDimensionMismatch(indices.size(), values.size()));
1265  add(indices.size(), indices.data(), values.data());
1266 }
1267 
1268 
1269 
1270 template <typename Number>
1271 template <typename OtherNumber>
1272 inline void
1273 Vector<Number>::add(const std::vector<size_type> &indices,
1274  const Vector<OtherNumber> & values)
1275 {
1276  Assert(indices.size() == values.size(),
1277  ExcDimensionMismatch(indices.size(), values.size()));
1278  add(indices.size(), indices.data(), values.values.begin());
1279 }
1280 
1281 
1282 
1283 template <typename Number>
1284 template <typename OtherNumber>
1285 inline void
1286 Vector<Number>::add(const size_type n_indices,
1287  const size_type * indices,
1288  const OtherNumber *values)
1289 {
1290  for (size_type i = 0; i < n_indices; ++i)
1291  {
1292  AssertIndexRange(indices[i], size());
1293  Assert(
1294  numbers::is_finite(values[i]),
1295  ExcMessage(
1296  "The given value is not finite but either infinite or Not A Number (NaN)"));
1297 
1298  this->values[indices[i]] += values[i];
1299  }
1300 }
1301 
1302 
1303 
1304 template <typename Number>
1305 template <typename Number2>
1306 inline bool
1307 Vector<Number>::operator!=(const Vector<Number2> &v) const
1308 {
1309  return !(*this == v);
1310 }
1311 
1312 
1313 
1314 template <typename Number>
1316 {}
1317 
1318 
1319 template <typename Number>
1320 inline bool
1322 {
1323  return false;
1324 }
1325 
1326 template <typename Number>
1327 inline void
1329 {}
1330 
1331 
1332 
1333 // Moved from vector.templates.h as an inline function by Luca Heltai
1334 // on 2009/04/12 to prevent strange compiling errors, after making
1335 // swap virtual.
1336 template <typename Number>
1337 inline void
1338 Vector<Number>::swap(Vector<Number> &v)
1339 {
1340  values.swap(v.values);
1341  std::swap(thread_loop_partitioner, v.thread_loop_partitioner);
1342 }
1343 
1344 
1345 
1346 template <typename Number>
1347 template <class Archive>
1348 inline void
1349 Vector<Number>::save(Archive &ar, const unsigned int) const
1350 {
1351  // forward to serialization function in the base class.
1352  ar &static_cast<const Subscriptor &>(*this);
1353  ar &values;
1354 }
1355 
1356 
1357 
1358 template <typename Number>
1359 template <class Archive>
1360 inline void
1361 Vector<Number>::load(Archive &ar, const unsigned int)
1362 {
1363  // the load stuff again from the archive
1364  ar &static_cast<Subscriptor &>(*this);
1365  ar &values;
1367 }
1368 
1369 #endif
1370 
1371 
1384 template <typename Number>
1385 inline void
1386 swap(Vector<Number> &u, Vector<Number> &v)
1387 {
1388  u.swap(v);
1389 }
1390 
1391 
1399 template <typename number>
1400 inline std::ostream &
1401 operator<<(std::ostream &out, const Vector<number> &v)
1402 {
1403  Assert(v.size() != 0, ExcEmptyObject());
1404  AssertThrow(out, ExcIO());
1405 
1406  for (typename Vector<number>::size_type i = 0; i < v.size() - 1; ++i)
1407  out << v(i) << ' ';
1408  out << v(v.size() - 1);
1409 
1410  AssertThrow(out, ExcIO());
1411 
1412  return out;
1413 }
1414 
1423 template <typename Number>
1424 struct is_serial_vector<Vector<Number>> : std::true_type
1425 {};
1426 
1427 
1429 
1430 #endif
void update_ghost_values() const
void reinit(MatrixBlock< MatrixType > &v, const BlockSparsityPattern &p)
Definition: matrix_block.h:618
Number operator[](const size_type i) const
bool operator!=(const Vector< Number2 > &v) const
void load(Archive &ar, const unsigned int version)
void extract_subvector_to(const std::vector< size_type > &indices, std::vector< OtherNumber > &values) const
Number operator()(const size_type i) const
static ::ExceptionBase & ExcIO()
typename numbers::NumberTraits< typename VectorType::value_type >::real_type real_type
Definition: vector.h:140
bool operator!=(const AlignedVector< T > &lhs, const AlignedVector< T > &rhs)
std::vector< value_type > l2_norm(const typename ::Triangulation< dim, spacedim >::cell_iterator &parent, const value_type parent_value)
void scale(const double scaling_factor, Triangulation< dim, spacedim > &triangulation)
Definition: grid_tools.cc:2047
#define AssertIndexRange(index, range)
Definition: exceptions.h:1720
pointer data()
void swap(VectorBase &u, VectorBase &v)
__global__ void add_and_dot(Number *res, Number *v1, const Number *v2, const Number *v3, const Number a, const size_type N)
static const char V
#define AssertThrow(cond, exc)
Definition: exceptions.h:1571
pointer data()
iterator end()
bool operator==(const AlignedVector< T > &lhs, const AlignedVector< T > &rhs)
bool is_finite(const double x)
Definition: numbers.h:537
void add(const std::vector< size_type > &indices, const std::vector< OtherNumber > &values)
static ::ExceptionBase & ExcMessage(std::string arg1)
bool has_ghost_elements() const
std::string compress(const std::string &input)
Definition: utilities.cc:392
typename VectorType::value_type value_type
Definition: vector.h:122
bool in_local_range(const size_type global_index) const
void maybe_reset_thread_partitioner()
#define Assert(cond, exc)
Definition: exceptions.h:1461
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
real_type lp_norm(const real_type p) const
Number linfty_norm(const Tensor< 2, dim, Number > &t)
Definition: tensor.h:3025
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:401
VectorType::value_type * end(VectorType &V)
Vector< Number > & operator*=(const Number factor)
iterator begin()
__global__ void equ(Number *val, const Number a, const Number *V_val, const size_type N)
__global__ void sadd(const Number s, Number *val, const Number a, const Number *V_val, const size_type N)
SymmetricTensor< 2, dim, Number > b(const Tensor< 2, dim, Number > &F)
Point< 2 > first
Definition: grid_out.cc:4586
void swap(Vector< Number > &u, Vector< Number > &v)
Definition: vector.h:1386
unsigned int global_dof_index
Definition: types.h:76
virtual void swap(Vector< Number > &v)
Vector< Number > & operator/=(const Number factor)
void save(Archive &ar, const unsigned int version) const
iterator end()
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:400
VectorType::value_type * begin(VectorType &V)
iterator begin()
size_type locally_owned_size() const
Number l1_norm(const Tensor< 2, dim, Number > &t)
Definition: tensor.h:2999
static const char N
static ::ExceptionBase & ExcEmptyObject()
void compress(::VectorOperation::values operation=::VectorOperation::unknown) const
size_type size() const
virtual void reinit(const size_type N, const bool omit_zeroing_entries=false)
AlignedVector< Number > values
Definition: vector.h:1009
size_type size() const
std::shared_ptr< parallel::internal::TBBPartitioner > thread_loop_partitioner
Definition: vector.h:1032
TrilinosWrappers::types::int_type global_index(const Epetra_BlockMap &map, const ::types::global_dof_index i)
std::enable_if< std::is_floating_point< T >::value &&std::is_floating_point< U >::value, typename ProductType< std::complex< T >, std::complex< U > >::type >::type operator*(const std::complex< T > &left, const std::complex< U > &right)
types::global_dof_index size_type
Definition: vector.h:129
static ::ExceptionBase & ExcZero()
void copy(const T *begin, const T *end, U *dest)
#define AssertIsFinite(number)
Definition: exceptions.h:1746
std::enable_if< std::is_fundamental< T >::value, std::size_t >::type memory_consumption(const T &t)