Line data Source code
1 0 : // Distributed under the MIT License.
2 : // See LICENSE.txt for details.
3 :
4 : #pragma once
5 :
6 : #include <array>
7 : #include <cstddef>
8 : #include <limits>
9 : #include <ostream>
10 :
11 : #include "Domain/Structure/Side.hpp"
12 : #include "Utilities/ConstantExpressions.hpp"
13 : #include "Utilities/Gsl.hpp"
14 : #include "Utilities/Requires.hpp"
15 : #include "Utilities/StdHelpers.hpp"
16 :
17 : /*!
18 : * \brief An element of dimension `ElementDim` on the boundary of a hypercube of
19 : * dimension `HypercubeDim`
20 : *
21 : * A hypercube of dimension \f$n\f$ (`HypercubeDim`) is composed of
22 : * \f$2^{n-k}\binom{n}{k}\f$ elements of dimension \f$k \leq n\f$
23 : * (`ElementDim`). For example, a 3D cube has 8 vertices (\f$k=0\f$), 12 edges
24 : * (\f$k=1\f$), 6 faces (\f$k=2\f$) and 1 cell (\f$k=3\f$). Each element is
25 : * identified by the \f$k\f$ dimensions it shares with the parent hypercube
26 : * and \f$n - k\f$ indices that specify whether it is located on the lower or
27 : * upper side of the parent hypercube's remaining dimensions.
28 : */
29 : template <size_t ElementDim, size_t HypercubeDim>
30 1 : struct HypercubeElement {
31 : static_assert(ElementDim <= HypercubeDim);
32 :
33 0 : HypercubeElement(std::array<size_t, ElementDim> dimensions_in_parent,
34 : std::array<Side, HypercubeDim - ElementDim> index);
35 :
36 : template <size_t NumIndices = HypercubeDim - ElementDim,
37 : Requires<NumIndices == 0> = nullptr>
38 0 : HypercubeElement() {
39 : for (size_t d = 0; d < ElementDim; ++d) {
40 : gsl::at(dimensions_in_parent_, d) = d;
41 : }
42 : }
43 :
44 : template <size_t LocalElementDim = ElementDim,
45 : Requires<LocalElementDim == 0> = nullptr>
46 0 : explicit HypercubeElement(std::array<Side, HypercubeDim> index);
47 :
48 : template <typename... Indices, size_t LocalElementDim = ElementDim,
49 : size_t LocalHypercubeDim = HypercubeDim,
50 : Requires<(LocalElementDim == 0 and LocalHypercubeDim > 0 and
51 : sizeof...(Indices) == LocalHypercubeDim)> = nullptr>
52 0 : explicit HypercubeElement(Indices... indices)
53 : : index_{{static_cast<Side>(indices)...}} {}
54 :
55 : template <size_t LocalElementDim = ElementDim,
56 : Requires<LocalElementDim == 1> = nullptr>
57 0 : HypercubeElement(size_t dim_in_parent,
58 : std::array<Side, HypercubeDim - 1> index);
59 :
60 : /// @{
61 : /// The parent hypercube's dimensions that this element shares
62 1 : const std::array<size_t, ElementDim>& dimensions_in_parent() const;
63 :
64 : template <size_t LocalElementDim = ElementDim,
65 : Requires<LocalElementDim == 1> = nullptr>
66 1 : size_t dimension_in_parent() const;
67 : /// @}
68 :
69 : /// @{
70 : /// Whether this element is located on the lower or upper side in those
71 : /// dimensions that it does not share with its parent hypercube
72 1 : const std::array<Side, HypercubeDim - ElementDim>& index() const;
73 :
74 1 : const Side& side_in_parent_dimension(size_t d) const;
75 :
76 : template <size_t NumIndices = HypercubeDim - ElementDim,
77 : Requires<NumIndices == 1> = nullptr>
78 1 : const Side& side() const {
79 : return index_[0];
80 : }
81 : /// @}
82 :
83 0 : bool operator==(const HypercubeElement& rhs) const {
84 : return dimensions_in_parent_ == rhs.dimensions_in_parent_ and
85 : index_ == rhs.index_;
86 : }
87 :
88 0 : bool operator!=(const HypercubeElement& rhs) const {
89 : return not(*this == rhs);
90 : }
91 :
92 : private:
93 0 : std::array<size_t, ElementDim> dimensions_in_parent_{};
94 0 : std::array<Side, HypercubeDim - ElementDim> index_{};
95 : };
96 :
97 : template <size_t ElementDim, size_t HypercubeDim>
98 0 : std::ostream& operator<<(
99 : std::ostream& os,
100 : const HypercubeElement<ElementDim, HypercubeDim>& element);
101 :
102 : /// A vertex in a `Dim`-dimensional hypercube
103 : template <size_t Dim>
104 1 : using Vertex = HypercubeElement<0, Dim>;
105 :
106 : /// An edge in a `Dim`-dimensional hypercube
107 : template <size_t Dim>
108 1 : using Edge = HypercubeElement<1, Dim>;
109 :
110 : /// A face in a `Dim`-dimensional hypercube
111 : template <size_t Dim>
112 1 : using Face = HypercubeElement<2, Dim>;
113 :
114 : /// A cell in a `Dim`-dimensional hypercube
115 : template <size_t Dim>
116 1 : using Cell = HypercubeElement<3, Dim>;
117 :
118 : /*!
119 : * \brief Iterator over all `ElementDim`-dimensional elements on the boundary of
120 : * a `HypercubeDim`-dimensional hypercube.
121 : *
122 : * \see `HypercubeElement`
123 : */
124 : template <size_t ElementDim, size_t HypercubeDim>
125 1 : struct HypercubeElementsIterator {
126 : static_assert(
127 : ElementDim <= HypercubeDim,
128 : "Hypercube element dimension must not exceed hypercube dimension.");
129 :
130 : public:
131 0 : static constexpr size_t num_indices = HypercubeDim - ElementDim;
132 :
133 : /// The number of `ElementDim`-dimensional elements on the boundary of a
134 : /// `HypercubeDim`-dimensional hypercube.
135 1 : static constexpr size_t size() {
136 : return two_to_the(num_indices) * factorial(HypercubeDim) /
137 : factorial(ElementDim) / factorial(num_indices);
138 : }
139 :
140 0 : HypercubeElementsIterator();
141 :
142 0 : static HypercubeElementsIterator begin();
143 :
144 0 : static HypercubeElementsIterator end();
145 :
146 0 : HypercubeElementsIterator& operator++();
147 :
148 0 : HypercubeElementsIterator operator++(int);
149 :
150 0 : HypercubeElement<ElementDim, HypercubeDim> operator*() const;
151 :
152 : private:
153 : template <size_t LocalElementDim = ElementDim,
154 : Requires<(LocalElementDim > 0)> = nullptr>
155 0 : void increment_dimension_in_parent(size_t d);
156 :
157 : template <size_t LocalElementDim, size_t LocalHypercubeDim>
158 : // NOLINTNEXTLINE(readability-redundant-declaration)
159 0 : friend bool operator==(
160 : const HypercubeElementsIterator<LocalElementDim, LocalHypercubeDim>& lhs,
161 : const HypercubeElementsIterator<LocalElementDim, LocalHypercubeDim>& rhs);
162 :
163 0 : std::array<size_t, ElementDim> dimensions_in_parent_{};
164 0 : size_t index_ = std::numeric_limits<size_t>::max();
165 : };
166 :
167 : template <size_t ElementDim, size_t HypercubeDim>
168 0 : bool operator!=(const HypercubeElementsIterator<ElementDim, HypercubeDim>& lhs,
169 : const HypercubeElementsIterator<ElementDim, HypercubeDim>& rhs);
170 :
171 : /// Iterate over all vertices in a `Dim`-dimensional hypercube
172 : template <size_t Dim>
173 1 : using VertexIterator = HypercubeElementsIterator<0, Dim>;
174 :
175 : /// Iterate over all edges in a `Dim`-dimensional hypercube
176 : template <size_t Dim>
177 1 : using EdgeIterator = HypercubeElementsIterator<1, Dim>;
178 :
179 : /// Iterate over all faces in a `Dim`-dimensional hypercube
180 : template <size_t Dim>
181 1 : using FaceIterator = HypercubeElementsIterator<2, Dim>;
182 :
183 : /// Iterate over all cells in a `Dim`-dimensional hypercube
184 : template <size_t Dim>
185 1 : using CellIterator = HypercubeElementsIterator<3, Dim>;
|