Reference documentation for deal.II version Git 662679e6de 2020-10-27 13:02:35 -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 - 2020 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 operator*(const Vector<Number2> &V) const;
478 
486  real_type
487  norm_sqr() const;
488 
496  Number
497  mean_value() const;
498 
506  real_type
507  l1_norm() const;
508 
517  real_type
518  l2_norm() const;
519 
528  real_type
529  lp_norm(const real_type p) const;
530 
534  real_type
535  linfty_norm() const;
536 
561  Number
562  add_and_dot(const Number a, const Vector<Number> &V, const Vector<Number> &W);
563 
565 
566 
571 
575  pointer
576  data();
577 
582  data() const;
583 
589  iterator
590  begin();
591 
596  begin() const;
597 
601  iterator
602  end();
603 
609  end() const;
610 
614  Number
615  operator()(const size_type i) const;
616 
620  Number &
621  operator()(const size_type i);
622 
628  Number operator[](const size_type i) const;
629 
635  Number &operator[](const size_type i);
636 
652  template <typename OtherNumber>
653  void
654  extract_subvector_to(const std::vector<size_type> &indices,
655  std::vector<OtherNumber> & values) const;
656 
684  template <typename ForwardIterator, typename OutputIterator>
685  void
686  extract_subvector_to(ForwardIterator indices_begin,
687  const ForwardIterator indices_end,
688  OutputIterator values_begin) const;
690 
691 
696 
702  Vector<Number> &
703  operator+=(const Vector<Number> &V);
704 
710  Vector<Number> &
711  operator-=(const Vector<Number> &V);
712 
717  template <typename OtherNumber>
718  void
719  add(const std::vector<size_type> & indices,
720  const std::vector<OtherNumber> &values);
721 
726  template <typename OtherNumber>
727  void
728  add(const std::vector<size_type> &indices, const Vector<OtherNumber> &values);
729 
735  template <typename OtherNumber>
736  void
737  add(const size_type n_elements,
738  const size_type * indices,
739  const OtherNumber *values);
740 
747  void
748  add(const Number s);
749 
755  void
756  add(const Number a,
757  const Vector<Number> &V,
758  const Number b,
759  const Vector<Number> &W);
760 
766  void
767  add(const Number a, const Vector<Number> &V);
768 
774  void
775  sadd(const Number s, const Vector<Number> &V);
776 
782  void
783  sadd(const Number s, const Number a, const Vector<Number> &V);
784 
790  Vector<Number> &
791  operator*=(const Number factor);
792 
798  Vector<Number> &
799  operator/=(const Number factor);
800 
808  void
809  scale(const Vector<Number> &scaling_factors);
810 
816  template <typename Number2>
817  void
818  scale(const Vector<Number2> &scaling_factors);
819 
825  void
826  equ(const Number a, const Vector<Number> &u);
827 
831  template <typename Number2>
832  void
833  equ(const Number a, const Vector<Number2> &u);
834 
839  void
840  update_ghost_values() const;
842 
843 
854  void
855  print(std::ostream & out,
856  const unsigned int precision = 3,
857  const bool scientific = true,
858  const bool across = true) const;
859 
865  void
866  block_write(std::ostream &out) const;
867 
879  void
880  block_read(std::istream &in);
881 
886  template <class Archive>
887  void
888  save(Archive &ar, const unsigned int version) const;
889 
894  template <class Archive>
895  void
896  load(Archive &ar, const unsigned int version);
897 
898 #ifdef DOXYGEN
899 
903  template <class Archive>
904  void
905  serialize(Archive &archive, const unsigned int version);
906 #else
907  // This macro defines the serialize() method that is compatible with
908  // the templated save() and load() method that have been implemented.
909  BOOST_SERIALIZATION_SPLIT_MEMBER()
910 #endif
911 
920 
926  bool
927  in_local_range(const size_type global_index) const;
928 
944  IndexSet
945  locally_owned_elements() const;
946 
950  size_type
951  size() const;
952 
958  bool
959  all_zero() const;
960 
970  bool
971  is_non_negative() const;
972 
977  std::size_t
978  memory_consumption() const;
979 
985  bool
986  has_ghost_elements() const;
988 
989 private:
994 
1000  void
1001  maybe_reset_thread_partitioner();
1002 
1006  void
1007  do_reinit(const size_type new_size,
1008  const bool omit_zeroing_entries,
1009  const bool reset_partitioner);
1010 
1015  mutable std::shared_ptr<parallel::internal::TBBPartitioner>
1017 
1018  // Make all other vector types friends.
1019  template <typename Number2>
1020  friend class Vector;
1021 };
1022 
1024 /*----------------------- Inline functions ----------------------------------*/
1025 
1026 
1027 #ifndef DOXYGEN
1028 
1029 
1030 //------------------------ declarations for explicit specializations
1031 template <>
1033 Vector<int>::lp_norm(const real_type) const;
1034 
1035 
1036 //------------------------ inline functions
1037 
1038 template <typename Number>
1039 inline Vector<Number>::Vector()
1040 {
1041  // virtual functions called in constructors and destructors never use the
1042  // override in a derived class
1043  // for clarity be explicit on which function is called
1045 }
1046 
1047 
1048 
1049 template <typename Number>
1050 template <typename OtherNumber>
1051 Vector<Number>::Vector(const std::initializer_list<OtherNumber> &v)
1052  : Vector(v.begin(), v.end())
1053 {}
1054 
1055 
1056 
1057 template <typename Number>
1058 template <typename InputIterator>
1059 Vector<Number>::Vector(const InputIterator first, const InputIterator last)
1060 {
1061  // allocate memory. do not initialize it, as we will copy over to it in a
1062  // second
1063  reinit(std::distance(first, last), true);
1064  std::copy(first, last, begin());
1065 }
1066 
1067 
1068 
1069 template <typename Number>
1070 inline Vector<Number>::Vector(const size_type n)
1071 {
1072  // virtual functions called in constructors and destructors never use the
1073  // override in a derived class
1074  // for clarity be explicit on which function is called
1075  Vector<Number>::reinit(n, false);
1076 }
1077 
1078 
1079 
1080 template <typename Number>
1081 inline typename Vector<Number>::size_type
1082 Vector<Number>::size() const
1083 {
1084  return values.size();
1085 }
1086 
1087 
1088 template <typename Number>
1089 inline bool
1091 {
1092  return true;
1093 }
1094 
1095 
1096 
1097 template <typename Number>
1098 inline typename Vector<Number>::pointer
1100 {
1101  return values.data();
1102 }
1103 
1104 
1105 
1106 template <typename Number>
1107 inline typename Vector<Number>::const_pointer
1108 Vector<Number>::data() const
1109 {
1110  return values.data();
1111 }
1112 
1113 
1114 
1115 template <typename Number>
1116 inline typename Vector<Number>::iterator
1118 {
1119  return values.begin();
1120 }
1121 
1122 
1123 
1124 template <typename Number>
1125 inline typename Vector<Number>::const_iterator
1126 Vector<Number>::begin() const
1127 {
1128  return values.begin();
1129 }
1130 
1131 
1132 
1133 template <typename Number>
1134 inline typename Vector<Number>::iterator
1136 {
1137  return values.end();
1138 }
1139 
1140 
1141 
1142 template <typename Number>
1143 inline typename Vector<Number>::const_iterator
1144 Vector<Number>::end() const
1145 {
1146  return values.end();
1147 }
1148 
1149 
1150 
1151 template <typename Number>
1152 inline Number
1154 {
1155  AssertIndexRange(i, size());
1156  return values[i];
1157 }
1158 
1159 
1160 
1161 template <typename Number>
1162 inline Number &
1164 {
1165  AssertIndexRange(i, size());
1166  return values[i];
1167 }
1168 
1169 
1170 
1171 template <typename Number>
1172 inline Number Vector<Number>::operator[](const size_type i) const
1173 {
1174  return operator()(i);
1175 }
1176 
1177 
1178 
1179 template <typename Number>
1180 inline Number &Vector<Number>::operator[](const size_type i)
1181 {
1182  return operator()(i);
1183 }
1184 
1185 
1186 
1187 template <typename Number>
1188 template <typename OtherNumber>
1189 inline void
1190 Vector<Number>::extract_subvector_to(const std::vector<size_type> &indices,
1191  std::vector<OtherNumber> & values) const
1192 {
1193  for (size_type i = 0; i < indices.size(); ++i)
1194  values[i] = operator()(indices[i]);
1195 }
1196 
1197 
1198 
1199 template <typename Number>
1200 template <typename ForwardIterator, typename OutputIterator>
1201 inline void
1202 Vector<Number>::extract_subvector_to(ForwardIterator indices_begin,
1203  const ForwardIterator indices_end,
1204  OutputIterator values_begin) const
1205 {
1206  while (indices_begin != indices_end)
1207  {
1208  *values_begin = operator()(*indices_begin);
1209  indices_begin++;
1210  values_begin++;
1211  }
1212 }
1213 
1214 
1215 
1216 template <typename Number>
1217 inline Vector<Number> &
1218 Vector<Number>::operator/=(const Number factor)
1219 {
1220  AssertIsFinite(factor);
1221  Assert(factor != Number(0.), ExcZero());
1222 
1223  this->operator*=(Number(1.) / factor);
1224  return *this;
1225 }
1226 
1227 
1228 
1229 template <typename Number>
1230 template <typename OtherNumber>
1231 inline void
1232 Vector<Number>::add(const std::vector<size_type> & indices,
1233  const std::vector<OtherNumber> &values)
1234 {
1235  Assert(indices.size() == values.size(),
1236  ExcDimensionMismatch(indices.size(), values.size()));
1237  add(indices.size(), indices.data(), values.data());
1238 }
1239 
1240 
1241 
1242 template <typename Number>
1243 template <typename OtherNumber>
1244 inline void
1245 Vector<Number>::add(const std::vector<size_type> &indices,
1246  const Vector<OtherNumber> & values)
1247 {
1248  Assert(indices.size() == values.size(),
1249  ExcDimensionMismatch(indices.size(), values.size()));
1250  add(indices.size(), indices.data(), values.values.begin());
1251 }
1252 
1253 
1254 
1255 template <typename Number>
1256 template <typename OtherNumber>
1257 inline void
1258 Vector<Number>::add(const size_type n_indices,
1259  const size_type * indices,
1260  const OtherNumber *values)
1261 {
1262  for (size_type i = 0; i < n_indices; ++i)
1263  {
1264  AssertIndexRange(indices[i], size());
1265  Assert(
1266  numbers::is_finite(values[i]),
1267  ExcMessage(
1268  "The given value is not finite but either infinite or Not A Number (NaN)"));
1269 
1270  this->values[indices[i]] += values[i];
1271  }
1272 }
1273 
1274 
1275 
1276 template <typename Number>
1277 template <typename Number2>
1278 inline bool
1279 Vector<Number>::operator!=(const Vector<Number2> &v) const
1280 {
1281  return !(*this == v);
1282 }
1283 
1284 
1285 
1286 template <typename Number>
1288 {}
1289 
1290 
1291 template <typename Number>
1292 inline bool
1294 {
1295  return false;
1296 }
1297 
1298 template <typename Number>
1299 inline void
1301 {}
1302 
1303 
1304 
1305 // Moved from vector.templates.h as an inline function by Luca Heltai
1306 // on 2009/04/12 to prevent strange compiling errors, after making
1307 // swap virtual.
1308 template <typename Number>
1309 inline void
1310 Vector<Number>::swap(Vector<Number> &v)
1311 {
1312  values.swap(v.values);
1313  std::swap(thread_loop_partitioner, v.thread_loop_partitioner);
1314 }
1315 
1316 
1317 
1318 template <typename Number>
1319 template <class Archive>
1320 inline void
1321 Vector<Number>::save(Archive &ar, const unsigned int) const
1322 {
1323  // forward to serialization function in the base class.
1324  ar &static_cast<const Subscriptor &>(*this);
1325  ar &values;
1326 }
1327 
1328 
1329 
1330 template <typename Number>
1331 template <class Archive>
1332 inline void
1333 Vector<Number>::load(Archive &ar, const unsigned int)
1334 {
1335  // the load stuff again from the archive
1336  ar &static_cast<Subscriptor &>(*this);
1337  ar &values;
1339 }
1340 
1341 #endif
1342 
1343 
1356 template <typename Number>
1357 inline void
1358 swap(Vector<Number> &u, Vector<Number> &v)
1359 {
1360  u.swap(v);
1361 }
1362 
1363 
1371 template <typename number>
1372 inline std::ostream &
1373 operator<<(std::ostream &out, const Vector<number> &v)
1374 {
1375  Assert(v.size() != 0, ExcEmptyObject());
1376  AssertThrow(out, ExcIO());
1377 
1378  for (typename Vector<number>::size_type i = 0; i < v.size() - 1; ++i)
1379  out << v(i) << ' ';
1380  out << v(v.size() - 1);
1381 
1382  AssertThrow(out, ExcIO());
1383 
1384  return out;
1385 }
1386 
1395 template <typename Number>
1396 struct is_serial_vector<Vector<Number>> : std::true_type
1397 {};
1398 
1399 
1401 
1402 #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:951
#define AssertIndexRange(index, range)
Definition: exceptions.h:1648
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:1533
pointer data()
iterator end()
bool operator==(const AlignedVector< T > &lhs, const AlignedVector< T > &rhs)
bool is_finite(const double x)
Definition: numbers.h:527
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:1423
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:2816
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:369
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:4334
void swap(Vector< Number > &u, Vector< Number > &v)
Definition: vector.h:1358
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:368
VectorType::value_type * begin(VectorType &V)
iterator begin()
Number l1_norm(const Tensor< 2, dim, Number > &t)
Definition: tensor.h:2790
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:993
size_type size() const
std::shared_ptr< parallel::internal::TBBPartitioner > thread_loop_partitioner
Definition: vector.h:1016
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:1679
std::enable_if< std::is_fundamental< T >::value, std::size_t >::type memory_consumption(const T &t)