Line data Source code
1 1 : // Distributed under the MIT License.
2 : // See LICENSE.txt for details.
3 :
4 : /// \file
5 : /// Defines class SegmentId.
6 :
7 : #pragma once
8 :
9 : #include <cstddef>
10 : #include <iosfwd>
11 : #include <pup.h>
12 :
13 : #include "Domain/Structure/Side.hpp"
14 : #include "Utilities/ConstantExpressions.hpp"
15 : #include "Utilities/ErrorHandling/Assert.hpp"
16 :
17 : /*!
18 : * \ingroup ComputationalDomainGroup
19 : * \brief A SegmentId labels a segment of the interval \f$[-1,1]\f$ and is used
20 : * to identify the bounds of an Element in a Block in each dimension.
21 : *
22 : * In \f$d\f$ dimensions, \f$d\f$ SegmentId%s are used to identify an Element.
23 : * In each dimension, a segment spans the subinterval \f$[-1 + 2 \frac{i}{N}, -1
24 : * + 2 \frac{i+1}{N}]\f$ of the logical coordinates of a Block, where \f$i \f$=
25 : * `index` and \f$N = 2^L\f$ where \f$L\f$ = `refinement_level`.
26 : *
27 : * \image html SegmentId.svg "SegmentIds"
28 : *
29 : * In the figure, The `index` of segments increase from the `lower side` to the
30 : * `upper side` in each dimension of a Block, while the `refinement level`
31 : * increases as the segments are subdivided. For example, let the segment
32 : * labeled `self` be on `refinement level` \f$L\f$, with `index` \f$i\f$. Its
33 : * `parent` segment is on `refinement level` \f$L-1\f$ with `index`
34 : * \f$\frac{i-1}{2}\f$. The `children` of `self` are on `refinement level`
35 : * \f$L+1\f$, and have `index` \f$2i\f$ and \f$2i+1\f$ for the lower and upper
36 : * child respectively. Also labeled on the figure are the `sibling` and
37 : * `abutting nibling` (child of sibling) of `self`. These relationships between
38 : * segments are important for h-refinement, since in each dimension an Element
39 : * can be flagged to split into its two `children` segments, or join with its
40 : * `sibling` segment to form its `parent` segment. As refinement levels of
41 : * neighboring elements are kept within one, in the direction of its `sibling`,
42 : * a segment can only abut its `sibling` or `abutting nibling`, while on the
43 : * opposite side, it can abut a segment on its level, the next-lower, or the
44 : * next-higher level.
45 : */
46 1 : class SegmentId {
47 : public:
48 : /// Default constructor needed for Charm++ serialization.
49 1 : SegmentId() = default;
50 0 : SegmentId(const SegmentId& segment_id) = default;
51 0 : SegmentId(SegmentId&& segment_id) = default;
52 0 : ~SegmentId() = default;
53 0 : SegmentId& operator=(const SegmentId& segment_id) = default;
54 0 : SegmentId& operator=(SegmentId&& segment_id) = default;
55 :
56 0 : SegmentId(size_t refinement_level, size_t index);
57 :
58 0 : constexpr size_t refinement_level() const { return refinement_level_; }
59 :
60 0 : constexpr size_t index() const { return index_; }
61 :
62 0 : SegmentId id_of_parent() const;
63 :
64 0 : SegmentId id_of_child(Side side) const;
65 :
66 : /// The other child of the parent of this segment
67 1 : SegmentId id_of_sibling() const;
68 :
69 : /// The child of the sibling of this segment that shares an endpoint with it
70 1 : SegmentId id_of_abutting_nibling() const;
71 :
72 : /// The side on which this segment shares an endpoint with its sibling
73 1 : Side side_of_sibling() const;
74 :
75 : /// The id this segment would have if the coordinate axis were flipped.
76 1 : SegmentId id_if_flipped() const;
77 :
78 : /// The block logical coordinate of the endpoint of the segment on the given
79 : /// Side.
80 1 : double endpoint(Side side) const;
81 :
82 : /// The block logical coordinate of the midpoint of the segment
83 1 : double midpoint() const {
84 : return -1.0 + (1.0 + 2.0 * index_) / two_to_the(refinement_level_);
85 : }
86 :
87 : // NOLINTNEXTLINE
88 0 : void pup(PUP::er& p);
89 :
90 : private:
91 0 : size_t refinement_level_;
92 0 : size_t index_;
93 : };
94 :
95 : /// Output operator for SegmentId.
96 1 : std::ostream& operator<<(std::ostream& os, const SegmentId& id);
97 :
98 : /// Equivalence operator for SegmentId.
99 : bool operator==(const SegmentId& lhs, const SegmentId& rhs);
100 :
101 : /// Inequivalence operator for SegmentId.
102 : bool operator!=(const SegmentId& lhs, const SegmentId& rhs);
103 :
104 : /// \ingroup ComputationalDomainGroup
105 : /// Check if two segments overlap.
106 1 : bool overlapping(const SegmentId& a, const SegmentId& b);
107 :
108 : // #############################################################################
109 : // INLINE DEFINITIONS
110 : // #############################################################################
111 :
112 : inline SegmentId SegmentId::id_of_parent() const {
113 : ASSERT(0 != refinement_level_,
114 : "Cannot call id_of_parent() on root refinement level!");
115 : // The parent has half as many segments as the child.
116 : return {refinement_level() - 1, index() / 2};
117 : }
118 :
119 : inline SegmentId SegmentId::id_of_child(Side side) const {
120 : // We cannot ASSERT on the maximum level because it's only known at runtime
121 : // and only known elsewhere in the code, not by SegmentId. I.e. SegmentId is
122 : // too low-level to know about this.
123 : // The child has twice as many segments as the parent, so for a particular
124 : // parent segment, there is both an upper and lower child segment.
125 : if (Side::Lower == side) {
126 : return {refinement_level() + 1, index() * 2};
127 : }
128 : return {refinement_level() + 1, 1 + index() * 2};
129 : }
130 :
131 : inline SegmentId SegmentId::id_of_sibling() const {
132 : ASSERT(0 != refinement_level(),
133 : "The segment on the root refinement level has no sibling");
134 : return {refinement_level(), (0 == index() % 2 ? index() + 1 : index() - 1)};
135 : }
136 :
137 : inline SegmentId SegmentId::id_of_abutting_nibling() const {
138 : ASSERT(0 != refinement_level(),
139 : "The segment on the root refinement level has no abutting nibling");
140 : return {refinement_level() + 1,
141 : (0 == index() % 2 ? 2 * index() + 2 : 2 * index() - 1)};
142 : }
143 :
144 : inline Side SegmentId::side_of_sibling() const {
145 : ASSERT(0 != refinement_level(),
146 : "The segment on the root refinement level has no sibling");
147 : return 0 == index() % 2 ? Side::Upper : Side::Lower;
148 : }
149 :
150 : inline SegmentId SegmentId::id_if_flipped() const {
151 : return {refinement_level(), two_to_the(refinement_level()) - 1 - index()};
152 : }
153 :
154 : inline double SegmentId::endpoint(Side side) const {
155 : if (Side::Lower == side) {
156 : return -1.0 + (2.0 * index()) / two_to_the(refinement_level());
157 : }
158 : return -1.0 + (2.0 * index() + 2.0) / two_to_the(refinement_level());
159 : }
160 :
161 : // These are defined so that a SegmentId can be used as part of a key of an
162 : // unordered_set or unordered_map.
163 :
164 : // hash_value is called by boost::hash and related functions.
165 0 : size_t hash_value(const SegmentId& s);
166 :
167 : namespace std {
168 : template <>
169 : struct hash<SegmentId> {
170 : size_t operator()(const SegmentId& segment_id) const;
171 : };
172 : } // namespace std
173 :
174 1 : inline bool operator==(const SegmentId& lhs, const SegmentId& rhs) {
175 : return (lhs.refinement_level() == rhs.refinement_level() and
176 : lhs.index() == rhs.index());
177 : }
178 :
179 1 : inline bool operator!=(const SegmentId& lhs, const SegmentId& rhs) {
180 : return not(lhs == rhs);
181 : }
|
