1 #ifndef _MARRAY_UTILITY_HPP_
2 #define _MARRAY_UTILITY_HPP_
3
4 #include <type_traits>
5 #include <array>
6 #include <vector>
7 #include <utility>
8 #include <iterator>
9 #include <cassert>
10 #include <algorithm>
11 #include <numeric>
12 #include <tuple>
13 #include <cmath>
14 #include <cstddef>
15 #include <string>
16 #include <functional>
17 #include <ostream>
18
19 #ifdef MARRAY_ENABLE_ASSERTS
20 #define MARRAY_ASSERT(e) assert(e)
21 #else
22 #define MARRAY_ASSERT(e) ((void)0)
23 #endif
24
25 #include "short_vector.hpp"
26
27 namespace MArray
28 {
29
30 /*
31 * The type all_t specifies a range [0,len_i) for an array
32 * dimension i of length len_i (i.e. it selects all of the data along
33 * that dimension).
34 */
all_tMArray::all_t35 struct all_t { constexpr all_t() {} };
bcast_tMArray::bcast_t36 struct bcast_t { constexpr bcast_t() {} };
37 namespace slice
38 {
39 constexpr all_t all;
40 constexpr bcast_t bcast;
41 }
42
43 #ifndef MARRAY_LEN_TYPE
44 #define MARRAY_LEN_TYPE ptrdiff_t
45 #endif
46
47 typedef MARRAY_LEN_TYPE len_type;
48
49 #ifndef MARRAY_STRIDE_TYPE
50 #define MARRAY_STRIDE_TYPE ptrdiff_t
51 #endif
52
53 typedef MARRAY_STRIDE_TYPE stride_type;
54
55 #ifndef MARRAY_OPT_NDIM
56 #define MARRAY_OPT_NDIM 6
57 #endif
58
59 typedef short_vector<len_type,MARRAY_OPT_NDIM> len_vector;
60 typedef short_vector<stride_type,MARRAY_OPT_NDIM> stride_vector;
61 typedef short_vector<unsigned,MARRAY_OPT_NDIM> dim_vector;
62 typedef short_vector<len_type,MARRAY_OPT_NDIM> index_vector;
63 typedef short_vector<unsigned,MARRAY_OPT_NDIM> irrep_vector;
64 template <typename T>
65 using ptr_vector = short_vector<T*,MARRAY_OPT_NDIM>;
66
67 #ifndef MARRAY_DEFAULT_LAYOUT
68 #define MARRAY_DEFAULT_LAYOUT ROW_MAJOR
69 #endif
70
71 #define MARRAY_PASTE_(x,y) x##y
72 #define MARRAY_PASTE(x,y) MARRAY_PASTE_(x,y)
73
74 #define MARRAY_DEFAULT_DPD_LAYOUT_(type) \
75 MARRAY_PASTE(MARRAY_PASTE(type,_),MARRAY_DEFAULT_LAYOUT)
76
77 #ifndef MARRAY_DEFAULT_DPD_LAYOUT
78 #define MARRAY_DEFAULT_DPD_LAYOUT PREFIX
79 #endif
80
81 /*
82 * The special value uninitialized is used to construct an array which
83 * does not default- or value-initialize its elements (useful for avoiding
84 * redundant memory operations for scalar types).
85 */
uninitialized_tMArray::uninitialized_t86 struct uninitialized_t { constexpr uninitialized_t() {} };
87 constexpr uninitialized_t uninitialized;
88
89 /*
90 * Specifies the layout of the array data.
91 */
92 struct layout
93 {
94 int type;
95
layoutMArray::layout96 constexpr explicit layout(int type) : type(type) {}
97
operator ==MArray::layout98 bool operator==(layout other) const { return type == other.type; }
operator !=MArray::layout99 bool operator!=(layout other) const { return type != other.type; }
100 };
101
column_major_layoutMArray::column_major_layout102 struct column_major_layout : layout { constexpr column_major_layout() : layout(0) {} };
103 constexpr column_major_layout COLUMN_MAJOR;
104
row_major_layoutMArray::row_major_layout105 struct row_major_layout : layout { constexpr row_major_layout() : layout(1) {} };
106 constexpr row_major_layout ROW_MAJOR;
107
108 constexpr decltype(MARRAY_DEFAULT_LAYOUT) DEFAULT;
109
110 struct dpd_layout
111 {
112 int type;
113
dpd_layoutMArray::dpd_layout114 constexpr explicit dpd_layout(int type) : type(type) {}
115
116 dpd_layout(layout layout);
117
118 layout base() const;
119
operator ==MArray::dpd_layout120 bool operator==(dpd_layout other) const { return type == other.type; }
operator !=MArray::dpd_layout121 bool operator!=(dpd_layout other) const { return type != other.type; }
122 };
123
124 struct balanced_column_major_layout : dpd_layout
balanced_column_major_layoutMArray::balanced_column_major_layout125 { constexpr balanced_column_major_layout() : dpd_layout(0) {} };
126 constexpr balanced_column_major_layout BALANCED_COLUMN_MAJOR;
127
128 struct balanced_row_major_layout : dpd_layout
balanced_row_major_layoutMArray::balanced_row_major_layout129 { constexpr balanced_row_major_layout() : dpd_layout(1) {} };
130 constexpr balanced_row_major_layout BALANCED_ROW_MAJOR;
131
132 struct blocked_column_major_layout : dpd_layout
blocked_column_major_layoutMArray::blocked_column_major_layout133 { constexpr blocked_column_major_layout() : dpd_layout(2) {} };
134 constexpr blocked_column_major_layout BLOCKED_COLUMN_MAJOR;
135
136 struct blocked_row_major_layout : dpd_layout
blocked_row_major_layoutMArray::blocked_row_major_layout137 { constexpr blocked_row_major_layout() : dpd_layout(3) {} };
138 constexpr blocked_row_major_layout BLOCKED_ROW_MAJOR;
139
140 struct prefix_column_major_layout : dpd_layout
prefix_column_major_layoutMArray::prefix_column_major_layout141 { constexpr prefix_column_major_layout() : dpd_layout(4) {} };
142 constexpr prefix_column_major_layout PREFIX_COLUMN_MAJOR;
143
144 struct prefix_row_major_layout : dpd_layout
prefix_row_major_layoutMArray::prefix_row_major_layout145 { constexpr prefix_row_major_layout() : dpd_layout(5) {} };
146 constexpr prefix_row_major_layout PREFIX_ROW_MAJOR;
147
148 constexpr decltype(MARRAY_DEFAULT_DPD_LAYOUT_(BALANCED)) BALANCED;
149 constexpr decltype(MARRAY_DEFAULT_DPD_LAYOUT_(BLOCKED)) BLOCKED;
150 constexpr decltype(MARRAY_DEFAULT_DPD_LAYOUT_(PREFIX)) PREFIX;
151
dpd_layout(layout layout)152 inline dpd_layout::dpd_layout(layout layout)
153 : type(layout == DEFAULT ? MARRAY_DEFAULT_DPD_LAYOUT.type :
154 layout == ROW_MAJOR ? MARRAY_PASTE(MARRAY_DEFAULT_DPD_LAYOUT,_ROW_MAJOR).type
155 : MARRAY_PASTE(MARRAY_DEFAULT_DPD_LAYOUT,_COLUMN_MAJOR).type) {}
156
base() const157 inline layout dpd_layout::base() const
158 {
159 if (*this == BALANCED_COLUMN_MAJOR ||
160 *this == BLOCKED_COLUMN_MAJOR ||
161 *this == PREFIX_COLUMN_MAJOR)
162 {
163 return COLUMN_MAJOR;
164 }
165 else
166 {
167 return ROW_MAJOR;
168 }
169 }
170
171 template <typename I> class range_t;
172
173 namespace detail
174 {
175 template <size_t N>
inverse_permutation(const std::array<unsigned,N> & p)176 std::array<unsigned, N> inverse_permutation(const std::array<unsigned, N>& p)
177 {
178 std::array<unsigned, N> ip;
179 for (unsigned i = 0;i < N;i++) ip[p[i]] = i;
180 return ip;
181 }
182
inverse_permutation(const dim_vector & p)183 inline dim_vector inverse_permutation(const dim_vector& p)
184 {
185 dim_vector ip(p.size());
186 for (unsigned i = 0;i < p.size();i++) ip[p[i]] = i;
187 return ip;
188 }
189
190 template <typename...>
191 struct exists {};
192
193 template <typename T>
194 using decay_t = typename std::decay<T>::type;
195
196 template <typename T>
197 using remove_cv_t = typename std::remove_cv<T>::type;
198
199 template <bool Cond, typename T=void>
200 using enable_if_t = typename std::enable_if<Cond,T>::type;
201
202 template <bool Cond, typename T, typename U>
203 using conditional_t = typename std::conditional<Cond,T,U>::type;
204
205 template <typename T, typename U=void>
206 using enable_if_integral_t = enable_if_t<std::is_integral<T>::value,U>;
207
208 template <typename T, typename U=void>
209 using enable_if_not_integral_t = enable_if_t<!std::is_integral<T>::value,U>;
210
211 template <typename T, typename U=void>
212 using enable_if_const_t = enable_if_t<std::is_const<T>::value,U>;
213
214 template <typename T, typename U=void>
215 using enable_if_not_const_t = enable_if_t<!std::is_const<T>::value,U>;
216
217 template <typename T, typename U, typename V=void>
218 using enable_if_assignable_t = enable_if_t<std::is_assignable<T,U>::value,V>;
219
220 template <typename T, typename U, typename V=void>
221 using enable_if_convertible_t = enable_if_t<std::is_convertible<T,U>::value,V>;
222
223 template <typename T, typename... Args>
224 struct are_convertible;
225
226 template <typename T>
227 struct are_convertible<T> : std::true_type {};
228
229 template <typename T, typename Arg, typename... Args>
230 struct are_convertible<T, Arg, Args...> :
231 conditional_t<std::is_convertible<Arg, T>::value,
232 are_convertible<T, Args...>,
233 std::false_type> {};
234
235 template <typename T, typename... Args>
236 struct are_assignable;
237
238 template <typename T>
239 struct are_assignable<T> : std::true_type {};
240
241 template <typename T, typename Arg, typename... Args>
242 struct are_assignable<T, Arg, Args...> :
243 conditional_t<std::is_assignable<T, Arg>::value,
244 are_assignable<T, Args...>,
245 std::false_type> {};
246
247 template <typename T, typename=void>
248 struct is_index_or_slice_helper : std::false_type {};
249
250 template <typename T>
251 struct is_index_or_slice_helper<T, enable_if_convertible_t<T, int>> : std::true_type {};
252
253 template <typename I>
254 struct is_index_or_slice_helper<range_t<I>, enable_if_integral_t<I>> : std::true_type {};
255
256 template <>
257 struct is_index_or_slice_helper<all_t> : std::true_type {};
258
259 template <>
260 struct is_index_or_slice_helper<bcast_t> : std::true_type {};
261
262 template <typename T>
263 struct is_index_or_slice : is_index_or_slice_helper<typename std::decay<T>::type> {};
264
265 template <typename... Args>
266 struct are_indices_or_slices;
267
268 template<>
269 struct are_indices_or_slices<> : std::true_type {};
270
271 template <typename Arg, typename... Args>
272 struct are_indices_or_slices<Arg, Args...> :
273 conditional_t<is_index_or_slice<Arg>::value,
274 are_indices_or_slices<Args...>,
275 std::false_type> {};
276
277 template <typename T, typename=void>
278 struct is_container : std::false_type {};
279
280 template <typename T>
281 struct is_container<T,
282 conditional_t<false,
283 exists<typename T::value_type,
284 decltype(std::declval<T>().size()),
285 decltype(std::declval<T>().begin()),
286 decltype(std::declval<T>().end())>,
287 void>>
288 : std::true_type {};
289
290 template <typename C, typename T, typename=void>
291 struct is_container_of : std::false_type {};
292
293 template <typename C, typename T>
294 struct is_container_of<C, T, typename std::enable_if<is_container<C>::value>::type>
295 : std::is_assignable<T&, typename C::value_type> {};
296
297 template <typename C, typename T, typename U=void>
298 using enable_if_container_of_t = typename std::enable_if<is_container_of<C, T>::value,U>::type;
299
300 template <typename C, typename=void>
301 struct is_container_of_containers : std::false_type {};
302
303 template <typename C>
304 struct is_container_of_containers<C, enable_if_t<is_container<C>::value>> :
305 is_container<typename C::value_type> {};
306
307 template <typename C, typename T, typename=void>
308 struct is_container_of_containers_of : std::false_type {};
309
310 template <typename C, typename T>
311 struct is_container_of_containers_of<C, T, enable_if_t<is_container<C>::value>> :
312 is_container_of<typename C::value_type, T> {};
313
314 template <typename C, typename T, typename U=void>
315 using enable_if_container_of_containers_of_t =
316 typename std::enable_if<is_container_of_containers_of<C, T>::value,U>::type;
317
318 template <typename T, typename... Ts>
319 struct are_containers_helper;
320
321 template <typename T>
322 struct are_containers_helper<T> : is_container<T> {};
323
324 template <typename T, typename... Ts>
325 struct are_containers_helper
326 : std::conditional<is_container<T>::value,
327 are_containers_helper<Ts...>,
328 std::false_type>::type {};
329
330 template <typename... Ts>
331 struct are_containers;
332
333 template <>
334 struct are_containers<> : std::true_type {};
335
336 template <typename... Ts>
337 struct are_containers : are_containers_helper<Ts...> {};
338
339 template <typename T, typename C, typename... Cs>
340 struct are_containers_of_helper;
341
342 template <typename T, typename C>
343 struct are_containers_of_helper<T, C> : is_container_of<C, T> {};
344
345 template <typename T, typename C, typename... Cs>
346 struct are_containers_of_helper
347 : std::conditional<is_container_of<C, T>::value,
348 are_containers_of_helper<T, Cs...>,
349 std::false_type>::type {};
350
351 template <typename T, typename... Cs>
352 struct are_containers_of;
353
354 template <typename T>
355 struct are_containers_of<T> : std::true_type {};
356
357 template <typename T, typename... Cs>
358 struct are_containers_of : are_containers_of_helper<T, Cs...> {};
359
360 template <typename Iterator>
inc_offsets_helper(unsigned i,Iterator)361 void inc_offsets_helper(unsigned i, Iterator) {}
362
363 template <typename Iterator, typename Offset, typename... Offsets>
inc_offsets_helper(unsigned i,Iterator it,Offset & off0,Offsets &...off)364 void inc_offsets_helper(unsigned i, Iterator it, Offset& off0,
365 Offsets&... off)
366 {
367 off0 += (*it)[i];
368 inc_offsets_helper(i, ++it, off...);
369 }
370
371 template <typename Strides, typename... Offsets>
inc_offsets(unsigned i,const Strides & strides,Offsets &...off)372 void inc_offsets(unsigned i, const Strides& strides, Offsets&... off)
373 {
374 inc_offsets_helper(i, strides.begin(), off...);
375 }
376
377 template <typename Pos, typename Iterator>
dec_offsets_helper(unsigned i,const Pos &,Iterator)378 void dec_offsets_helper(unsigned i, const Pos&, Iterator) {}
379
380 template <typename Pos, typename Iterator, typename Offset, typename... Offsets>
dec_offsets_helper(unsigned i,const Pos & pos,Iterator it,Offset & off0,Offsets &...off)381 void dec_offsets_helper(unsigned i, const Pos& pos, Iterator it,
382 Offset& off0, Offsets&... off)
383 {
384 off0 -= pos[i]*(*it)[i];
385 dec_offsets_helper(i, pos, ++it, off...);
386 }
387
388 template <typename Pos, typename Strides, typename... Offsets>
dec_offsets(unsigned i,const Pos & pos,const Strides & strides,Offsets &...off)389 void dec_offsets(unsigned i, const Pos& pos, const Strides& strides,
390 Offsets&... off)
391 {
392 dec_offsets_helper(i, pos, strides.begin(), off...);
393 }
394
395 template <typename Pos, typename Iterator>
move_offsets_helper(const Pos &,Iterator)396 void move_offsets_helper(const Pos&, Iterator) {}
397
398 template <typename Pos, typename Iterator, typename Offset, typename... Offsets>
move_offsets_helper(const Pos & pos,Iterator it,Offset & off0,Offsets &...off)399 void move_offsets_helper(const Pos& pos, Iterator it,
400 Offset& off0, Offsets&... off)
401 {
402 for (unsigned i = 0;i < pos.size();i++) off0 += pos[i]*(*it)[i];
403 move_offsets_helper(pos, ++it, off...);
404 }
405
406 template <typename Pos, typename Strides, typename... Offsets>
move_offsets(const Pos & pos,const Strides & strides,Offsets &...off)407 void move_offsets(const Pos& pos, const Strides& strides,
408 Offsets&... off)
409 {
410 move_offsets_helper(pos, strides.begin(), off...);
411 }
412
413 template <typename Iterator>
set_strides_helper(Iterator)414 void set_strides_helper(Iterator) {}
415
416 template <typename Iterator, typename Stride, typename... Strides>
set_strides_helper(Iterator it,const Stride & stride,const Strides &...strides)417 void set_strides_helper(Iterator it, const Stride& stride, const Strides&... strides)
418 {
419 std::copy(stride.begin(), stride.end(), it->begin());
420 set_strides_helper(++it, strides...);
421 }
422
423 template <typename Strides_, typename... Strides>
set_strides(Strides_ & strides_,const Strides &...strides)424 void set_strides(Strides_& strides_, const Strides&... strides)
425 {
426 set_strides_helper(strides_.begin(), strides...);
427 }
428
429 template <typename T, T... S> struct integer_sequence {};
430
431 template <typename T, typename U, typename V> struct concat_sequences;
432 template <typename T, T... S, T... R>
433 struct concat_sequences<T, integer_sequence<T, S...>, integer_sequence<T, R...>>
434 {
435 typedef integer_sequence<T, S..., (R+sizeof...(S))...> type;
436 };
437
438 template <typename T, T N, typename=void> struct static_range_helper;
439
440 template <typename T, T N> struct static_range_helper<T, N, enable_if_t<N==0>>
441 {
442 typedef integer_sequence<T> type;
443 };
444
445 template <typename T, T N> struct static_range_helper<T, N, enable_if_t<N==1>>
446 {
447 typedef integer_sequence<T, 0> type;
448 };
449
450 template <typename T, T N> struct static_range_helper<T, N, enable_if_t<(N>1)>>
451 {
452 typedef typename concat_sequences<T, typename static_range_helper<T, (N+1)/2>::type,
453 typename static_range_helper<T, N/2>::type>::type type;
454 };
455
456 template <typename T, T N>
457 using static_range = typename static_range_helper<T, N>::type;
458 }
459
460 }
461
462 #endif
463