Branch data Line data Source code
1 : : //
2 : : // Class Index
3 : : // Define a slice in an array.
4 : : //
5 : : // This essentially defines a list of evenly spaced numbers.
6 : : // Most commonly this list will be increasing (positive stride)
7 : : // but it can also have negative stride and be decreasing.
8 : : //
9 : : // Index() --> A null interval with no elements.
10 : : // Index(n) --> make an Index on [0..n-1]
11 : : // Index(a,b) --> make an Index on [a..b]
12 : : // Index(a,b,s) --> make an Index on [a..b] with stride s
13 : : //
14 : : // Example1:
15 : : // --------
16 : : // Index I(10); --> Index on [0..9]
17 : : // Index Low(5); --> Index on [0..4]
18 : : // Index High(5,9); --> Index on [5..9]
19 : : // Index IOdd(1,9,2); --> Index on [1..9] stride 2
20 : : // Index IEven(0,9,2); --> Index on [0..9] stride 2
21 : : //
22 : : // Given an Index I(a,n,s), and an integer j you can do the following:
23 : : //
24 : : // I+j : a+j+i*s for i in [0..n-1]
25 : : // j-I : j-a-i*s
26 : : // j*I : j*a + i*j*s
27 : : // I/j : a/j + i*s/j
28 : : //
29 : : // j/I we don't do because it is not a uniform stride, and we don't
30 : : // allow strides that are fractions.
31 : : //
32 : : #include "Utility/PAssert.h"
33 : :
34 : : namespace ippl {
35 : :
36 : 14996 : KOKKOS_INLINE_FUNCTION Index::Index()
37 : 14996 : : first_m(0)
38 : 14996 : , stride_m(0)
39 : 14996 : , length_m(0) {}
40 : :
41 : 2478 : KOKKOS_INLINE_FUNCTION Index::Index(size_t n)
42 : 2478 : : first_m(0)
43 : 2478 : , stride_m(1)
44 : 2478 : , length_m(n) {}
45 : :
46 : 1252 : KOKKOS_INLINE_FUNCTION Index::Index(int f, int l)
47 : 1252 : : first_m(f)
48 : 1252 : , stride_m(1)
49 : 1252 : , length_m(l - f + 1) {
50 [ - + ]: 1252 : PAssert_GE(l - f + 1, 0);
51 : 1252 : }
52 : :
53 : 0 : KOKKOS_INLINE_FUNCTION Index::Index(int f, int l, int s)
54 : 0 : : first_m(f)
55 : 0 : , stride_m(s) {
56 [ # # ]: 0 : PAssert_NE(s, 0);
57 [ # # ]: 0 : if (f == l) {
58 : 0 : length_m = 1;
59 [ # # ]: 0 : } else if ((l > f) ^ (s < 0)) {
60 : 0 : length_m = (l - f) / s + 1;
61 : : } else {
62 : 0 : length_m = 0;
63 : : }
64 : 0 : }
65 : :
66 : 0 : KOKKOS_INLINE_FUNCTION Index::Index(int m, int a, const Index& b)
67 : 0 : : first_m(b.first_m * m + a)
68 : 0 : , stride_m(b.stride_m * m)
69 : 0 : , length_m(b.length_m) {}
70 : :
71 : : KOKKOS_INLINE_FUNCTION Index::Index(int f, int s, const Index* b)
72 : : : first_m(f)
73 : : , stride_m(s)
74 : : , length_m(b->length_m) {}
75 : :
76 : 28011048 : KOKKOS_INLINE_FUNCTION int Index::first() const noexcept {
77 : 28011048 : return first_m;
78 : : }
79 : :
80 : 42 : KOKKOS_INLINE_FUNCTION int Index::stride() const noexcept {
81 : 42 : return stride_m;
82 : : }
83 : :
84 : : KOKKOS_INLINE_FUNCTION bool Index::empty() const noexcept {
85 : : return length_m == 0;
86 : : }
87 : :
88 : 13926230 : KOKKOS_INLINE_FUNCTION size_t Index::length() const noexcept {
89 : 13926230 : return length_m;
90 : : }
91 : :
92 : 1044 : KOKKOS_INLINE_FUNCTION int Index::last() const noexcept {
93 [ - + ]: 1044 : return (length_m == 0) ? first_m : first_m + stride_m * (length_m - 1);
94 : : }
95 : :
96 : 840 : KOKKOS_INLINE_FUNCTION int Index::min() const noexcept {
97 [ + - ]: 840 : return (stride_m >= 0) ? first_m : first_m + stride_m * (length_m - 1);
98 : : }
99 : :
100 : 840 : KOKKOS_INLINE_FUNCTION int Index::max() const noexcept {
101 [ + - ]: 840 : return (stride_m >= 0) ? first_m + stride_m * (length_m - 1) : first_m;
102 : : }
103 : :
104 : 0 : KOKKOS_INLINE_FUNCTION Index& Index::operator+=(int off) {
105 : 0 : first_m += off;
106 : 0 : return *this;
107 : : }
108 : :
109 : 0 : KOKKOS_INLINE_FUNCTION Index& Index::operator-=(int off) {
110 : 0 : first_m -= off;
111 : 0 : return *this;
112 : : }
113 : :
114 : 0 : KOKKOS_INLINE_FUNCTION Index operator+(const Index& i, int off) {
115 : 0 : return Index(1, off, i);
116 : : }
117 : :
118 : : KOKKOS_INLINE_FUNCTION Index operator+(int off, const Index& i) {
119 : : return Index(1, off, i);
120 : : }
121 : :
122 : 0 : KOKKOS_INLINE_FUNCTION Index operator-(const Index& i, int off) {
123 : 0 : return Index(1, -off, i);
124 : : }
125 : :
126 : : KOKKOS_INLINE_FUNCTION Index operator-(int off, const Index& i) {
127 : : return Index(-1, off, i);
128 : : }
129 : :
130 : : KOKKOS_INLINE_FUNCTION Index operator-(const Index& i) {
131 : : return Index(-1, 0, i);
132 : : }
133 : :
134 : : KOKKOS_INLINE_FUNCTION Index operator*(const Index& i, int m) {
135 : : return Index(m, 0, i);
136 : : }
137 : :
138 : : KOKKOS_INLINE_FUNCTION Index operator*(int m, const Index& i) {
139 : : return Index(m, 0, i);
140 : : }
141 : :
142 : : KOKKOS_INLINE_FUNCTION Index operator/(const Index& i, int d) {
143 : : return Index(i.first_m / d, i.stride_m / d, &i);
144 : : }
145 : :
146 : 0 : KOKKOS_INLINE_FUNCTION Index Index::reverse() const {
147 : 0 : Index j;
148 : 0 : j.first_m = last();
149 : 0 : j.length_m = length_m;
150 : 0 : j.stride_m = -stride_m;
151 : 0 : return j;
152 : : }
153 : :
154 : 0 : KOKKOS_INLINE_FUNCTION bool Index::touches(const Index& a) const {
155 [ # # # # ]: 0 : return (min() <= a.max()) && (max() >= a.min());
156 : : }
157 : :
158 : : KOKKOS_INLINE_FUNCTION bool Index::contains(const Index& a) const {
159 : : return (min() <= a.min()) && (max() >= a.max());
160 : : }
161 : :
162 : 0 : KOKKOS_INLINE_FUNCTION bool Index::split(Index& l, Index& r) const {
163 [ # # ]: 0 : PAssert_EQ(stride_m, 1);
164 [ # # ]: 0 : PAssert_GT(length_m, 1);
165 : 0 : auto mid = first_m + length_m / 2 - 1;
166 [ # # ]: 0 : l = Index(first_m, mid);
167 [ # # ]: 0 : r = Index(mid + 1, first_m + length_m - 1);
168 : 0 : return true;
169 : : }
170 : :
171 : 0 : KOKKOS_INLINE_FUNCTION bool Index::split(Index& l, Index& r, int i) const {
172 [ # # ]: 0 : PAssert_EQ(stride_m, 1);
173 [ # # ]: 0 : PAssert_GT(length_m, 1);
174 [ # # ]: 0 : if (i >= first_m + static_cast<int>(length_m)) {
175 : 0 : return false;
176 : : }
177 [ # # ]: 0 : l = Index(first_m, i);
178 [ # # ]: 0 : r = Index(i + 1, first_m + length_m - 1);
179 : 0 : return true;
180 : : }
181 : :
182 : 0 : KOKKOS_INLINE_FUNCTION bool Index::split(Index& l, Index& r, double a) const {
183 [ # # ]: 0 : PAssert_EQ(stride_m, 1);
184 [ # # ]: 0 : PAssert_GT(length_m, 1);
185 [ # # ]: 0 : PAssert_LT(a, 1.0);
186 [ # # ]: 0 : PAssert_GT(a, 0.0);
187 : 0 : int mid = first_m + static_cast<int>(length_m * a + 0.5) - 1;
188 [ # # ]: 0 : l = Index(first_m, mid);
189 [ # # ]: 0 : r = Index(mid + 1, first_m + length_m - 1);
190 : 0 : return true;
191 : : }
192 : :
193 : : //////////////////////////////////////////////////////////////////////
194 : : // Calculate the lowest common multiple of s1 and s2.
195 : : // put the result in s.
196 : : // Also calculate m1 = s/s1 and m2 = s/s2.
197 : : // This version is optimized for small s1 and s2 and
198 : : // just uses an exhaustive search.
199 : : //////////////////////////////////////////////////////////////////////
200 : 0 : KOKKOS_INLINE_FUNCTION void lcm(int s1, int s2, int& s, int& m1, int& m2) {
201 [ # # ]: 0 : PAssert_GT(s1, 0); // For simplicity, make some assumptions.
202 [ # # ]: 0 : PAssert_GT(s2, 0);
203 : 0 : int i1 = s1;
204 : 0 : int i2 = s2;
205 : 0 : int _m1 = 1;
206 : 0 : int _m2 = 1;
207 [ # # ]: 0 : if (i2 < i1) {
208 : : while (true) {
209 [ # # ]: 0 : while (i2 < i1) {
210 : 0 : i2 += s2;
211 : 0 : ++_m2;
212 : : }
213 [ # # ]: 0 : if (i1 == i2) {
214 : 0 : m1 = _m1;
215 : 0 : m2 = _m2;
216 : 0 : s = i1;
217 : 0 : return;
218 : : }
219 : 0 : i1 += s1;
220 : 0 : ++_m1;
221 : : }
222 : : } else {
223 : : while (true) {
224 [ # # ]: 0 : while (i1 < i2) {
225 : 0 : i1 += s1;
226 : 0 : ++_m1;
227 : : }
228 [ # # ]: 0 : if (i1 == i2) {
229 : 0 : m1 = _m1;
230 : 0 : m2 = _m2;
231 : 0 : s = i1;
232 : 0 : return;
233 : : }
234 : 0 : i2 += s2;
235 : 0 : ++_m2;
236 : : }
237 : : }
238 : : }
239 : :
240 : : //
241 : : // Intersect, with the code for the common case of
242 : : // both strides equal to one.
243 : : //
244 : 0 : KOKKOS_INLINE_FUNCTION Index Index::intersect(const Index& rhs) const {
245 : 0 : Index ret;
246 [ # # # # : 0 : if ((stride() == 1) && (rhs.stride() == 1)) {
# # ]
247 : 0 : int lf = first();
248 : 0 : int rf = rhs.first();
249 : 0 : int ll = last();
250 : 0 : int rl = rhs.last();
251 [ # # ]: 0 : int f = lf > rf ? lf : rf;
252 [ # # ]: 0 : int l = ll < rl ? ll : rl;
253 : 0 : ret.first_m = f;
254 [ # # ]: 0 : ret.length_m = ((l >= f) ? l - f + 1 : 0);
255 : 0 : ret.stride_m = 1;
256 : : } else
257 [ # # ]: 0 : ret = general_intersect(rhs);
258 : 0 : return ret;
259 : : }
260 : :
261 : 84 : KOKKOS_INLINE_FUNCTION Index Index::grow(int ncells) const {
262 : 84 : Index index;
263 : :
264 : 84 : index.first_m = this->first_m - ncells;
265 : 84 : index.length_m = this->length_m + 2 * ncells;
266 : 84 : index.stride_m = this->stride_m;
267 : :
268 : 84 : return index;
269 : : }
270 : :
271 : 0 : KOKKOS_INLINE_FUNCTION static Index do_intersect(const Index& a, const Index& b) {
272 [ # # # # ]: 0 : PAssert_GT(a.stride(), 0); // This should be assured by the
273 [ # # # # ]: 0 : PAssert_GT(b.stride(), 0); // caller of this function.
274 : :
275 : : int newStride; // The stride for the new index is
276 : : int a_mul, b_mul; // a_mul=newStride/a.stride() ...
277 [ # # ]: 0 : lcm(a.stride(), b.stride(), // The input strides...
278 : : newStride, a_mul, b_mul); // the lcm of the strides of a and b.
279 : :
280 : : // Find the offset from a.first() in units of newStride
281 : : // that puts the ranges close together.
282 : 0 : int a_i = (b.first() - a.first()) / a.stride();
283 : 0 : int a_off = a.first() + a_i * a.stride();
284 [ # # ]: 0 : if (a_off < b.first()) {
285 : 0 : a_i++;
286 : 0 : a_off += a.stride();
287 : : }
288 : :
289 [ # # # # ]: 0 : PAssert_GE(a_off, b.first()); // make sure I'm understanding this right...
290 : :
291 : : // Now do an exhaustive search for the first point in common.
292 : : // Count over all possible offsets for a.
293 [ # # # # : 0 : for (int a_m = 0; (a_m < a_mul) && (a_i < (int)a.length());
# # ]
294 : 0 : a_m++, a_i++, a_off += a.stride()) {
295 : 0 : int b_off = b.first();
296 : : // Count over all possible offsets for b.
297 [ # # # # : 0 : for (int b_m = 0; (b_m < b_mul) && (b_m < (int)b.length()); b_m++, b_off += b.stride())
# # ]
298 [ # # ]: 0 : if (a_off == b_off) { // If the offsets are the same, we found it!
299 : 0 : int am = a.max(); // Find the minimum maximum of a and b...
300 : 0 : int bm = b.max();
301 [ # # ]: 0 : int m = am < bm ? am : bm;
302 [ # # ]: 0 : return Index(a_off, m, newStride);
303 : : }
304 : : }
305 : 0 : return Index(0); // If we get to here there is no intersection.
306 : : }
307 : :
308 : 0 : KOKKOS_INLINE_FUNCTION Index Index::general_intersect(const Index& that) const {
309 : : // If they just don't overlap, return null indexes.
310 [ # # # # : 0 : if ((min() > that.max()) || (that.min() > max()))
# # ]
311 : 0 : return Index(0);
312 [ # # # # ]: 0 : if ((stride_m == 0) || (that.stride_m == 0))
313 : 0 : return Index(0);
314 : :
315 : : // If one or the other counts -ve, reverse it and intersect result.
316 [ # # ]: 0 : if (that.stride_m < 0)
317 [ # # ]: 0 : return intersect(that.reverse());
318 [ # # ]: 0 : if (stride_m < 0) {
319 : 0 : Index r;
320 [ # # ]: 0 : r = reverse().intersect(that).reverse();
321 : 0 : return r;
322 : : }
323 : :
324 : : // Getting closer to the real thing: intersect them.
325 : : // Pass the one that starts lower as the first argument.
326 : 0 : Index r;
327 [ # # ]: 0 : if (first_m < that.first_m)
328 [ # # ]: 0 : r = do_intersect(*this, that);
329 : : else
330 [ # # ]: 0 : r = do_intersect(that, *this);
331 : :
332 : 0 : return r;
333 : : }
334 : : } // namespace ippl
|
