Reference documentation for deal.II version 9.3.3
\(\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\}}\)
block_vector_base.h
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2004 - 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_block_vector_base_h
17#define dealii_block_vector_base_h
18
19
20#include <deal.II/base/config.h>
21
26
29#include <deal.II/lac/vector.h>
31
32#include <cmath>
33#include <cstddef>
34#include <iterator>
35#include <type_traits>
36#include <vector>
37
39
40
45// Forward declaration
46#ifndef DOXYGEN
47template <typename>
48class BlockVectorBase;
49#endif
50
65template <typename VectorType>
67{
68private:
73 template <typename T>
74 static std::true_type
76
81 static std::false_type
83
84public:
90 static const bool value =
91 std::is_same<decltype(check_for_block_vector(std::declval<VectorType *>())),
92 std::true_type>::value;
93};
94
95
96// instantiation of the static member
97template <typename VectorType>
99
100
101
102namespace internal
103{
107 namespace BlockVectorIterators
108 {
124 template <class BlockVectorType, bool Constness>
126 {
127 public:
132
139 typename std::conditional<Constness,
140 const typename BlockVectorType::value_type,
141 typename BlockVectorType::value_type>::type;
142
148 using iterator_category = std::random_access_iterator_tag;
149 using difference_type = std::ptrdiff_t;
150 using reference = typename BlockVectorType::reference;
152
153 using dereference_type = typename std::conditional<
154 Constness,
156 typename BlockVectorType::BlockType::reference>::type;
157
162 using BlockVector = typename std::
163 conditional<Constness, const BlockVectorType, BlockVectorType>::type;
164
174
183
184
189
190 private:
201
202 public:
206 Iterator &
208
215
221
226 Iterator &
228
236
241 Iterator &
243
251
256 template <bool OtherConstness>
257 bool
259
264 template <bool OtherConstness>
265 bool
267
273 template <bool OtherConstness>
274 bool
276
280 template <bool OtherConstness>
281 bool
283
287 template <bool OtherConstness>
288 bool
290
294 template <bool OtherConstness>
295 bool
297
301 template <bool OtherConstness>
304
311
318
323 Iterator &
325
330 Iterator &
332
343 "Your program tried to compare iterators pointing to "
344 "different block vectors. There is no reasonable way "
345 "to do this.");
346
348 private:
355
360
365 unsigned int current_block;
367
376
380 void
382
386 void
388
389 // Mark all other instances of this template as friends.
390 template <typename, bool>
391 friend class Iterator;
392 };
393 } // namespace BlockVectorIterators
394} // namespace internal
395
396
434template <class VectorType>
436{
437public:
441 using BlockType = VectorType;
442
443 /*
444 * Declare standard types used in
445 * all containers. These types
446 * parallel those in the
447 * <tt>C++</tt> standard
448 * libraries
449 * <tt>std::vector<...></tt>
450 * class. This includes iterator
451 * types.
452 */
453 using value_type = typename BlockType::value_type;
455 using const_pointer = const value_type *;
456 using iterator =
460 using reference = typename BlockType::reference;
461 using const_reference = typename BlockType::const_reference;
463
473 using real_type = typename BlockType::real_type;
474
478 BlockVectorBase() = default;
479
483 BlockVectorBase(const BlockVectorBase & /*V*/) = default;
484
490 BlockVectorBase(BlockVectorBase && /*V*/) noexcept = default;
491
498 void
500
511 void
513
517 BlockType &
518 block(const unsigned int i);
519
523 const BlockType &
524 block(const unsigned int i) const;
525
531 const BlockIndices &
533
537 unsigned int
538 n_blocks() const;
539
544 std::size_t
545 size() const;
546
552 std::size_t
554
573
579
585 begin() const;
586
592
598 end() const;
599
604 operator()(const size_type i) const;
605
610 operator()(const size_type i);
611
617 value_type operator[](const size_type i) const;
618
624 reference operator[](const size_type i);
625
641 template <typename OtherNumber>
642 void
643 extract_subvector_to(const std::vector<size_type> &indices,
644 std::vector<OtherNumber> & values) const;
645
673 template <typename ForwardIterator, typename OutputIterator>
674 void
675 extract_subvector_to(ForwardIterator indices_begin,
676 const ForwardIterator indices_end,
677 OutputIterator values_begin) const;
678
684 operator=(const value_type s);
685
690 operator=(const BlockVectorBase &V);
691
698 operator=(BlockVectorBase && /*V*/) = default; // NOLINT
699
703 template <class VectorType2>
705 operator=(const BlockVectorBase<VectorType2> &V);
706
711 operator=(const VectorType &v);
712
717 template <class VectorType2>
718 bool
719 operator==(const BlockVectorBase<VectorType2> &v) const;
720
724 value_type operator*(const BlockVectorBase &V) const;
725
730 norm_sqr() const;
731
736 mean_value() const;
737
742 l1_norm() const;
743
749 l2_norm() const;
750
756 linfty_norm() const;
757
782 const BlockVectorBase &V,
783 const BlockVectorBase &W);
784
789 bool
791
797 bool
798 all_zero() const;
799
805 bool
807
812 operator+=(const BlockVectorBase &V);
813
818 operator-=(const BlockVectorBase &V);
819
820
825 template <typename Number>
826 void
827 add(const std::vector<size_type> &indices, const std::vector<Number> &values);
828
833 template <typename Number>
834 void
835 add(const std::vector<size_type> &indices, const Vector<Number> &values);
836
842 template <typename Number>
843 void
844 add(const size_type n_elements,
845 const size_type *indices,
846 const Number * values);
847
852 void
853 add(const value_type s);
854
858 void
859 add(const value_type a, const BlockVectorBase &V);
860
864 void
866 const BlockVectorBase &V,
867 const value_type b,
868 const BlockVectorBase &W);
869
873 void
875
879 void
880 sadd(const value_type s, const value_type a, const BlockVectorBase &V);
881
885 void
887 const value_type a,
888 const BlockVectorBase &V,
889 const value_type b,
890 const BlockVectorBase &W);
891
895 void
897 const value_type a,
898 const BlockVectorBase &V,
899 const value_type b,
900 const BlockVectorBase &W,
901 const value_type c,
902 const BlockVectorBase &X);
903
908 operator*=(const value_type factor);
909
914 operator/=(const value_type factor);
915
920 template <class BlockVector2>
921 void
922 scale(const BlockVector2 &v);
923
927 template <class BlockVector2>
928 void
929 equ(const value_type a, const BlockVector2 &V);
930
935 void
937
942 std::size_t
944
945protected:
949 std::vector<VectorType> components;
950
956
957 // Make the iterator class a friend.
958 template <typename N, bool C>
959 friend class ::internal::BlockVectorIterators::Iterator;
960
961 template <typename>
962 friend class BlockVectorBase;
963};
964
965
968/*----------------------- Inline functions ----------------------------------*/
969
970
971#ifndef DOXYGEN
972namespace internal
973{
974 namespace BlockVectorIterators
975 {
976 template <class BlockVectorType, bool Constness>
977 inline Iterator<BlockVectorType, Constness>::Iterator(
978 const Iterator<BlockVectorType, Constness> &c)
979 : parent(c.parent)
980 , global_index(c.global_index)
981 , current_block(c.current_block)
982 , index_within_block(c.index_within_block)
983 , next_break_forward(c.next_break_forward)
984 , next_break_backward(c.next_break_backward)
985 {}
986
987
988
989 template <class BlockVectorType, bool Constness>
990 inline Iterator<BlockVectorType, Constness>::Iterator(
991 const Iterator<BlockVectorType, !Constness> &c)
992 : parent(c.parent)
993 , global_index(c.global_index)
994 , current_block(c.current_block)
995 , index_within_block(c.index_within_block)
996 , next_break_forward(c.next_break_forward)
997 , next_break_backward(c.next_break_backward)
998 {
999 // Only permit copy-constructing const iterators from non-const
1000 // iterators, and not vice versa (i.e., Constness must always be
1001 // true).
1002 static_assert(Constness == true,
1003 "Constructing a non-const iterator from a const iterator "
1004 "does not make sense.");
1005 }
1006
1007
1008
1009 template <class BlockVectorType, bool Constness>
1010 inline Iterator<BlockVectorType, Constness>::Iterator(
1011 BlockVector & parent,
1012 const size_type global_index,
1013 const size_type current_block,
1014 const size_type index_within_block,
1015 const size_type next_break_forward,
1016 const size_type next_break_backward)
1017 : parent(&parent)
1019 , current_block(current_block)
1020 , index_within_block(index_within_block)
1021 , next_break_forward(next_break_forward)
1022 , next_break_backward(next_break_backward)
1023 {}
1024
1025
1026
1027 template <class BlockVectorType, bool Constness>
1028 inline Iterator<BlockVectorType, Constness> &
1029 Iterator<BlockVectorType, Constness>::operator=(const Iterator &c)
1030 {
1031 parent = c.parent;
1032 global_index = c.global_index;
1033 index_within_block = c.index_within_block;
1034 current_block = c.current_block;
1035 next_break_forward = c.next_break_forward;
1036 next_break_backward = c.next_break_backward;
1037
1038 return *this;
1039 }
1040
1041
1042
1043 template <class BlockVectorType, bool Constness>
1044 inline typename Iterator<BlockVectorType, Constness>::dereference_type
1046 {
1047 return parent->block(current_block)(index_within_block);
1048 }
1049
1050
1051
1052 template <class BlockVectorType, bool Constness>
1053 inline typename Iterator<BlockVectorType, Constness>::dereference_type
1054 Iterator<BlockVectorType, Constness>::
1055 operator[](const difference_type d) const
1056 {
1057 // if the index pointed to is
1058 // still within the block we
1059 // currently point into, then we
1060 // can save the computation of
1061 // the block
1062 if ((global_index + d >= next_break_backward) &&
1063 (global_index + d <= next_break_forward))
1064 return parent->block(current_block)(index_within_block + d);
1065
1066 // if the index is not within the
1067 // block of the block vector into
1068 // which we presently point, then
1069 // there is no way: we have to
1070 // search for the block. this can
1071 // be done through the parent
1072 // class as well.
1073 return (*parent)(global_index + d);
1074 }
1075
1076
1077
1078 template <class BlockVectorType, bool Constness>
1079 inline Iterator<BlockVectorType, Constness> &
1081 {
1082 move_forward();
1083 return *this;
1084 }
1085
1086
1087
1088 template <class BlockVectorType, bool Constness>
1089 inline Iterator<BlockVectorType, Constness>
1091 {
1092 const Iterator old_value = *this;
1093 move_forward();
1094 return old_value;
1095 }
1096
1097
1098
1099 template <class BlockVectorType, bool Constness>
1100 inline Iterator<BlockVectorType, Constness> &
1101 Iterator<BlockVectorType, Constness>::operator--()
1102 {
1103 move_backward();
1104 return *this;
1105 }
1106
1107
1108
1109 template <class BlockVectorType, bool Constness>
1110 inline Iterator<BlockVectorType, Constness>
1111 Iterator<BlockVectorType, Constness>::operator--(int)
1112 {
1113 const Iterator old_value = *this;
1114 move_backward();
1115 return old_value;
1116 }
1117
1118
1119
1120 template <class BlockVectorType, bool Constness>
1121 template <bool OtherConstness>
1122 inline bool
1124 operator==(const Iterator<BlockVectorType, OtherConstness> &i) const
1125 {
1126 Assert(parent == i.parent, ExcPointerToDifferentVectors());
1127
1128 return (global_index == i.global_index);
1129 }
1130
1131
1132
1133 template <class BlockVectorType, bool Constness>
1134 template <bool OtherConstness>
1135 inline bool
1137 operator!=(const Iterator<BlockVectorType, OtherConstness> &i) const
1138 {
1139 Assert(parent == i.parent, ExcPointerToDifferentVectors());
1140
1141 return (global_index != i.global_index);
1142 }
1143
1144
1145
1146 template <class BlockVectorType, bool Constness>
1147 template <bool OtherConstness>
1148 inline bool
1150 operator<(const Iterator<BlockVectorType, OtherConstness> &i) const
1151 {
1152 Assert(parent == i.parent, ExcPointerToDifferentVectors());
1153
1154 return (global_index < i.global_index);
1155 }
1156
1157
1158
1159 template <class BlockVectorType, bool Constness>
1160 template <bool OtherConstness>
1161 inline bool
1163 operator<=(const Iterator<BlockVectorType, OtherConstness> &i) const
1164 {
1165 Assert(parent == i.parent, ExcPointerToDifferentVectors());
1166
1167 return (global_index <= i.global_index);
1168 }
1169
1170
1171
1172 template <class BlockVectorType, bool Constness>
1173 template <bool OtherConstness>
1174 inline bool
1176 operator>(const Iterator<BlockVectorType, OtherConstness> &i) const
1177 {
1178 Assert(parent == i.parent, ExcPointerToDifferentVectors());
1179
1180 return (global_index > i.global_index);
1181 }
1182
1183
1184
1185 template <class BlockVectorType, bool Constness>
1186 template <bool OtherConstness>
1187 inline bool
1189 operator>=(const Iterator<BlockVectorType, OtherConstness> &i) const
1190 {
1191 Assert(parent == i.parent, ExcPointerToDifferentVectors());
1192
1193 return (global_index >= i.global_index);
1194 }
1195
1196
1197
1198 template <class BlockVectorType, bool Constness>
1199 template <bool OtherConstness>
1200 inline typename Iterator<BlockVectorType, Constness>::difference_type
1202 operator-(const Iterator<BlockVectorType, OtherConstness> &i) const
1203 {
1204 Assert(parent == i.parent, ExcPointerToDifferentVectors());
1205
1206 return (static_cast<signed int>(global_index) -
1207 static_cast<signed int>(i.global_index));
1208 }
1209
1210
1211
1212 template <class BlockVectorType, bool Constness>
1213 inline Iterator<BlockVectorType, Constness>
1215 operator+(const difference_type &d) const
1216 {
1217 // if the index pointed to is
1218 // still within the block we
1219 // currently point into, then we
1220 // can save the computation of
1221 // the block
1222 if ((global_index + d >= next_break_backward) &&
1223 (global_index + d <= next_break_forward))
1224 return Iterator(*parent,
1225 global_index + d,
1226 current_block,
1227 index_within_block + d,
1228 next_break_forward,
1229 next_break_backward);
1230 else
1231 // outside present block, so
1232 // have to seek new block
1233 // anyway
1234 return Iterator(*parent, global_index + d);
1235 }
1236
1237
1238
1239 template <class BlockVectorType, bool Constness>
1240 inline Iterator<BlockVectorType, Constness>
1242 operator-(const difference_type &d) const
1243 {
1244 // if the index pointed to is
1245 // still within the block we
1246 // currently point into, then we
1247 // can save the computation of
1248 // the block
1249 if ((global_index - d >= next_break_backward) &&
1250 (global_index - d <= next_break_forward))
1251 return Iterator(*parent,
1252 global_index - d,
1253 current_block,
1254 index_within_block - d,
1255 next_break_forward,
1256 next_break_backward);
1257 else
1258 // outside present block, so
1259 // have to seek new block
1260 // anyway
1261 return Iterator(*parent, global_index - d);
1262 }
1263
1264
1265
1266 template <class BlockVectorType, bool Constness>
1267 inline Iterator<BlockVectorType, Constness> &
1268 Iterator<BlockVectorType, Constness>::operator+=(const difference_type &d)
1269 {
1270 // if the index pointed to is
1271 // still within the block we
1272 // currently point into, then we
1273 // can save the computation of
1274 // the block
1275 if ((global_index + d >= next_break_backward) &&
1276 (global_index + d <= next_break_forward))
1277 {
1278 global_index += d;
1279 index_within_block += d;
1280 }
1281 else
1282 // outside present block, so
1283 // have to seek new block
1284 // anyway
1285 *this = Iterator(*parent, global_index + d);
1286
1287 return *this;
1288 }
1289
1290
1291
1292 template <class BlockVectorType, bool Constness>
1293 inline Iterator<BlockVectorType, Constness> &
1294 Iterator<BlockVectorType, Constness>::operator-=(const difference_type &d)
1295 {
1296 // if the index pointed to is
1297 // still within the block we
1298 // currently point into, then we
1299 // can save the computation of
1300 // the block
1301 if ((global_index - d >= next_break_backward) &&
1302 (global_index - d <= next_break_forward))
1303 {
1304 global_index -= d;
1305 index_within_block -= d;
1306 }
1307 else
1308 // outside present block, so
1309 // have to seek new block
1310 // anyway
1311 *this = Iterator(*parent, global_index - d);
1312
1313 return *this;
1314 }
1315
1316
1317 template <class BlockVectorType, bool Constness>
1318 Iterator<BlockVectorType, Constness>::Iterator(BlockVector & parent,
1319 const size_type global_index)
1320 : parent(&parent)
1322 {
1323 // find which block we are
1324 // in. for this, take into
1325 // account that it happens at
1326 // times that people want to
1327 // initialize iterators
1328 // past-the-end
1329 if (global_index < parent.size())
1330 {
1331 const std::pair<size_type, size_type> indices =
1333 current_block = indices.first;
1334 index_within_block = indices.second;
1335
1336 next_break_backward =
1337 parent.block_indices.local_to_global(current_block, 0);
1338 next_break_forward =
1339 (parent.block_indices.local_to_global(current_block, 0) +
1340 parent.block_indices.block_size(current_block) - 1);
1341 }
1342 else
1343 // past the end. only have one
1344 // value for this
1345 {
1346 this->global_index = parent.size();
1347 current_block = parent.n_blocks();
1348 index_within_block = 0;
1349 next_break_backward = global_index;
1350 next_break_forward = numbers::invalid_size_type;
1351 };
1352 }
1353
1354
1355
1356 template <class BlockVectorType, bool Constness>
1357 void
1358 Iterator<BlockVectorType, Constness>::move_forward()
1359 {
1360 if (global_index != next_break_forward)
1361 ++index_within_block;
1362 else
1363 {
1364 // ok, we traverse a boundary
1365 // between blocks:
1366 index_within_block = 0;
1367 ++current_block;
1368
1369 // break backwards is now old
1370 // break forward
1371 next_break_backward = next_break_forward + 1;
1372
1373 // compute new break forward
1374 if (current_block < parent->block_indices.size())
1375 next_break_forward +=
1376 parent->block_indices.block_size(current_block);
1377 else
1378 // if we are beyond the end,
1379 // then move the next
1380 // boundary arbitrarily far
1381 // away
1382 next_break_forward = numbers::invalid_size_type;
1383 };
1384
1385 ++global_index;
1386 }
1387
1388
1389
1390 template <class BlockVectorType, bool Constness>
1391 void
1392 Iterator<BlockVectorType, Constness>::move_backward()
1393 {
1394 if (global_index != next_break_backward)
1395 --index_within_block;
1396 else if (current_block != 0)
1397 {
1398 // ok, we traverse a boundary
1399 // between blocks:
1400 --current_block;
1401 index_within_block =
1402 parent->block_indices.block_size(current_block) - 1;
1403
1404 // break forwards is now old
1405 // break backward
1406 next_break_forward = next_break_backward - 1;
1407
1408 // compute new break forward
1409 next_break_backward -=
1410 parent->block_indices.block_size(current_block);
1411 }
1412 else
1413 // current block was 0, we now
1414 // get into unspecified terrain
1415 {
1416 --current_block;
1417 index_within_block = numbers::invalid_size_type;
1418 next_break_forward = 0;
1419 next_break_backward = 0;
1420 };
1421
1422 --global_index;
1423 }
1424
1425
1426 } // namespace BlockVectorIterators
1427
1428} // namespace internal
1429
1430
1431
1432template <class VectorType>
1433inline std::size_t
1435{
1436 return block_indices.total_size();
1437}
1438
1439
1440
1441template <class VectorType>
1442inline std::size_t
1444{
1445 std::size_t local_size = 0;
1446 for (unsigned int b = 0; b < n_blocks(); ++b)
1447 local_size += block(b).locally_owned_size();
1448 return local_size;
1449}
1450
1451
1452
1453template <class VectorType>
1454inline IndexSet
1456{
1457 IndexSet is(size());
1458
1459 // copy index sets from blocks into the global one, shifted
1460 // by the appropriate amount for each block
1461 for (unsigned int b = 0; b < n_blocks(); ++b)
1462 {
1463 IndexSet x = block(b).locally_owned_elements();
1464 is.add_indices(x, block_indices.block_start(b));
1465 }
1466
1467 is.compress();
1468
1469 return is;
1470}
1471
1472
1473
1474template <class VectorType>
1475inline unsigned int
1477{
1478 return block_indices.size();
1479}
1480
1481
1482template <class VectorType>
1484BlockVectorBase<VectorType>::block(const unsigned int i)
1485{
1487
1488 return components[i];
1489}
1490
1491
1492
1493template <class VectorType>
1494inline const typename BlockVectorBase<VectorType>::BlockType &
1495BlockVectorBase<VectorType>::block(const unsigned int i) const
1496{
1498
1499 return components[i];
1500}
1501
1502
1503
1504template <class VectorType>
1505inline const BlockIndices &
1507{
1508 return block_indices;
1509}
1510
1511
1512template <class VectorType>
1513inline void
1515{
1516 std::vector<size_type> sizes(n_blocks());
1517
1518 for (size_type i = 0; i < n_blocks(); ++i)
1519 sizes[i] = block(i).size();
1520
1521 block_indices.reinit(sizes);
1522}
1523
1524
1525
1526template <class VectorType>
1527inline void
1529 ::VectorOperation::values operation)
1530{
1531 for (unsigned int i = 0; i < n_blocks(); ++i)
1532 block(i).compress(operation);
1533}
1534
1535
1536
1537template <class VectorType>
1540{
1541 return iterator(*this, 0U);
1542}
1543
1544
1545
1546template <class VectorType>
1549{
1550 return const_iterator(*this, 0U);
1551}
1552
1553
1554template <class VectorType>
1557{
1558 return iterator(*this, size());
1559}
1560
1561
1562
1563template <class VectorType>
1566{
1567 return const_iterator(*this, size());
1568}
1569
1570
1571template <class VectorType>
1572inline bool
1574{
1575 const std::pair<size_type, size_type> local_index =
1577
1578 return components[local_index.first].in_local_range(global_index);
1579}
1580
1581
1582template <class VectorType>
1583bool
1585{
1586 for (size_type i = 0; i < n_blocks(); ++i)
1587 if (components[i].all_zero() == false)
1588 return false;
1589
1590 return true;
1591}
1592
1593
1594
1595template <class VectorType>
1596bool
1598{
1599 for (size_type i = 0; i < n_blocks(); ++i)
1600 if (components[i].is_non_negative() == false)
1601 return false;
1602
1603 return true;
1604}
1605
1606
1607
1608template <class VectorType>
1611{
1612 Assert(n_blocks() == v.n_blocks(),
1614
1615 value_type sum = 0.;
1616 for (size_type i = 0; i < n_blocks(); ++i)
1617 sum += components[i] * v.components[i];
1618
1619 return sum;
1620}
1621
1622
1623template <class VectorType>
1626{
1627 real_type sum = 0.;
1628 for (size_type i = 0; i < n_blocks(); ++i)
1629 sum += components[i].norm_sqr();
1630
1631 return sum;
1632}
1633
1634
1635
1636template <class VectorType>
1639{
1640 value_type sum = 0.;
1641 // need to do static_cast as otherwise it won't work with
1642 // value_type=complex<T>
1643 for (size_type i = 0; i < n_blocks(); ++i)
1644 sum += components[i].mean_value() *
1646 components[i].size()));
1647
1649}
1650
1651
1652
1653template <class VectorType>
1656{
1657 real_type sum = 0.;
1658 for (size_type i = 0; i < n_blocks(); ++i)
1659 sum += components[i].l1_norm();
1660
1661 return sum;
1662}
1663
1664
1665
1666template <class VectorType>
1669{
1670 return std::sqrt(norm_sqr());
1671}
1672
1673
1674
1675template <class VectorType>
1678{
1679 real_type sum = 0.;
1680 for (size_type i = 0; i < n_blocks(); ++i)
1681 {
1682 value_type newval = components[i].linfty_norm();
1683 if (sum < newval)
1684 sum = newval;
1685 }
1686 return sum;
1687}
1688
1689
1690
1691template <class VectorType>
1697{
1698 AssertDimension(n_blocks(), V.n_blocks());
1700
1701 value_type sum = 0.;
1702 for (size_type i = 0; i < n_blocks(); ++i)
1703 sum += components[i].add_and_dot(a, V.components[i], W.components[i]);
1704
1705 return sum;
1706}
1707
1708
1709
1710template <class VectorType>
1713{
1714 Assert(n_blocks() == v.n_blocks(),
1716
1717 for (size_type i = 0; i < n_blocks(); ++i)
1718 {
1719 components[i] += v.components[i];
1720 }
1721
1722 return *this;
1723}
1724
1725
1726
1727template <class VectorType>
1730{
1731 Assert(n_blocks() == v.n_blocks(),
1733
1734 for (size_type i = 0; i < n_blocks(); ++i)
1735 {
1736 components[i] -= v.components[i];
1737 }
1738 return *this;
1739}
1740
1741
1742
1743template <class VectorType>
1744template <typename Number>
1745inline void
1746BlockVectorBase<VectorType>::add(const std::vector<size_type> &indices,
1747 const std::vector<Number> & values)
1748{
1749 Assert(indices.size() == values.size(),
1750 ExcDimensionMismatch(indices.size(), values.size()));
1751 add(indices.size(), indices.data(), values.data());
1752}
1753
1754
1755
1756template <class VectorType>
1757template <typename Number>
1758inline void
1759BlockVectorBase<VectorType>::add(const std::vector<size_type> &indices,
1760 const Vector<Number> & values)
1761{
1762 Assert(indices.size() == values.size(),
1763 ExcDimensionMismatch(indices.size(), values.size()));
1764 const size_type n_indices = indices.size();
1765 for (size_type i = 0; i < n_indices; ++i)
1766 (*this)(indices[i]) += values(i);
1767}
1768
1769
1770
1771template <class VectorType>
1772template <typename Number>
1773inline void
1775 const size_type *indices,
1776 const Number * values)
1777{
1778 for (size_type i = 0; i < n_indices; ++i)
1779 (*this)(indices[i]) += values[i];
1780}
1781
1782
1783
1784template <class VectorType>
1785void
1787{
1788 AssertIsFinite(a);
1789
1790 for (size_type i = 0; i < n_blocks(); ++i)
1791 {
1792 components[i].add(a);
1793 }
1794}
1795
1796
1797
1798template <class VectorType>
1799void
1802{
1803 AssertIsFinite(a);
1804
1805 Assert(n_blocks() == v.n_blocks(),
1807
1808 for (size_type i = 0; i < n_blocks(); ++i)
1809 {
1810 components[i].add(a, v.components[i]);
1811 }
1812}
1813
1814
1815
1816template <class VectorType>
1817void
1820 const value_type b,
1822{
1823 AssertIsFinite(a);
1825
1826 Assert(n_blocks() == v.n_blocks(),
1828 Assert(n_blocks() == w.n_blocks(),
1829 ExcDimensionMismatch(n_blocks(), w.n_blocks()));
1830
1831
1832 for (size_type i = 0; i < n_blocks(); ++i)
1833 {
1834 components[i].add(a, v.components[i], b, w.components[i]);
1835 }
1836}
1837
1838
1839
1840template <class VectorType>
1841void
1844{
1845 AssertIsFinite(x);
1846
1847 Assert(n_blocks() == v.n_blocks(),
1849
1850 for (size_type i = 0; i < n_blocks(); ++i)
1851 {
1852 components[i].sadd(x, v.components[i]);
1853 }
1854}
1855
1856
1857
1858template <class VectorType>
1859void
1861 const value_type a,
1863{
1864 AssertIsFinite(x);
1865 AssertIsFinite(a);
1866
1867 Assert(n_blocks() == v.n_blocks(),
1869
1870 for (size_type i = 0; i < n_blocks(); ++i)
1871 {
1872 components[i].sadd(x, a, v.components[i]);
1873 }
1874}
1875
1876
1877
1878template <class VectorType>
1879void
1881 const value_type a,
1883 const value_type b,
1885{
1886 AssertIsFinite(x);
1887 AssertIsFinite(a);
1889
1890 Assert(n_blocks() == v.n_blocks(),
1892 Assert(n_blocks() == w.n_blocks(),
1893 ExcDimensionMismatch(n_blocks(), w.n_blocks()));
1894
1895 for (size_type i = 0; i < n_blocks(); ++i)
1896 {
1897 components[i].sadd(x, a, v.components[i], b, w.components[i]);
1898 }
1899}
1900
1901
1902
1903template <class VectorType>
1904void
1906 const value_type a,
1908 const value_type b,
1910 const value_type c,
1912{
1913 AssertIsFinite(x);
1914 AssertIsFinite(a);
1916 AssertIsFinite(c);
1917
1918 Assert(n_blocks() == v.n_blocks(),
1920 Assert(n_blocks() == w.n_blocks(),
1921 ExcDimensionMismatch(n_blocks(), w.n_blocks()));
1922 Assert(n_blocks() == y.n_blocks(),
1924
1925 for (size_type i = 0; i < n_blocks(); ++i)
1926 {
1927 components[i].sadd(
1928 x, a, v.components[i], b, w.components[i], c, y.components[i]);
1929 }
1930}
1931
1932
1933
1934template <class VectorType>
1935template <class BlockVector2>
1936void
1937BlockVectorBase<VectorType>::scale(const BlockVector2 &v)
1938{
1939 Assert(n_blocks() == v.n_blocks(),
1940 ExcDimensionMismatch(n_blocks(), v.n_blocks()));
1941 for (size_type i = 0; i < n_blocks(); ++i)
1942 components[i].scale(v.block(i));
1943}
1944
1945
1946
1947template <class VectorType>
1948std::size_t
1950{
1951 return (MemoryConsumption::memory_consumption(this->block_indices) +
1952 MemoryConsumption::memory_consumption(this->components));
1953}
1954
1955
1956
1957template <class VectorType>
1958template <class BlockVector2>
1959void
1960BlockVectorBase<VectorType>::equ(const value_type a, const BlockVector2 &v)
1961{
1962 AssertIsFinite(a);
1963
1964 Assert(n_blocks() == v.n_blocks(),
1965 ExcDimensionMismatch(n_blocks(), v.n_blocks()));
1966
1967 for (size_type i = 0; i < n_blocks(); ++i)
1968 components[i].equ(a, v.components[i]);
1969}
1970
1971
1972
1973template <class VectorType>
1974void
1976{
1977 for (size_type i = 0; i < n_blocks(); ++i)
1978 block(i).update_ghost_values();
1979}
1980
1981
1982
1983template <class VectorType>
1986{
1987 AssertIsFinite(s);
1988
1989 for (size_type i = 0; i < n_blocks(); ++i)
1990 components[i] = s;
1991
1992 return *this;
1993}
1994
1995
1996template <class VectorType>
1999{
2001
2002 for (size_type i = 0; i < n_blocks(); ++i)
2003 components[i] = v.components[i];
2004
2005 return *this;
2006}
2007
2008
2009template <class VectorType>
2010template <class VectorType2>
2013{
2015
2016 for (size_type i = 0; i < n_blocks(); ++i)
2017 components[i] = v.components[i];
2018
2019 return *this;
2020}
2021
2022
2023
2024template <class VectorType>
2026BlockVectorBase<VectorType>::operator=(const VectorType &v)
2027{
2028 Assert(size() == v.size(), ExcDimensionMismatch(size(), v.size()));
2029
2030 size_type index_v = 0;
2031 for (size_type b = 0; b < n_blocks(); ++b)
2032 for (size_type i = 0; i < block(b).size(); ++i, ++index_v)
2033 block(b)(i) = v(index_v);
2034
2035 return *this;
2036}
2037
2038
2039
2040template <class VectorType>
2041template <class VectorType2>
2042inline bool
2045{
2047
2048 for (size_type i = 0; i < n_blocks(); ++i)
2049 if (!(components[i] == v.components[i]))
2050 return false;
2051
2052 return true;
2053}
2054
2055
2056
2057template <class VectorType>
2060{
2061 AssertIsFinite(factor);
2062
2063 for (size_type i = 0; i < n_blocks(); ++i)
2064 components[i] *= factor;
2065
2066 return *this;
2067}
2068
2069
2070
2071template <class VectorType>
2074{
2075 AssertIsFinite(factor);
2076 Assert(factor != 0., ExcDivideByZero());
2077
2078 const value_type inverse_factor = value_type(1.) / factor;
2079
2080 for (size_type i = 0; i < n_blocks(); ++i)
2081 components[i] *= inverse_factor;
2082
2083 return *this;
2084}
2085
2086
2087template <class VectorType>
2090{
2091 const std::pair<unsigned int, size_type> local_index =
2093 return components[local_index.first](local_index.second);
2094}
2095
2096
2097
2098template <class VectorType>
2101{
2102 const std::pair<unsigned int, size_type> local_index =
2104 return components[local_index.first](local_index.second);
2105}
2106
2107
2108
2109template <class VectorType>
2112{
2113 return operator()(i);
2114}
2115
2116
2117
2118template <class VectorType>
2121{
2122 return operator()(i);
2123}
2124
2125
2126
2127template <typename VectorType>
2128template <typename OtherNumber>
2129inline void
2131 const std::vector<size_type> &indices,
2132 std::vector<OtherNumber> & values) const
2133{
2134 for (size_type i = 0; i < indices.size(); ++i)
2135 values[i] = operator()(indices[i]);
2136}
2137
2138
2139
2140template <typename VectorType>
2141template <typename ForwardIterator, typename OutputIterator>
2142inline void
2144 ForwardIterator indices_begin,
2145 const ForwardIterator indices_end,
2146 OutputIterator values_begin) const
2147{
2148 while (indices_begin != indices_end)
2149 {
2150 *values_begin = operator()(*indices_begin);
2151 indices_begin++;
2152 values_begin++;
2153 }
2154}
2155
2156#endif // DOXYGEN
2157
2159
2160#endif
bool operator!=(const AlignedVector< T > &lhs, const AlignedVector< T > &rhs)
bool operator==(const AlignedVector< T > &lhs, const AlignedVector< T > &rhs)
size_type block_size(const unsigned int i) const
size_type total_size() const
void reinit(const unsigned int n_blocks, const size_type n_elements_per_block)
unsigned int size() const
size_type local_to_global(const unsigned int block, const size_type index) const
size_type block_start(const unsigned int i) const
std::pair< unsigned int, size_type > global_to_local(const size_type i) const
Definition: vector.h:110
std::enable_if< std::is_floating_point< T >::value &&std::is_floating_point< U >::value, typenameProductType< std::complex< T >, std::complex< U > >::type >::type operator*(const std::complex< T > &left, const std::complex< U > &right)
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:402
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:403
BlockVectorBase & operator+=(const BlockVectorBase &V)
BlockVectorBase()=default
void extract_subvector_to(const std::vector< size_type > &indices, std::vector< OtherNumber > &values) const
value_type mean_value() const
void add(const std::vector< size_type > &indices, const std::vector< Number > &values)
::internal::BlockVectorIterators::Iterator< BlockVectorBase, false > iterator
BlockVectorBase & operator/=(const value_type factor)
unsigned int n_blocks() const
real_type norm_sqr() const
const value_type * const_pointer
typename BlockType::const_reference const_reference
void update_ghost_values() const
std::size_t memory_consumption() const
real_type l1_norm() const
types::global_dof_index size_type
real_type linfty_norm() const
value_type add_and_dot(const value_type a, const BlockVectorBase &V, const BlockVectorBase &W)
typename BlockType::value_type value_type
value_type operator*(const BlockVectorBase &V) const
bool all_zero() const
BlockVectorBase & operator=(const value_type s)
std::size_t size() const
void collect_sizes()
#define Assert(cond, exc)
Definition: exceptions.h:1465
BlockVectorBase & operator*=(const value_type factor)
void sadd(const value_type s, const BlockVectorBase &V)
std::vector< VectorType > components
static ::ExceptionBase & ExcPointerToDifferentVectors()
iterator begin()
Expression operator>=(const Expression &lhs, const Expression &rhs)
bool in_local_range(const size_type global_index) const
#define AssertIsFinite(number)
Definition: exceptions.h:1721
::internal::BlockVectorIterators::Iterator< BlockVectorBase, true > const_iterator
BlockVectorBase(BlockVectorBase &&) noexcept=default
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:1622
BlockVectorBase(const BlockVectorBase &)=default
typename BlockType::reference reference
iterator end()
value_type operator[](const size_type i) const
BlockIndices block_indices
#define AssertIndexRange(index, range)
Definition: exceptions.h:1690
#define DeclExceptionMsg(Exception, defaulttext)
Definition: exceptions.h:493
typename BlockType::real_type real_type
Expression operator>(const Expression &lhs, const Expression &rhs)
void scale(const BlockVector2 &v)
void compress(::VectorOperation::values operation)
BlockType & block(const unsigned int i)
static ::ExceptionBase & ExcDivideByZero()
static ::ExceptionBase & ExcDifferentBlockIndices()
Expression operator<=(const Expression &lhs, const Expression &rhs)
real_type l2_norm() const
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
value_type operator()(const size_type i) const
IndexSet locally_owned_elements() const
bool is_non_negative() const
BlockVectorBase & operator-=(const BlockVectorBase &V)
bool operator==(const BlockVectorBase< VectorType2 > &v) const
std::size_t locally_owned_size() const
const BlockIndices & get_block_indices() const
value_type * pointer
void equ(const value_type a, const BlockVector2 &V)
bool operator==(const Iterator< BlockVectorType, OtherConstness > &i) const
static std::true_type check_for_block_vector(const BlockVectorBase< T > *)
dereference_type operator[](const difference_type d) const
Iterator & operator=(const Iterator &c)
bool operator<(const Iterator< BlockVectorType, OtherConstness > &i) const
Iterator(BlockVector &parent, const size_type global_index, const size_type current_block, const size_type index_within_block, const size_type next_break_forward, const size_type next_break_backward)
Iterator(const Iterator< BlockVectorType, !Constness > &c)
dereference_type operator*() const
Iterator operator-(const difference_type &d) const
difference_type operator-(const Iterator< BlockVectorType, OtherConstness > &i) const
std::random_access_iterator_tag iterator_category
typename std::conditional< Constness, value_type, typename BlockVectorType::BlockType::reference >::type dereference_type
Iterator(BlockVector &parent, const size_type global_index)
bool operator>=(const Iterator< BlockVectorType, OtherConstness > &i) const
static std::false_type check_for_block_vector(...)
typename std::conditional< Constness, const typename BlockVectorType::value_type, typename BlockVectorType::value_type >::type value_type
Iterator & operator-=(const difference_type &d)
bool operator<=(const Iterator< BlockVectorType, OtherConstness > &i) const
typename BlockVectorType::reference reference
static const bool value
Iterator operator+(const difference_type &d) const
Iterator & operator+=(const difference_type &d)
bool operator!=(const Iterator< BlockVectorType, OtherConstness > &i) const
bool operator>(const Iterator< BlockVectorType, OtherConstness > &i) const
static const char U
static const char N
static const char V
std::enable_if< std::is_fundamental< T >::value, std::size_t >::type memory_consumption(const T &t)
SymmetricTensor< 2, dim, Number > C(const Tensor< 2, dim, Number > &F)
Tensor< 2, dim, Number > w(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
SymmetricTensor< 2, dim, Number > b(const Tensor< 2, dim, Number > &F)
SymmetricTensor< 2, dim, Number > d(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
TrilinosWrappers::types::int_type global_index(const Epetra_BlockMap &map, const ::types::global_dof_index i)
T sum(const T &t, const MPI_Comm &mpi_communicator)
const types::global_dof_index invalid_size_type
Definition: types.h:205
STL namespace.
::VectorizedArray< Number, width > sqrt(const ::VectorizedArray< Number, width > &)
unsigned int global_dof_index
Definition: types.h:76
BarycentricPolynomial< dim, Number1 > operator-(const Number2 &a, const BarycentricPolynomial< dim, Number1 > &bp)
BarycentricPolynomial< dim, Number1 > operator+(const Number2 &a, const BarycentricPolynomial< dim, Number1 > &bp)
SynchronousIterators< Iterators > operator++(SynchronousIterators< Iterators > &a)
bool operator<(const SynchronousIterators< Iterators > &a, const SynchronousIterators< Iterators > &b)