xtensor
 
Loading...
Searching...
No Matches
xfixed.hpp
1/***************************************************************************
2 * Copyright (c) Johan Mabille, Sylvain Corlay and Wolf Vollprecht *
3 * Copyright (c) QuantStack *
4 * *
5 * Distributed under the terms of the BSD 3-Clause License. *
6 * *
7 * The full license is in the file LICENSE, distributed with this software. *
8 ****************************************************************************/
9
10#ifndef XTENSOR_FIXED_HPP
11#define XTENSOR_FIXED_HPP
12
13#include <algorithm>
14#include <array>
15#include <cstddef>
16#include <utility>
17#include <vector>
18
19#include <xtl/xsequence.hpp>
20
21#include "../containers/xcontainer.hpp"
22#include "../containers/xstorage.hpp"
23#include "../core/xsemantic.hpp"
24#include "../core/xstrides.hpp"
25#include "../core/xtensor_config.hpp"
26
27namespace xtl
28{
29 namespace detail
30 {
31 template <class T, std::size_t N>
32 struct sequence_builder<xt::const_array<T, N>>
33 {
34 using sequence_type = xt::const_array<T, N>;
35 using value_type = typename sequence_type::value_type;
36 using size_type = typename sequence_type::size_type;
37
38 inline static sequence_type make(size_type /*size*/, value_type /*v*/)
39 {
40 return sequence_type();
41 }
42 };
43 }
44}
45
46namespace xt
47{
48
49 /**********************
50 * xfixed declaration *
51 **********************/
52
53 template <class ET, class S, layout_type L, bool SH, class Tag>
54 class xfixed_container;
55
56 namespace detail
57 {
58 /**************************************************************************************
59 The following is something we can currently only dream about -- for when we drop
60 support for a lot of the old compilers (e.g. GCC 4.9, MSVC 2017 ;)
61
62 template <class T>
63 constexpr std::size_t calculate_stride(T& shape, std::size_t idx, layout_type L)
64 {
65 if (shape[idx] == 1)
66 {
67 return std::size_t(0);
68 }
69
70 std::size_t data_size = 1;
71 std::size_t stride = 1;
72 if (L == layout_type::row_major)
73 {
74 // because we have a integer sequence that counts
75 // from 0 to sz - 1, we need to "invert" idx here
76 idx = shape.size() - idx;
77 for (std::size_t i = idx; i != 0; --i)
78 {
79 stride = data_size;
80 data_size = stride * shape[i - 1];
81 }
82 }
83 else
84 {
85 for (std::size_t i = 0; i < idx + 1; ++i)
86 {
87 stride = data_size;
88 data_size = stride * shape[i];
89 }
90 }
91 return stride;
92 }
93
94 *****************************************************************************************/
95
96 template <layout_type L, std::size_t I, std::size_t... X>
97 struct calculate_stride;
98
99 template <std::size_t I, std::size_t Y, std::size_t... X>
100 struct calculate_stride<layout_type::column_major, I, Y, X...>
101 {
102 static constexpr std::ptrdiff_t value = Y
103 * calculate_stride<layout_type::column_major, I - 1, X...>::value;
104 };
105
106 template <std::size_t Y, std::size_t... X>
107 struct calculate_stride<layout_type::column_major, 0, Y, X...>
108 {
109 static constexpr std::ptrdiff_t value = 1;
110 };
111
112 template <std::size_t I, std::size_t... X>
113 struct calculate_stride_row_major
114 {
115 static constexpr std::ptrdiff_t value = at<sizeof...(X) - I, X...>::value
116 * calculate_stride_row_major<I - 1, X...>::value;
117 };
118
119 template <std::size_t... X>
120 struct calculate_stride_row_major<0, X...>
121 {
122 static constexpr std::ptrdiff_t value = 1;
123 };
124
125 template <std::size_t I, std::size_t... X>
126 struct calculate_stride<layout_type::row_major, I, X...>
127 {
128 static constexpr std::ptrdiff_t value = calculate_stride_row_major<sizeof...(X) - I - 1, X...>::value;
129 };
130
131 namespace workaround
132 {
133 template <layout_type L, size_t I, class SEQ>
134 struct computed_strides;
135
136 template <layout_type L, size_t I, size_t... X>
137 struct computed_strides<L, I, std::index_sequence<X...>>
138 {
139 static constexpr std::ptrdiff_t value = calculate_stride<L, I, X...>::value;
140 };
141
142 template <layout_type L, size_t I, class SEQ>
143 constexpr std::ptrdiff_t get_computed_strides(bool cond)
144 {
145 return cond ? 0 : computed_strides<L, I, SEQ>::value;
146 }
147 }
148
149 template <layout_type L, class R, std::size_t... X, std::size_t... I>
150 constexpr R get_strides_impl(const xt::fixed_shape<X...>& shape, std::index_sequence<I...>)
151 {
152 static_assert(
154 "Layout not supported for fixed array"
155 );
156#if (_MSC_VER >= 1910)
157 using temp_type = std::index_sequence<X...>;
158 return R({workaround::get_computed_strides<L, I, temp_type>(shape[I] == 1)...});
159#else
160 return R({shape[I] == 1 ? 0 : calculate_stride<L, I, X...>::value...});
161#endif
162 }
163
164 template <class S, class T, std::size_t... I>
165 constexpr T get_backstrides_impl(const S& shape, const T& strides, std::index_sequence<I...>)
166 {
167 return T({(strides[I] * std::ptrdiff_t(shape[I] - 1))...});
168 }
169
170 template <std::size_t... X>
171 struct fixed_compute_size_impl;
172
173 template <std::size_t Y, std::size_t... X>
174 struct fixed_compute_size_impl<Y, X...>
175 {
176 static constexpr std::size_t value = Y * fixed_compute_size_impl<X...>::value;
177 };
178
179 template <std::size_t X>
180 struct fixed_compute_size_impl<X>
181 {
182 static constexpr std::size_t value = X;
183 };
184
185 template <>
186 struct fixed_compute_size_impl<>
187 {
188 // support for 0D xtensor fixed (empty shape = xshape<>)
189 static constexpr std::size_t value = 1;
190 };
191
192 // TODO unify with constexpr compute_size when dropping MSVC 2015
193 template <class T>
194 struct fixed_compute_size;
195
196 template <std::size_t... X>
197 struct fixed_compute_size<xt::fixed_shape<X...>>
198 {
199 static constexpr std::size_t value = fixed_compute_size_impl<X...>::value;
200 };
201
202 template <class V, std::size_t... X>
203 struct get_init_type_impl;
204
205 template <class V, std::size_t Y>
206 struct get_init_type_impl<V, Y>
207 {
208 using type = V[Y];
209 };
210
211 template <class V>
212 struct get_init_type_impl<V>
213 {
214 using type = V[1];
215 };
216
217 template <class V, std::size_t Y, std::size_t... X>
218 struct get_init_type_impl<V, Y, X...>
219 {
220 using tmp_type = typename get_init_type_impl<V, X...>::type;
221 using type = tmp_type[Y];
222 };
223 }
224
225 template <layout_type L, class R, std::size_t... X>
226 constexpr R get_strides(const fixed_shape<X...>& shape) noexcept
227 {
228 return detail::get_strides_impl<L, R>(shape, std::make_index_sequence<sizeof...(X)>{});
229 }
230
231 template <class S, class T>
232 constexpr T get_backstrides(const S& shape, const T& strides) noexcept
233 {
234 return detail::get_backstrides_impl(shape, strides, std::make_index_sequence<std::tuple_size<T>::value>{});
235 }
236
237 template <class V, class S>
239
240 template <class V, std::size_t... X>
241 struct get_init_type<V, fixed_shape<X...>>
242 {
243 using type = typename detail::get_init_type_impl<V, X...>::type;
244 };
245
246 template <class V, class S>
247 using get_init_type_t = typename get_init_type<V, S>::type;
248
249 template <class ET, class S, layout_type L, bool SH, class Tag>
250 struct xcontainer_inner_types<xfixed_container<ET, S, L, SH, Tag>>
251 {
252 using shape_type = S;
253 using inner_shape_type = typename S::cast_type;
254 using strides_type = get_strides_t<inner_shape_type>;
255 using inner_strides_type = strides_type;
256 using backstrides_type = inner_strides_type;
257 using inner_backstrides_type = backstrides_type;
258
259 // NOTE: 0D (S::size() == 0) results in storage for 1 element (scalar)
260#if defined(_MSC_VER) && _MSC_VER < 1910 && !defined(_WIN64)
261 // WORKAROUND FOR MSVC 2015 32 bit, fallback to unaligned container for 0D scalar case
262 using storage_type = std::array<ET, detail::fixed_compute_size<S>::value>;
263#else
265#endif
266
267 using reference = typename storage_type::reference;
268 using const_reference = typename storage_type::const_reference;
269 using size_type = typename storage_type::size_type;
270 using temporary_type = xfixed_container<ET, S, L, SH, Tag>;
271 static constexpr layout_type layout = L;
272 };
273
274 template <class ET, class S, layout_type L, bool SH, class Tag>
275 struct xiterable_inner_types<xfixed_container<ET, S, L, SH, Tag>>
276 : xcontainer_iterable_types<xfixed_container<ET, S, L, SH, Tag>>
277 {
278 };
279
295 template <class ET, class S, layout_type L, bool SH, class Tag>
296 class xfixed_container : public xcontainer<xfixed_container<ET, S, L, SH, Tag>>,
297 public xcontainer_semantic<xfixed_container<ET, S, L, SH, Tag>>
298 {
299 public:
300
301 using self_type = xfixed_container<ET, S, L, SH, Tag>;
302 using base_type = xcontainer<self_type>;
303 using semantic_base = xcontainer_semantic<self_type>;
304
305 using storage_type = typename base_type::storage_type;
306 using value_type = typename base_type::value_type;
307 using reference = typename base_type::reference;
308 using const_reference = typename base_type::const_reference;
309 using pointer = typename base_type::pointer;
310 using const_pointer = typename base_type::const_pointer;
311 using shape_type = typename base_type::shape_type;
312 using inner_shape_type = typename base_type::inner_shape_type;
313 using strides_type = typename base_type::strides_type;
314 using backstrides_type = typename base_type::backstrides_type;
315 using inner_backstrides_type = typename base_type::inner_backstrides_type;
316 using inner_strides_type = typename base_type::inner_strides_type;
317 using temporary_type = typename semantic_base::temporary_type;
318 using expression_tag = Tag;
319
320 static constexpr std::size_t N = std::tuple_size<shape_type>::value;
321 static constexpr std::size_t rank = N;
322
323 xfixed_container() = default;
324 xfixed_container(const value_type& v);
325 explicit xfixed_container(const inner_shape_type& shape, layout_type l = L);
326 explicit xfixed_container(const inner_shape_type& shape, value_type v, layout_type l = L);
327
328 template <class IX = std::integral_constant<std::size_t, N>>
329 xfixed_container(nested_initializer_list_t<value_type, N> t)
330 requires(IX::value != 0);
331
332 ~xfixed_container() = default;
333
334 xfixed_container(const xfixed_container&) = default;
335 xfixed_container& operator=(const xfixed_container&) = default;
336
337 xfixed_container(xfixed_container&&) = default;
338 xfixed_container& operator=(xfixed_container&&) = default;
339
340 template <class E>
342
343 template <class E>
344 xfixed_container& operator=(const xexpression<E>& e);
345
346 template <class ST = std::array<std::size_t, N>>
347 static xfixed_container from_shape(ST&& /*s*/);
348
349 template <class ST = std::array<std::size_t, N>>
350 void resize(ST&& shape, bool force = false) const;
351 template <class ST = shape_type>
352 void resize(ST&& shape, layout_type l) const;
353 template <class ST = shape_type>
354 void resize(ST&& shape, const strides_type& strides) const;
355
356 template <class ST = std::array<std::size_t, N>>
357 const auto& reshape(ST&& shape, layout_type layout = L) const;
358
359 template <class ST>
360 bool broadcast_shape(ST& s, bool reuse_cache = false) const;
361
362 constexpr layout_type layout() const noexcept;
363 bool is_contiguous() const noexcept;
364
365 private:
366
367 storage_type m_storage;
368
369 XTENSOR_CONSTEXPR_ENHANCED_STATIC inner_shape_type m_shape = S();
370 XTENSOR_CONSTEXPR_ENHANCED_STATIC inner_strides_type m_strides = get_strides<L, inner_strides_type>(S());
371 XTENSOR_CONSTEXPR_ENHANCED_STATIC inner_backstrides_type
372 m_backstrides = get_backstrides(m_shape, m_strides);
373
374 storage_type& storage_impl() noexcept;
375 const storage_type& storage_impl() const noexcept;
376
377 XTENSOR_CONSTEXPR_RETURN const inner_shape_type& shape_impl() const noexcept;
378 XTENSOR_CONSTEXPR_RETURN const inner_strides_type& strides_impl() const noexcept;
379 XTENSOR_CONSTEXPR_RETURN const inner_backstrides_type& backstrides_impl() const noexcept;
380
381 friend class xcontainer<xfixed_container<ET, S, L, SH, Tag>>;
382 };
383
384 /****************************************
385 * xfixed_container_adaptor declaration *
386 ****************************************/
387
388 template <class EC, class S, layout_type L, bool SH, class Tag>
389 class xfixed_adaptor;
390
391 template <class EC, class S, layout_type L, bool SH, class Tag>
392 struct xcontainer_inner_types<xfixed_adaptor<EC, S, L, SH, Tag>>
393 {
394 using storage_type = std::remove_reference_t<EC>;
395 using reference = typename storage_type::reference;
396 using const_reference = typename storage_type::const_reference;
397 using size_type = typename storage_type::size_type;
398 using shape_type = S;
399 using inner_shape_type = typename S::cast_type;
400 using strides_type = get_strides_t<inner_shape_type>;
401 using backstrides_type = strides_type;
402 using inner_strides_type = strides_type;
403 using inner_backstrides_type = backstrides_type;
405 static constexpr layout_type layout = L;
406 };
407
408 template <class EC, class S, layout_type L, bool SH, class Tag>
409 struct xiterable_inner_types<xfixed_adaptor<EC, S, L, SH, Tag>>
410 : xcontainer_iterable_types<xfixed_adaptor<EC, S, L, SH, Tag>>
411 {
412 };
413
430 template <class EC, class S, layout_type L, bool SH, class Tag>
431 class xfixed_adaptor : public xcontainer<xfixed_adaptor<EC, S, L, SH, Tag>>,
432 public xcontainer_semantic<xfixed_adaptor<EC, S, L, SH, Tag>>
433 {
434 public:
435
436 using container_closure_type = EC;
437
438 using self_type = xfixed_adaptor<EC, S, L, SH, Tag>;
439 using base_type = xcontainer<self_type>;
440 using semantic_base = xcontainer_semantic<self_type>;
441 using storage_type = typename base_type::storage_type;
442 using shape_type = typename base_type::shape_type;
443 using strides_type = typename base_type::strides_type;
444 using backstrides_type = typename base_type::backstrides_type;
445 using inner_shape_type = typename base_type::inner_shape_type;
446 using inner_strides_type = typename base_type::inner_strides_type;
447 using inner_backstrides_type = typename base_type::inner_backstrides_type;
448 using temporary_type = typename semantic_base::temporary_type;
449 using expression_tag = Tag;
450
451 static constexpr std::size_t N = S::size();
452
453 xfixed_adaptor(storage_type&& data);
454 xfixed_adaptor(const storage_type& data);
455
456 template <class D>
457 xfixed_adaptor(D&& data);
458
459 ~xfixed_adaptor() = default;
460
461 xfixed_adaptor(const xfixed_adaptor&) = default;
462 xfixed_adaptor& operator=(const xfixed_adaptor&);
463
464 xfixed_adaptor(xfixed_adaptor&&) = default;
465 xfixed_adaptor& operator=(xfixed_adaptor&&);
466 xfixed_adaptor& operator=(temporary_type&&);
467
468 template <class E>
469 xfixed_adaptor& operator=(const xexpression<E>& e);
470
471 template <class ST = std::array<std::size_t, N>>
472 void resize(ST&& shape, bool force = false) const;
473 template <class ST = shape_type>
474 void resize(ST&& shape, layout_type l) const;
475 template <class ST = shape_type>
476 void resize(ST&& shape, const strides_type& strides) const;
477
478 template <class ST = std::array<std::size_t, N>>
479 const auto& reshape(ST&& shape, layout_type layout = L) const;
480
481 template <class ST>
482 bool broadcast_shape(ST& s, bool reuse_cache = false) const;
483
484 constexpr layout_type layout() const noexcept;
485 bool is_contiguous() const noexcept;
486
487 private:
488
489 container_closure_type m_storage;
490
491 XTENSOR_CONSTEXPR_ENHANCED_STATIC inner_shape_type m_shape = S();
492 XTENSOR_CONSTEXPR_ENHANCED_STATIC inner_strides_type m_strides = get_strides<L, inner_strides_type>(S());
493 XTENSOR_CONSTEXPR_ENHANCED_STATIC inner_backstrides_type
494 m_backstrides = get_backstrides(m_shape, m_strides);
495
496 storage_type& storage_impl() noexcept;
497 const storage_type& storage_impl() const noexcept;
498
499 XTENSOR_CONSTEXPR_RETURN const inner_shape_type& shape_impl() const noexcept;
500 XTENSOR_CONSTEXPR_RETURN const inner_strides_type& strides_impl() const noexcept;
501 XTENSOR_CONSTEXPR_RETURN const inner_backstrides_type& backstrides_impl() const noexcept;
502
503 friend class xcontainer<xfixed_adaptor<EC, S, L, SH, Tag>>;
504 };
505
506 /************************************
507 * xfixed_container implementation *
508 ************************************/
509
514
523 template <class ET, class S, layout_type L, bool SH, class Tag>
524 inline xfixed_container<ET, S, L, SH, Tag>::xfixed_container(const inner_shape_type& shape, layout_type l)
525 {
526 (void) (shape);
527 (void) (l);
528 XTENSOR_ASSERT(shape.size() == N && std::equal(shape.begin(), shape.end(), m_shape.begin()));
529 XTENSOR_ASSERT(L == l);
530 }
531
532 template <class ET, class S, layout_type L, bool SH, class Tag>
533 inline xfixed_container<ET, S, L, SH, Tag>::xfixed_container(const value_type& v)
534 {
535 if (this->size() != 1)
536 {
537 XTENSOR_THROW(std::runtime_error, "wrong shape for scalar assignment (has to be xshape<>).");
538 }
539 m_storage[0] = v;
540 }
541
551 template <class ET, class S, layout_type L, bool SH, class Tag>
552 inline xfixed_container<ET, S, L, SH, Tag>::xfixed_container(
553 const inner_shape_type& shape,
554 value_type v,
556 )
557 {
558 (void) (shape);
559 (void) (l);
560 XTENSOR_ASSERT(shape.size() == N && std::equal(shape.begin(), shape.end(), m_shape.begin()));
561 XTENSOR_ASSERT(L == l);
562 std::fill(m_storage.begin(), m_storage.end(), v);
563 }
564
565 namespace detail
566 {
567 template <std::size_t X>
568 struct check_initializer_list_shape
569 {
570 template <class T, class S>
571 static bool run(const T& t, const S& shape)
572 {
573 std::size_t IX = shape.size() - X;
574 bool result = (shape[IX] == t.size());
575 for (std::size_t i = 0; i < shape[IX]; ++i)
576 {
577 result = result && check_initializer_list_shape<X - 1>::run(t.begin()[i], shape);
578 }
579 return result;
580 }
581 };
582
583 template <>
584 struct check_initializer_list_shape<0>
585 {
586 template <class T, class S>
587 static bool run(const T& /*t*/, const S& /*shape*/)
588 {
589 return true;
590 }
591 };
592 }
593
594 template <class ET, class S, layout_type L, bool SH, class Tag>
595 template <class ST>
596 inline xfixed_container<ET, S, L, SH, Tag> xfixed_container<ET, S, L, SH, Tag>::from_shape(ST&& shape)
597 {
598 (void) shape;
599 self_type tmp;
600 XTENSOR_ASSERT(shape.size() == N && std::equal(shape.begin(), shape.end(), tmp.shape().begin()));
601 return tmp;
602 }
603
611 template <class ET, class S, layout_type L, bool SH, class Tag>
612 template <class IX>
613 inline xfixed_container<ET, S, L, SH, Tag>::xfixed_container(nested_initializer_list_t<value_type, N> t)
614 requires(IX::value != 0)
615 {
616 XTENSOR_ASSERT_MSG(
617 detail::check_initializer_list_shape<N>::run(t, this->shape()) == true,
618 "initializer list shape does not match fixed shape"
619 );
620 constexpr auto tmp = layout_type::row_major;
621 L == tmp ? nested_copy(m_storage.begin(), t) : nested_copy(this->template begin<tmp>(), t);
622 }
623
625
630
633 template <class ET, class S, layout_type L, bool SH, class Tag>
634 template <class E>
635 inline xfixed_container<ET, S, L, SH, Tag>::xfixed_container(const xexpression<E>& e)
636 {
637 semantic_base::assign(e);
638 }
639
643 template <class ET, class S, layout_type L, bool SH, class Tag>
644 template <class E>
645 inline auto xfixed_container<ET, S, L, SH, Tag>::operator=(const xexpression<E>& e) -> self_type&
646 {
647 return semantic_base::operator=(e);
648 }
649
651
656 template <class ET, class S, layout_type L, bool SH, class Tag>
657 template <class ST>
659 {
660 (void) (shape); // remove unused parameter warning if XTENSOR_ASSERT undefined
661 XTENSOR_ASSERT(std::equal(shape.begin(), shape.end(), m_shape.begin()) && shape.size() == m_shape.size());
662 }
663
668 template <class ET, class S, layout_type L, bool SH, class Tag>
669 template <class ST>
671 {
672 (void) (shape); // remove unused parameter warning if XTENSOR_ASSERT undefined
673 (void) (l);
674 XTENSOR_ASSERT(
675 std::equal(shape.begin(), shape.end(), m_shape.begin()) && shape.size() == m_shape.size() && L == l
676 );
677 }
678
683 template <class ET, class S, layout_type L, bool SH, class Tag>
684 template <class ST>
685 inline void xfixed_container<ET, S, L, SH, Tag>::resize(ST&& shape, const strides_type& strides) const
686 {
687 (void) (shape); // remove unused parameter warning if XTENSOR_ASSERT undefined
688 (void) (strides);
689 XTENSOR_ASSERT(std::equal(shape.begin(), shape.end(), m_shape.begin()) && shape.size() == m_shape.size());
690 XTENSOR_ASSERT(
691 std::equal(strides.begin(), strides.end(), m_strides.begin()) && strides.size() == m_strides.size()
692 );
693 }
694
698 template <class ET, class S, layout_type L, bool SH, class Tag>
699 template <class ST>
701 {
702 if (!(std::equal(shape.begin(), shape.end(), m_shape.begin()) && shape.size() == m_shape.size()
703 && layout == L))
704 {
705 XTENSOR_THROW(std::runtime_error, "Trying to reshape xtensor_fixed with different shape or layout.");
706 }
707 return *this;
708 }
709
710 template <class ET, class S, layout_type L, bool SH, class Tag>
711 template <class ST>
712 inline bool xfixed_container<ET, S, L, SH, Tag>::broadcast_shape(ST& shape, bool) const
713 {
714 return xt::broadcast_shape(m_shape, shape);
715 }
716
717 template <class ET, class S, layout_type L, bool SH, class Tag>
718 constexpr layout_type xfixed_container<ET, S, L, SH, Tag>::layout() const noexcept
719 {
720 return base_type::static_layout;
721 }
722
723 template <class ET, class S, layout_type L, bool SH, class Tag>
724 inline bool xfixed_container<ET, S, L, SH, Tag>::is_contiguous() const noexcept
725 {
726 using str_type = typename inner_strides_type::value_type;
727 return m_strides.empty() || (layout() == layout_type::row_major && m_strides.back() == str_type(1))
728 || (layout() == layout_type::column_major && m_strides.front() == str_type(1));
729 }
730
731 template <class ET, class S, layout_type L, bool SH, class Tag>
732 inline auto xfixed_container<ET, S, L, SH, Tag>::storage_impl() noexcept -> storage_type&
733 {
734 return m_storage;
735 }
736
737 template <class ET, class S, layout_type L, bool SH, class Tag>
738 inline auto xfixed_container<ET, S, L, SH, Tag>::storage_impl() const noexcept -> const storage_type&
739 {
740 return m_storage;
741 }
742
743 template <class ET, class S, layout_type L, bool SH, class Tag>
744 XTENSOR_CONSTEXPR_RETURN auto xfixed_container<ET, S, L, SH, Tag>::shape_impl() const noexcept
745 -> const inner_shape_type&
746 {
747 return m_shape;
748 }
749
750 template <class ET, class S, layout_type L, bool SH, class Tag>
751 XTENSOR_CONSTEXPR_RETURN auto xfixed_container<ET, S, L, SH, Tag>::strides_impl() const noexcept
752 -> const inner_strides_type&
753 {
754 return m_strides;
755 }
756
757 template <class ET, class S, layout_type L, bool SH, class Tag>
758 XTENSOR_CONSTEXPR_RETURN auto xfixed_container<ET, S, L, SH, Tag>::backstrides_impl() const noexcept
759 -> const inner_backstrides_type&
760 {
761 return m_backstrides;
762 }
763
764 /*******************
765 * xfixed_adaptor *
766 *******************/
767
772
776 template <class EC, class S, layout_type L, bool SH, class Tag>
778 : base_type()
779 , m_storage(std::move(data))
780 {
781 }
782
787 template <class EC, class S, layout_type L, bool SH, class Tag>
789 : base_type()
790 , m_storage(data)
791 {
792 }
793
799 template <class EC, class S, layout_type L, bool SH, class Tag>
800 template <class D>
802 : base_type()
803 , m_storage(std::forward<D>(data))
804 {
805 }
806
808
809 template <class EC, class S, layout_type L, bool SH, class Tag>
810 inline auto xfixed_adaptor<EC, S, L, SH, Tag>::operator=(const xfixed_adaptor& rhs) -> self_type&
811 {
812 base_type::operator=(rhs);
813 m_storage = rhs.m_storage;
814 return *this;
815 }
816
817 template <class EC, class S, layout_type L, bool SH, class Tag>
818 inline auto xfixed_adaptor<EC, S, L, SH, Tag>::operator=(xfixed_adaptor&& rhs) -> self_type&
819 {
820 base_type::operator=(std::move(rhs));
821 m_storage = rhs.m_storage;
822 return *this;
823 }
824
825 template <class EC, class S, layout_type L, bool SH, class Tag>
826 inline auto xfixed_adaptor<EC, S, L, SH, Tag>::operator=(temporary_type&& rhs) -> self_type&
827 {
828 m_storage.resize(rhs.storage().size());
829 std::copy(rhs.storage().cbegin(), rhs.storage().cend(), m_storage.begin());
830 return *this;
831 }
832
837
840 template <class EC, class S, layout_type L, bool SH, class Tag>
841 template <class E>
842 inline auto xfixed_adaptor<EC, S, L, SH, Tag>::operator=(const xexpression<E>& e) -> self_type&
843 {
844 return semantic_base::operator=(e);
845 }
846
848
853 template <class ET, class S, layout_type L, bool SH, class Tag>
854 template <class ST>
856 {
857 (void) (shape); // remove unused parameter warning if XTENSOR_ASSERT undefined
858 XTENSOR_ASSERT(std::equal(shape.begin(), shape.end(), m_shape.begin()) && shape.size() == m_shape.size());
859 }
860
865 template <class ET, class S, layout_type L, bool SH, class Tag>
866 template <class ST>
868 {
869 (void) (shape); // remove unused parameter warning if XTENSOR_ASSERT undefined
870 (void) (l);
871 XTENSOR_ASSERT(
872 std::equal(shape.begin(), shape.end(), m_shape.begin()) && shape.size() == m_shape.size() && L == l
873 );
874 }
875
880 template <class ET, class S, layout_type L, bool SH, class Tag>
881 template <class ST>
882 inline void xfixed_adaptor<ET, S, L, SH, Tag>::resize(ST&& shape, const strides_type& strides) const
883 {
884 (void) (shape); // remove unused parameter warning if XTENSOR_ASSERT undefined
885 (void) (strides);
886 XTENSOR_ASSERT(std::equal(shape.begin(), shape.end(), m_shape.begin()) && shape.size() == m_shape.size());
887 XTENSOR_ASSERT(
888 std::equal(strides.begin(), strides.end(), m_strides.begin()) && strides.size() == m_strides.size()
889 );
890 }
891
895 template <class ET, class S, layout_type L, bool SH, class Tag>
896 template <class ST>
897 inline const auto& xfixed_adaptor<ET, S, L, SH, Tag>::reshape(ST&& shape, layout_type layout) const
898 {
899 if (!(std::equal(shape.begin(), shape.end(), m_shape.begin()) && shape.size() == m_shape.size()
900 && layout == L))
901 {
902 XTENSOR_THROW(std::runtime_error, "Trying to reshape xtensor_fixed with different shape or layout.");
903 }
904 return *this;
905 }
906
907 template <class ET, class S, layout_type L, bool SH, class Tag>
908 template <class ST>
909 inline bool xfixed_adaptor<ET, S, L, SH, Tag>::broadcast_shape(ST& shape, bool) const
910 {
911 return xt::broadcast_shape(m_shape, shape);
912 }
913
914 template <class EC, class S, layout_type L, bool SH, class Tag>
915 inline auto xfixed_adaptor<EC, S, L, SH, Tag>::storage_impl() noexcept -> storage_type&
916 {
917 return m_storage;
918 }
919
920 template <class EC, class S, layout_type L, bool SH, class Tag>
921 inline auto xfixed_adaptor<EC, S, L, SH, Tag>::storage_impl() const noexcept -> const storage_type&
922 {
923 return m_storage;
924 }
925
926 template <class EC, class S, layout_type L, bool SH, class Tag>
927 constexpr layout_type xfixed_adaptor<EC, S, L, SH, Tag>::layout() const noexcept
928 {
929 return base_type::static_layout;
930 }
931
932 template <class EC, class S, layout_type L, bool SH, class Tag>
933 inline bool xfixed_adaptor<EC, S, L, SH, Tag>::is_contiguous() const noexcept
934 {
935 using str_type = typename inner_strides_type::value_type;
936 return m_strides.empty() || (layout() == layout_type::row_major && m_strides.back() == str_type(1))
937 || (layout() == layout_type::column_major && m_strides.front() == str_type(1));
938 }
939
940 template <class EC, class S, layout_type L, bool SH, class Tag>
941 XTENSOR_CONSTEXPR_RETURN auto xfixed_adaptor<EC, S, L, SH, Tag>::shape_impl() const noexcept
942 -> const inner_shape_type&
943 {
944 return m_shape;
945 }
946
947 template <class EC, class S, layout_type L, bool SH, class Tag>
948 XTENSOR_CONSTEXPR_RETURN auto xfixed_adaptor<EC, S, L, SH, Tag>::strides_impl() const noexcept
949 -> const inner_strides_type&
950 {
951 return m_strides;
952 }
953
954 template <class EC, class S, layout_type L, bool SH, class Tag>
955 XTENSOR_CONSTEXPR_RETURN auto xfixed_adaptor<EC, S, L, SH, Tag>::backstrides_impl() const noexcept
956 -> const inner_backstrides_type&
957 {
958 return m_backstrides;
959 }
960}
961
962#endif
This array class is modeled after std::array but adds optional alignment through a template parameter...
Fixed shape implementation for compile time defined arrays.
constexpr const inner_strides_type & strides() const noexcept
Returns the strides of the container.
constexpr const inner_shape_type & shape() const noexcept
Returns the shape of the container.
Base class for xexpressions.
Dense multidimensional container adaptor with tensor semantic and fixed dimension.
Definition xfixed.hpp:433
xfixed_adaptor(storage_type &&data)
Constructs an xfixed_adaptor of the given stl-like container.
Definition xfixed.hpp:777
const auto & reshape(ST &&shape, layout_type layout=L) const
Note that the xfixed_container cannot be reshaped to a shape different from S.
Definition xfixed.hpp:897
void resize(ST &&shape, bool force=false) const
Note that the xfixed_adaptor cannot be resized.
Definition xfixed.hpp:855
Dense multidimensional container with tensor semantic and fixed dimension.
Definition xfixed.hpp:298
xfixed_container(const inner_shape_type &shape, layout_type l=L)
Create an uninitialized xfixed_container.
Definition xfixed.hpp:524
xfixed_container(const xexpression< E > &e)
The extended copy constructor.
Definition xfixed.hpp:635
void resize(ST &&shape, bool force=false) const
Note that the xfixed_container cannot be resized.
Definition xfixed.hpp:658
void resize(ST &&shape, const strides_type &strides) const
Note that the xfixed_container cannot be resized.
Definition xfixed.hpp:685
xfixed_container(const inner_shape_type &shape, value_type v, layout_type l=L)
Create an xfixed_container, and initialize with the value of v.
Definition xfixed.hpp:552
void resize(ST &&shape, layout_type l) const
Note that the xfixed_container cannot be resized.
Definition xfixed.hpp:670
xfixed_container(nested_initializer_list_t< value_type, N > t)
Allocates an xfixed_container with shape S with values from a C array.
Definition xfixed.hpp:613
const auto & reshape(ST &&shape, layout_type layout=L) const
Note that the xfixed_container cannot be reshaped to a shape different from S.
Definition xfixed.hpp:700
auto strides(const E &e, stride_type type=stride_type::normal) noexcept
Get strides of an object.
Definition xstrides.hpp:248
standard mathematical functions for xexpressions
layout_type
Definition xlayout.hpp:24