Reference documentation for deal.II version 9.5.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\}}\)
Loading...
Searching...
No Matches
index_set.h
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2009 - 2023 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_index_set_h
17#define dealii_index_set_h
18
19#include <deal.II/base/config.h>
20
23#include <deal.II/base/mutex.h>
24
25#include <boost/container/small_vector.hpp>
26
27#include <algorithm>
28#include <vector>
29
30
31#ifdef DEAL_II_WITH_TRILINOS
32# include <Epetra_Map.h>
33# ifdef DEAL_II_TRILINOS_WITH_TPETRA
34# include <Tpetra_Map.hpp>
35# endif
36#endif
37
38#ifdef DEAL_II_WITH_PETSC
39# include <petscis.h>
40#endif
41
43
67{
68public:
69 // forward declarations:
70 class ElementIterator;
71 class IntervalIterator;
72
78
96 using value_type = signed int;
97
98
102 IndexSet();
103
107 explicit IndexSet(const size_type size);
108
112 IndexSet(const IndexSet &) = default;
113
117 IndexSet &
118 operator=(const IndexSet &) = default;
119
124 IndexSet(IndexSet &&is) noexcept;
125
130 IndexSet &
131 operator=(IndexSet &&is) noexcept;
132
133#ifdef DEAL_II_WITH_TRILINOS
137 explicit IndexSet(const Epetra_BlockMap &map);
138#endif
139
144 void
145 clear();
146
153 void
154 set_size(const size_type size);
155
164 size() const;
165
172 void
173 add_range(const size_type begin, const size_type end);
174
183 void
184 add_index(const size_type index);
185
203 template <typename ForwardIterator>
204 void
205 add_indices(const ForwardIterator &begin, const ForwardIterator &end);
206
222 void
223 add_indices(const IndexSet &other, const size_type offset = 0);
224
228 bool
229 is_element(const size_type index) const;
230
235 bool
236 is_contiguous() const;
237
242 bool
243 is_empty() const;
244
254 bool
255 is_ascending_and_one_to_one(const MPI_Comm communicator) const;
256
261 n_elements() const;
262
269 nth_index_in_set(const size_type local_index) const;
270
278 index_within_set(const size_type global_index) const;
279
288 unsigned int
289 n_intervals() const;
290
304
309 void
310 compress() const;
311
317 bool
318 operator==(const IndexSet &is) const;
319
325 bool
326 operator!=(const IndexSet &is) const;
327
335 operator&(const IndexSet &is) const;
336
351 get_view(const size_type begin, const size_type end) const;
352
358 std::vector<IndexSet>
360 const std::vector<types::global_dof_index> &n_indices_per_block) const;
361
367 void
368 subtract_set(const IndexSet &other);
369
387 tensor_product(const IndexSet &other) const;
388
394 pop_back();
395
401 pop_front();
402
411 std::vector<size_type>
412 get_index_vector() const;
413
420 DEAL_II_DEPRECATED_EARLY
421 void
422 fill_index_vector(std::vector<size_type> &indices) const;
423
436 template <typename VectorType>
437 void
438 fill_binary_vector(VectorType &vector) const;
439
444 template <class StreamType>
445 void
446 print(StreamType &out) const;
447
452 void
453 write(std::ostream &out) const;
454
459 void
460 read(std::istream &in);
461
466 void
467 block_write(std::ostream &out) const;
468
473 void
474 block_read(std::istream &in);
475
476#ifdef DEAL_II_WITH_TRILINOS
505 Epetra_Map
506 make_trilinos_map(const MPI_Comm communicator = MPI_COMM_WORLD,
507 const bool overlapping = false) const;
508
509# ifdef DEAL_II_TRILINOS_WITH_TPETRA
510 Tpetra::Map<int, types::signed_global_dof_index>
511 make_tpetra_map(const MPI_Comm communicator = MPI_COMM_WORLD,
512 const bool overlapping = false) const;
513# endif
514#endif
515
516#ifdef DEAL_II_WITH_PETSC
517 IS
518 make_petsc_is(const MPI_Comm communicator = MPI_COMM_WORLD) const;
519#endif
520
521
526 std::size_t
527 memory_consumption() const;
528
530 size_type,
531 << "The global index " << arg1
532 << " is not an element of this set.");
533
539 template <class Archive>
540 void
541 serialize(Archive &ar, const unsigned int version);
542
543
555 {
556 public:
561 IntervalAccessor(const IndexSet *idxset, const size_type range_idx);
562
566 explicit IntervalAccessor(const IndexSet *idxset);
567
572 n_elements() const;
573
577 bool
578 is_valid() const;
579
584 begin() const;
585
591 end() const;
592
597 last() const;
598
599 private:
608 operator=(const IntervalAccessor &other);
609
613 bool
614 operator==(const IntervalAccessor &other) const;
618 bool
619 operator<(const IntervalAccessor &other) const;
624 void
625 advance();
630
636
637 friend class IntervalIterator;
638 };
639
645 {
646 public:
651 IntervalIterator(const IndexSet *idxset, const size_type range_idx);
652
656 explicit IntervalIterator(const IndexSet *idxset);
657
662
666 IntervalIterator(const IntervalIterator &other) = default;
667
672 operator=(const IntervalIterator &other) = default;
673
678 operator++();
679
684 operator++(int);
685
689 const IntervalAccessor &
690 operator*() const;
691
695 const IntervalAccessor *
696 operator->() const;
697
701 bool
702 operator==(const IntervalIterator &) const;
703
707 bool
708 operator!=(const IntervalIterator &) const;
709
713 bool
714 operator<(const IntervalIterator &) const;
715
722 int
723 operator-(const IntervalIterator &p) const;
724
730 using iterator_category = std::forward_iterator_tag;
732 using difference_type = std::ptrdiff_t;
735
736 private:
741 };
742
748 {
749 public:
754 ElementIterator(const IndexSet *idxset,
755 const size_type range_idx,
756 const size_type index);
757
761 explicit ElementIterator(const IndexSet *idxset);
762
768 operator*() const;
769
773 bool
774 is_valid() const;
775
780 operator++();
781
786 operator++(int);
787
791 bool
792 operator==(const ElementIterator &) const;
793
797 bool
798 operator!=(const ElementIterator &) const;
799
803 bool
804 operator<(const ElementIterator &) const;
805
814 std::ptrdiff_t
815 operator-(const ElementIterator &p) const;
816
822 using iterator_category = std::forward_iterator_tag;
824 using difference_type = std::ptrdiff_t;
827
828 private:
832 void
833 advance();
834
847 };
848
854 begin() const;
855
871 at(const size_type global_index) const;
872
878 end() const;
879
884 begin_intervals() const;
885
891 end_intervals() const;
892
897private:
906 struct Range
907 {
910
912
921 Range();
922
930 Range(const size_type i1, const size_type i2);
931
932 friend inline bool
933 operator<(const Range &range_1, const Range &range_2)
934 {
935 return (
936 (range_1.begin < range_2.begin) ||
937 ((range_1.begin == range_2.begin) && (range_1.end < range_2.end)));
938 }
939
940 static bool
942 {
943 return x.end < y.end;
944 }
945
946 static bool
948 {
949 return (x.nth_index_in_set + (x.end - x.begin) <
950 y.nth_index_in_set + (y.end - y.begin));
951 }
952
953 friend inline bool
954 operator==(const Range &range_1, const Range &range_2)
955 {
956 return ((range_1.begin == range_2.begin) && (range_1.end == range_2.end));
957 }
958
959 static std::size_t
961 {
962 return sizeof(Range);
963 }
964
970 template <class Archive>
971 void
972 serialize(Archive &ar, const unsigned int version);
973 };
974
983 mutable std::vector<Range> ranges;
984
993 mutable bool is_compressed;
994
1000
1011
1017
1021 void
1022 do_compress() const;
1023
1030 bool
1031 is_element_binary_search(const size_type local_index) const;
1032
1039 size_type
1040 nth_index_in_set_binary_search(const size_type local_index) const;
1041
1048 size_type
1049 index_within_set_binary_search(const size_type global_index) const;
1050
1056 void
1057 add_range_lower_bound(const Range &range);
1058
1062 void
1064 boost::container::small_vector<std::pair<size_type, size_type>, 200>
1065 & tmp_ranges,
1066 const bool ranges_are_sorted);
1067};
1068
1069
1070
1088inline IndexSet
1090{
1091 IndexSet is(N);
1092 is.add_range(0, N);
1093 is.compress();
1094 return is;
1095}
1096
1097/* ------------------ inline functions ------------------ */
1098
1099
1100/* IntervalAccessor */
1101
1103 const IndexSet * idxset,
1104 const IndexSet::size_type range_idx)
1105 : index_set(idxset)
1106 , range_idx(range_idx)
1107{
1108 Assert(range_idx < idxset->n_intervals(),
1109 ExcInternalError("Invalid range index"));
1110}
1111
1112
1113
1115 : index_set(idxset)
1116 , range_idx(numbers::invalid_dof_index)
1117{}
1118
1119
1120
1122 const IndexSet::IntervalAccessor &other)
1123 : index_set(other.index_set)
1124 , range_idx(other.range_idx)
1125{
1127 ExcMessage("invalid iterator"));
1128}
1129
1130
1131
1134{
1135 Assert(is_valid(), ExcMessage("invalid iterator"));
1136 return index_set->ranges[range_idx].end - index_set->ranges[range_idx].begin;
1137}
1138
1139
1140
1141inline bool
1143{
1144 return index_set != nullptr && range_idx < index_set->n_intervals();
1145}
1146
1147
1148
1151{
1152 Assert(is_valid(), ExcMessage("invalid iterator"));
1153 return {index_set, range_idx, index_set->ranges[range_idx].begin};
1154}
1155
1156
1157
1160{
1161 Assert(is_valid(), ExcMessage("invalid iterator"));
1162
1163 // point to first index in next interval unless we are the last interval.
1164 if (range_idx < index_set->ranges.size() - 1)
1165 return {index_set, range_idx + 1, index_set->ranges[range_idx + 1].begin};
1166 else
1167 return index_set->end();
1168}
1169
1170
1171
1174{
1175 Assert(is_valid(), ExcMessage("invalid iterator"));
1176
1177 return index_set->ranges[range_idx].end - 1;
1178}
1179
1180
1181
1184{
1185 index_set = other.index_set;
1186 range_idx = other.range_idx;
1187 Assert(range_idx == numbers::invalid_dof_index || is_valid(),
1188 ExcMessage("invalid iterator"));
1189 return *this;
1190}
1191
1192
1193
1194inline bool
1196 const IndexSet::IntervalAccessor &other) const
1197{
1198 Assert(index_set == other.index_set,
1199 ExcMessage(
1200 "Can not compare accessors pointing to different IndexSets"));
1201 return range_idx == other.range_idx;
1202}
1203
1204
1205
1206inline bool
1208 const IndexSet::IntervalAccessor &other) const
1209{
1210 Assert(index_set == other.index_set,
1211 ExcMessage(
1212 "Can not compare accessors pointing to different IndexSets"));
1213 return range_idx < other.range_idx;
1214}
1215
1216
1217
1218inline void
1220{
1221 Assert(
1222 is_valid(),
1223 ExcMessage(
1224 "Impossible to advance an IndexSet::IntervalIterator that is invalid"));
1225 ++range_idx;
1226
1227 // set ourselves to invalid if we walk off the end
1228 if (range_idx >= index_set->ranges.size())
1229 range_idx = numbers::invalid_dof_index;
1230}
1231
1232
1233/* IntervalIterator */
1234
1236 const IndexSet * idxset,
1237 const IndexSet::size_type range_idx)
1238 : accessor(idxset, range_idx)
1239{}
1240
1241
1242
1244 : accessor(nullptr)
1245{}
1246
1247
1248
1250 : accessor(idxset)
1251{}
1252
1253
1254
1257{
1258 accessor.advance();
1259 return *this;
1260}
1261
1262
1263
1266{
1267 const IndexSet::IntervalIterator iter = *this;
1268 accessor.advance();
1269 return iter;
1270}
1271
1272
1273
1274inline const IndexSet::IntervalAccessor &
1276{
1277 return accessor;
1278}
1279
1280
1281
1282inline const IndexSet::IntervalAccessor *
1284{
1285 return &accessor;
1286}
1287
1288
1289
1290inline bool
1292 const IndexSet::IntervalIterator &other) const
1293{
1294 return accessor == other.accessor;
1295}
1296
1297
1298
1299inline bool
1301 const IndexSet::IntervalIterator &other) const
1302{
1303 return !(*this == other);
1304}
1305
1306
1307
1308inline bool
1310 const IndexSet::IntervalIterator &other) const
1311{
1312 return accessor < other.accessor;
1313}
1314
1315
1316
1317inline int
1319 const IndexSet::IntervalIterator &other) const
1320{
1321 Assert(accessor.index_set == other.accessor.index_set,
1322 ExcMessage(
1323 "Can not compare iterators belonging to different IndexSets"));
1324
1325 const size_type lhs = (accessor.range_idx == numbers::invalid_dof_index) ?
1326 accessor.index_set->ranges.size() :
1327 accessor.range_idx;
1328 const size_type rhs =
1330 accessor.index_set->ranges.size() :
1331 other.accessor.range_idx;
1332
1333 if (lhs > rhs)
1334 return static_cast<int>(lhs - rhs);
1335 else
1336 return -static_cast<int>(rhs - lhs);
1337}
1338
1339
1340
1341/* ElementIterator */
1342
1344 const IndexSet * idxset,
1345 const IndexSet::size_type range_idx,
1346 const IndexSet::size_type index)
1347 : index_set(idxset)
1348 , range_idx(range_idx)
1349 , idx(index)
1350{
1351 Assert(range_idx < index_set->ranges.size(),
1352 ExcMessage(
1353 "Invalid range index for IndexSet::ElementIterator constructor."));
1354 Assert(
1355 idx >= index_set->ranges[range_idx].begin &&
1356 idx < index_set->ranges[range_idx].end,
1358 "Invalid index argument for IndexSet::ElementIterator constructor."));
1359}
1360
1361
1362
1364 : index_set(idxset)
1365 , range_idx(numbers::invalid_dof_index)
1366 , idx(numbers::invalid_dof_index)
1367{}
1368
1369
1370
1371inline bool
1373{
1374 Assert((range_idx == numbers::invalid_dof_index &&
1376 (range_idx < index_set->ranges.size() &&
1377 idx < index_set->ranges[range_idx].end),
1378 ExcInternalError("Invalid ElementIterator state."));
1379
1380 return (range_idx < index_set->ranges.size() &&
1381 idx < index_set->ranges[range_idx].end);
1382}
1383
1384
1385
1388{
1389 Assert(
1390 is_valid(),
1391 ExcMessage(
1392 "Impossible to dereference an IndexSet::ElementIterator that is invalid"));
1393 return idx;
1394}
1395
1396
1397
1398inline bool
1400 const IndexSet::ElementIterator &other) const
1401{
1402 Assert(index_set == other.index_set,
1403 ExcMessage(
1404 "Can not compare iterators belonging to different IndexSets"));
1405 return range_idx == other.range_idx && idx == other.idx;
1406}
1407
1408
1409
1410inline void
1412{
1413 Assert(
1414 is_valid(),
1415 ExcMessage(
1416 "Impossible to advance an IndexSet::ElementIterator that is invalid"));
1417 if (idx < index_set->ranges[range_idx].end)
1418 ++idx;
1419 // end of this range?
1420 if (idx == index_set->ranges[range_idx].end)
1421 {
1422 // point to first element in next interval if possible
1423 if (range_idx < index_set->ranges.size() - 1)
1424 {
1425 ++range_idx;
1426 idx = index_set->ranges[range_idx].begin;
1427 }
1428 else
1429 {
1430 // we just fell off the end, set to invalid:
1431 range_idx = numbers::invalid_dof_index;
1433 }
1434 }
1435}
1436
1437
1438
1441{
1442 advance();
1443 return *this;
1444}
1445
1446
1447
1450{
1451 const IndexSet::ElementIterator it = *this;
1452 advance();
1453 return it;
1454}
1455
1456
1457
1458inline bool
1460 const IndexSet::ElementIterator &other) const
1461{
1462 return !(*this == other);
1463}
1464
1465
1466
1467inline bool
1469 const IndexSet::ElementIterator &other) const
1470{
1471 Assert(index_set == other.index_set,
1472 ExcMessage(
1473 "Can not compare iterators belonging to different IndexSets"));
1474 return range_idx < other.range_idx ||
1475 (range_idx == other.range_idx && idx < other.idx);
1476}
1477
1478
1479
1480inline std::ptrdiff_t
1482 const IndexSet::ElementIterator &other) const
1483{
1484 Assert(index_set == other.index_set,
1485 ExcMessage(
1486 "Can not compare iterators belonging to different IndexSets"));
1487 if (*this == other)
1488 return 0;
1489 if (!(*this < other))
1490 return -(other - *this);
1491
1492 // only other can be equal to end() because of the checks above.
1493 Assert(is_valid(), ExcInternalError());
1494
1495 // Note: we now compute how far advance *this in "*this < other" to get other,
1496 // so we need to return -c at the end.
1497
1498 // first finish the current range:
1499 std::ptrdiff_t c = index_set->ranges[range_idx].end - idx;
1500
1501 // now walk in steps of ranges (need to start one behind our current one):
1502 for (size_type range = range_idx + 1;
1503 range < index_set->ranges.size() && range <= other.range_idx;
1504 ++range)
1505 c += index_set->ranges[range].end - index_set->ranges[range].begin;
1506
1507 Assert(
1508 other.range_idx < index_set->ranges.size() ||
1510 ExcMessage(
1511 "Inconsistent iterator state. Did you invalidate iterators by modifying the IndexSet?"));
1512
1513 // We might have walked too far because we went until the end of
1514 // other.range_idx, so walk backwards to other.idx:
1516 c -= index_set->ranges[other.range_idx].end - other.idx;
1517
1518 return -c;
1519}
1520
1521
1522/* Range */
1523
1525 : begin(numbers::invalid_dof_index)
1526 , end(numbers::invalid_dof_index)
1527 , nth_index_in_set(numbers::invalid_dof_index)
1528{}
1529
1530
1531
1533 : begin(i1)
1534 , end(i2)
1535 , nth_index_in_set(numbers::invalid_dof_index)
1536{}
1537
1538
1539
1540/* IndexSet itself */
1541
1543 : is_compressed(true)
1544 , index_space_size(0)
1545 , largest_range(numbers::invalid_unsigned_int)
1546{}
1547
1548
1549
1551 : is_compressed(true)
1552 , index_space_size(size)
1553 , largest_range(numbers::invalid_unsigned_int)
1554{}
1555
1556
1557
1558inline IndexSet::IndexSet(IndexSet &&is) noexcept
1559 : ranges(std::move(is.ranges))
1560 , is_compressed(is.is_compressed)
1561 , index_space_size(is.index_space_size)
1562 , largest_range(is.largest_range)
1563{
1564 is.ranges.clear();
1565 is.is_compressed = true;
1566 is.index_space_size = 0;
1567 is.largest_range = numbers::invalid_unsigned_int;
1568
1569 compress();
1570}
1571
1572
1573
1574inline IndexSet &
1576{
1577 ranges = std::move(is.ranges);
1578 is_compressed = is.is_compressed;
1579 index_space_size = is.index_space_size;
1580 largest_range = is.largest_range;
1581
1582 is.ranges.clear();
1583 is.is_compressed = true;
1584 is.index_space_size = 0;
1585 is.largest_range = numbers::invalid_unsigned_int;
1586
1587 compress();
1588
1589 return *this;
1590}
1591
1592
1593
1596{
1597 compress();
1598 if (ranges.size() > 0)
1599 return {this, 0, ranges[0].begin};
1600 else
1601 return end();
1602}
1603
1604
1605
1608{
1609 compress();
1610 return IndexSet::ElementIterator(this);
1611}
1612
1613
1614
1617{
1618 compress();
1619 if (ranges.size() > 0)
1620 return IndexSet::IntervalIterator(this, 0);
1621 else
1622 return end_intervals();
1623}
1624
1625
1626
1629{
1630 compress();
1631 return IndexSet::IntervalIterator(this);
1632}
1633
1634
1635
1636inline void
1638{
1639 // reset so that there are no indices in the set any more; however,
1640 // as documented, the index set retains its size
1641 ranges.clear();
1642 is_compressed = true;
1644}
1645
1646
1647
1648inline void
1650{
1651 Assert(ranges.empty(),
1652 ExcMessage("This function can only be called if the current "
1653 "object does not yet contain any elements."));
1654 index_space_size = sz;
1655 is_compressed = true;
1656}
1657
1658
1659
1662{
1663 return index_space_size;
1664}
1665
1666
1667
1668inline void
1670{
1671 if (is_compressed == true)
1672 return;
1673
1674 do_compress();
1675}
1676
1677
1678
1679inline void
1681{
1683
1684 const Range new_range(index, index + 1);
1685 if (ranges.size() == 0 || index > ranges.back().end)
1686 ranges.push_back(new_range);
1687 else if (index == ranges.back().end)
1688 ranges.back().end++;
1689 else
1690 add_range_lower_bound(new_range);
1691 is_compressed = false;
1692}
1693
1694
1695
1696inline void
1698{
1701 ExcIndexRangeType<size_type>(begin, 0, index_space_size));
1703 ExcIndexRangeType<size_type>(end, 0, index_space_size + 1));
1705
1706 if (begin != end)
1707 {
1708 const Range new_range(begin, end);
1709
1710 // the new index might be larger than the last index present in the
1711 // ranges. Then we can skip the binary search
1712 if (ranges.size() == 0 || begin > ranges.back().end)
1713 ranges.push_back(new_range);
1714 else if (begin == ranges.back().end)
1715 ranges.back().end = end;
1716 else
1717 add_range_lower_bound(new_range);
1718
1719 is_compressed = false;
1720 }
1721}
1722
1723
1724
1725template <typename ForwardIterator>
1726inline void
1727IndexSet::add_indices(const ForwardIterator &begin, const ForwardIterator &end)
1728{
1729 if (begin == end)
1730 return;
1731
1732 // identify ranges in the given iterator range by checking whether some
1733 // indices happen to be consecutive. to avoid quadratic complexity when
1734 // calling add_range many times (as add_range() going into the middle of an
1735 // already existing range must shift entries around), we first collect a
1736 // vector of ranges.
1737 boost::container::small_vector<std::pair<size_type, size_type>, 200>
1738 tmp_ranges;
1739 bool ranges_are_sorted = true;
1740 for (ForwardIterator p = begin; p != end;)
1741 {
1742 // Starting with the current iterator 'p', find an iterator
1743 // 'q' so that the indices pointed to by the iterators in
1744 // the range [p,q) are consecutive. These indices then form
1745 // a range that is contiguous, and that can be added all
1746 // at once.
1747 const size_type begin_index = *p;
1748 size_type end_index = begin_index + 1;
1749 ForwardIterator q = p;
1750 ++q;
1751 while ((q != end) && (static_cast<size_type>(*q) == end_index))
1752 {
1753 ++end_index;
1754 ++q;
1755 }
1756
1757 // Add this range:
1758 tmp_ranges.emplace_back(begin_index, end_index);
1759
1760 // Then move on to the next element in the input range.
1761 // If the starting index of the next go-around of the for loop is less
1762 // than the end index of the one just identified, then we will have at
1763 // least one pair of ranges that are not sorted, and consequently the
1764 // whole collection of ranges is not sorted.
1765 p = q;
1766 if (p != end && static_cast<size_type>(*p) < end_index)
1767 ranges_are_sorted = false;
1768 }
1769
1770 add_ranges_internal(tmp_ranges, ranges_are_sorted);
1771}
1772
1773
1774
1775inline bool
1777{
1778 if (ranges.empty() == false)
1779 {
1780 compress();
1781
1782 // fast check whether the index is in the largest range
1784 if (index >= ranges[largest_range].begin &&
1785 index < ranges[largest_range].end)
1786 return true;
1787 else if (ranges.size() > 1)
1788 return is_element_binary_search(index);
1789 else
1790 return false;
1791 }
1792 else
1793 return false;
1794}
1795
1796
1797
1798inline bool
1800{
1801 compress();
1802 return (ranges.size() <= 1);
1803}
1804
1805
1806
1807inline bool
1809{
1810 return ranges.empty();
1811}
1812
1813
1814
1817{
1818 // make sure we have non-overlapping ranges
1819 compress();
1820
1821 size_type v = 0;
1822 if (!ranges.empty())
1823 {
1824 Range &r = ranges.back();
1825 v = r.nth_index_in_set + r.end - r.begin;
1826 }
1827
1828#ifdef DEBUG
1829 size_type s = 0;
1830 for (const auto &range : ranges)
1831 s += (range.end - range.begin);
1832 Assert(s == v, ExcInternalError());
1833#endif
1834
1835 return v;
1836}
1837
1838
1839
1840inline unsigned int
1842{
1843 compress();
1844 return ranges.size();
1845}
1846
1847
1848
1851{
1852 Assert(ranges.empty() == false, ExcMessage("IndexSet cannot be empty."));
1853
1854 compress();
1855 const std::vector<Range>::const_iterator main_range =
1856 ranges.begin() + largest_range;
1857
1858 return main_range->nth_index_in_set;
1859}
1860
1861
1862
1865{
1867
1868 compress();
1869
1870 // first check whether the index is in the largest range
1872 const auto main_range = ranges.begin() + largest_range;
1873 if (n >= main_range->nth_index_in_set &&
1874 n < main_range->nth_index_in_set + (main_range->end - main_range->begin))
1875 return main_range->begin + (n - main_range->nth_index_in_set);
1876 else
1878}
1879
1880
1881
1884{
1885 // to make this call thread-safe, compress() must not be called through this
1886 // function
1887 Assert(is_compressed == true, ExcMessage("IndexSet must be compressed."));
1888 AssertIndexRange(n, size());
1889
1890 // return immediately if the index set is empty
1891 if (is_empty())
1893
1894 // check whether the index is in the largest range. use the result to
1895 // perform a one-sided binary search afterward
1898 return (n - ranges[largest_range].begin) +
1899 ranges[largest_range].nth_index_in_set;
1900 else if (ranges.size() > 1)
1902 else
1904}
1905
1906
1907
1908inline bool
1910{
1911 Assert(size() == is.size(), ExcDimensionMismatch(size(), is.size()));
1912
1913 compress();
1914 is.compress();
1915
1916 return ranges == is.ranges;
1917}
1918
1919
1920
1921inline bool
1923{
1924 Assert(size() == is.size(), ExcDimensionMismatch(size(), is.size()));
1925
1926 compress();
1927 is.compress();
1928
1929 return ranges != is.ranges;
1930}
1931
1932
1933
1934template <typename Vector>
1935void
1937{
1938 Assert(vector.size() == size(), ExcDimensionMismatch(vector.size(), size()));
1939
1940 compress();
1941 // first fill all elements of the vector with zeroes.
1942 std::fill(vector.begin(), vector.end(), 0);
1943
1944 // then write ones into the elements whose indices are contained in the
1945 // index set
1946 for (const auto &range : ranges)
1947 for (size_type i = range.begin; i < range.end; ++i)
1948 vector[i] = 1;
1949}
1950
1951
1952
1953template <class StreamType>
1954inline void
1955IndexSet::print(StreamType &out) const
1956{
1957 compress();
1958 out << '{';
1959 std::vector<Range>::const_iterator p;
1960 for (p = ranges.begin(); p != ranges.end(); ++p)
1961 {
1962 if (p->end - p->begin == 1)
1963 out << p->begin;
1964 else
1965 out << '[' << p->begin << ',' << p->end - 1 << ']';
1966
1967 if (p != --ranges.end())
1968 out << ", ";
1969 }
1970 out << '}' << std::endl;
1971}
1972
1973
1974
1975template <class Archive>
1976inline void
1977IndexSet::Range::serialize(Archive &ar, const unsigned int)
1978{
1980}
1981
1982
1983
1984template <class Archive>
1985inline void
1986IndexSet::serialize(Archive &ar, const unsigned int)
1987{
1989}
1990
1992
1993#endif
size_type operator*() const
Definition index_set.h:1387
bool operator==(const ElementIterator &) const
Definition index_set.h:1399
ElementIterator(const IndexSet *idxset, const size_type range_idx, const size_type index)
Definition index_set.h:1343
std::ptrdiff_t difference_type
Definition index_set.h:824
ElementIterator & operator++()
Definition index_set.h:1440
bool operator<(const ElementIterator &) const
Definition index_set.h:1468
bool operator!=(const ElementIterator &) const
Definition index_set.h:1459
std::ptrdiff_t operator-(const ElementIterator &p) const
Definition index_set.h:1481
std::forward_iterator_tag iterator_category
Definition index_set.h:822
const IndexSet * index_set
Definition index_set.h:838
size_type last() const
Definition index_set.h:1173
ElementIterator end() const
Definition index_set.h:1159
size_type n_elements() const
Definition index_set.h:1133
bool operator==(const IntervalAccessor &other) const
Definition index_set.h:1195
IntervalAccessor & operator=(const IntervalAccessor &other)
Definition index_set.h:1183
bool operator<(const IntervalAccessor &other) const
Definition index_set.h:1207
IntervalAccessor(const IndexSet *idxset, const size_type range_idx)
Definition index_set.h:1102
const IndexSet * index_set
Definition index_set.h:629
ElementIterator begin() const
Definition index_set.h:1150
std::ptrdiff_t difference_type
Definition index_set.h:732
const IntervalAccessor & operator*() const
Definition index_set.h:1275
IntervalIterator(const IntervalIterator &other)=default
bool operator==(const IntervalIterator &) const
Definition index_set.h:1291
int operator-(const IntervalIterator &p) const
Definition index_set.h:1318
bool operator<(const IntervalIterator &) const
Definition index_set.h:1309
std::forward_iterator_tag iterator_category
Definition index_set.h:730
IntervalIterator & operator=(const IntervalIterator &other)=default
IntervalIterator & operator++()
Definition index_set.h:1256
const IntervalAccessor * operator->() const
Definition index_set.h:1283
IntervalAccessor accessor
Definition index_set.h:740
bool operator!=(const IntervalIterator &) const
Definition index_set.h:1300
bool is_element_binary_search(const size_type local_index) const
Definition index_set.cc:565
size_type index_within_set_binary_search(const size_type global_index) const
Definition index_set.cc:618
size_type largest_range
Definition index_set.h:1010
IS make_petsc_is(const MPI_Comm communicator=MPI_COMM_WORLD) const
Definition index_set.cc:858
bool is_ascending_and_one_to_one(const MPI_Comm communicator) const
Definition index_set.cc:878
bool is_contiguous() const
Definition index_set.h:1799
unsigned int n_intervals() const
Definition index_set.h:1841
IndexSet(const IndexSet &)=default
IntervalIterator end_intervals() const
Definition index_set.h:1628
void do_compress() const
Definition index_set.cc:100
ElementIterator at(const size_type global_index) const
Definition index_set.cc:640
size_type size() const
Definition index_set.h:1661
void fill_index_vector(std::vector< size_type > &indices) const
Definition index_set.cc:709
bool is_empty() const
Definition index_set.h:1808
size_type index_within_set(const size_type global_index) const
Definition index_set.h:1883
std::vector< IndexSet > split_by_block(const std::vector< types::global_dof_index > &n_indices_per_block) const
Definition index_set.cc:241
size_type n_elements() const
Definition index_set.h:1816
bool operator==(const IndexSet &is) const
Definition index_set.h:1909
void add_range_lower_bound(const Range &range)
Definition index_set.cc:389
bool is_element(const size_type index) const
Definition index_set.h:1776
void serialize(Archive &ar, const unsigned int version)
Definition index_set.h:1986
ElementIterator begin() const
Definition index_set.h:1595
size_type pop_front()
Definition index_set.cc:367
signed int value_type
Definition index_set.h:96
void set_size(const size_type size)
Definition index_set.h:1649
bool operator!=(const IndexSet &is) const
Definition index_set.h:1922
void read(std::istream &in)
Definition index_set.cc:509
Epetra_Map make_trilinos_map(const MPI_Comm communicator=MPI_COMM_WORLD, const bool overlapping=false) const
Definition index_set.cc:789
void clear()
Definition index_set.h:1637
size_type largest_range_starting_index() const
Definition index_set.h:1850
IndexSet tensor_product(const IndexSet &other) const
Definition index_set.cc:337
void add_index(const size_type index)
Definition index_set.h:1680
void write(std::ostream &out) const
Definition index_set.cc:494
void block_read(std::istream &in)
Definition index_set.cc:545
bool is_compressed
Definition index_set.h:993
void add_ranges_internal(boost::container::small_vector< std::pair< size_type, size_type >, 200 > &tmp_ranges, const bool ranges_are_sorted)
Definition index_set.cc:398
IntervalIterator begin_intervals() const
Definition index_set.h:1616
Tpetra::Map< int, types::signed_global_dof_index > make_tpetra_map(const MPI_Comm communicator=MPI_COMM_WORLD, const bool overlapping=false) const
Definition index_set.cc:720
IndexSet & operator=(const IndexSet &)=default
void fill_binary_vector(VectorType &vector) const
std::vector< Range > ranges
Definition index_set.h:983
void subtract_set(const IndexSet &other)
Definition index_set.cc:268
ElementIterator end() const
Definition index_set.h:1607
Threads::Mutex compress_mutex
Definition index_set.h:1016
IndexSet complete_index_set(const IndexSet::size_type N)
Definition index_set.h:1089
size_type index_space_size
Definition index_set.h:999
void block_write(std::ostream &out) const
Definition index_set.cc:531
IndexSet get_view(const size_type begin, const size_type end) const
Definition index_set.cc:214
void add_range(const size_type begin, const size_type end)
Definition index_set.h:1697
size_type nth_index_in_set(const size_type local_index) const
Definition index_set.h:1864
std::size_t memory_consumption() const
Definition index_set.cc:944
void print(StreamType &out) const
Definition index_set.h:1955
size_type nth_index_in_set_binary_search(const size_type local_index) const
Definition index_set.cc:601
void compress() const
Definition index_set.h:1669
size_type pop_back()
Definition index_set.cc:349
std::vector< size_type > get_index_vector() const
Definition index_set.cc:690
types::global_dof_index size_type
Definition index_set.h:77
void add_indices(const ForwardIterator &begin, const ForwardIterator &end)
Definition index_set.h:1727
IndexSet operator&(const IndexSet &is) const
iterator end()
size_type size() const
iterator begin()
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:472
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
#define Assert(cond, exc)
static ::ExceptionBase & ExcIndexNotPresent(size_type arg1)
#define AssertIndexRange(index, range)
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
#define DeclException1(Exception1, type1, outsequence)
Definition exceptions.h:510
static ::ExceptionBase & ExcMessage(std::string arg1)
static const unsigned int invalid_unsigned_int
Definition types.h:213
const types::global_dof_index invalid_dof_index
Definition types.h:233
unsigned int global_dof_index
Definition types.h:82
static bool end_compare(const IndexSet::Range &x, const IndexSet::Range &y)
Definition index_set.h:941
static bool nth_index_compare(const IndexSet::Range &x, const IndexSet::Range &y)
Definition index_set.h:947
friend bool operator==(const Range &range_1, const Range &range_2)
Definition index_set.h:954
size_type end
Definition index_set.h:909
size_type nth_index_in_set
Definition index_set.h:911
static std::size_t memory_consumption()
Definition index_set.h:960
size_type begin
Definition index_set.h:908
void serialize(Archive &ar, const unsigned int version)
Definition index_set.h:1977
friend bool operator<(const Range &range_1, const Range &range_2)
Definition index_set.h:933
void advance(std::tuple< I1, I2 > &t, const unsigned int n)