Line data Source code
1 0 : // Distributed under the MIT License.
2 : // See LICENSE.txt for details.
3 :
4 : #pragma once
5 :
6 : #include <algorithm>
7 : #include <array>
8 : #include <blaze/math/AlignmentFlag.h>
9 : #include <blaze/math/CustomVector.h>
10 : #include <blaze/math/DenseVector.h>
11 : #include <blaze/math/GroupTag.h>
12 : #include <blaze/math/PaddingFlag.h>
13 : #include <blaze/math/TransposeFlag.h>
14 : #include <cstddef>
15 : #include <cstring>
16 : #include <functional>
17 : #include <initializer_list>
18 : #include <limits>
19 : #include <memory>
20 : #include <ostream>
21 : #include <pup.h>
22 : #include <type_traits>
23 :
24 : #include "DataStructures/Blaze/StepFunction.hpp"
25 : #include "Utilities/ErrorHandling/Assert.hpp"
26 : #include "Utilities/ForceInline.hpp"
27 : #include "Utilities/Gsl.hpp"
28 : #include "Utilities/MakeString.hpp"
29 : #include "Utilities/MakeWithValue.hpp"
30 : #include "Utilities/MemoryHelpers.hpp"
31 : #include "Utilities/PrintHelpers.hpp"
32 : #include "Utilities/Requires.hpp"
33 : #include "Utilities/SetNumberOfGridPoints.hpp"
34 : #include "Utilities/StdArrayHelpers.hpp"
35 : #include "Utilities/TypeTraits/IsComplexOfFundamental.hpp"
36 : #include "Utilities/TypeTraits/IsStdArray.hpp"
37 :
38 : class ComplexDataVector;
39 : class ComplexModalVector;
40 : class DataVector;
41 : class ModalVector;
42 :
43 : namespace VectorImpl_detail {
44 : /// \brief Whether or not a given vector type is assignable to another
45 : ///
46 : /// \details
47 : /// This is used to define which types can be assigned to one another. For
48 : /// example, you can assign a `ComplexDataVector` to a `DataVector`, but not
49 : /// vice versa.
50 : ///
51 : /// To enable assignments between more types, modify a current template
52 : /// specialization or add a new one.
53 : ///
54 : /// \tparam LhsDataType the type being assigned
55 : /// \tparam RhsDataType the type to convert to `LhsDataType`
56 : template <typename LhsDataType, typename RhsDataType>
57 : struct is_assignable;
58 :
59 : /// No template specialization was matched, so LHS is not assignable to RHS
60 : template <typename LhsDataType, typename RhsDataType>
61 : struct is_assignable : std::false_type {};
62 : /// Can assign a type to itself
63 : template <typename RhsDataType>
64 : struct is_assignable<RhsDataType, RhsDataType> : std::true_type {};
65 : /// Can assign a `ComplexDataVector` to a `DataVector`
66 : template <>
67 : struct is_assignable<ComplexDataVector, DataVector> : std::true_type {};
68 : /// Can assign a `ComplexModalVector` to a `ModalVector`
69 : template <>
70 : struct is_assignable<ComplexModalVector, ModalVector> : std::true_type {};
71 :
72 : /// \brief Whether or not a given vector type is assignable to another
73 : ///
74 : /// \details
75 : /// See `is_assignable` for which assignments are permitted
76 : template <typename LhsDataType, typename RhsDataType>
77 : constexpr bool is_assignable_v = is_assignable<LhsDataType, RhsDataType>::value;
78 : } // namespace VectorImpl_detail
79 :
80 : /// \ingroup TensorExpressionsGroup
81 : /// \brief Marks a class as being a `VectorImpl`
82 : ///
83 : /// \details
84 : /// The empty base class provides a simple means for checking if a type is a
85 : /// `VectorImpl`
86 1 : struct MarkAsVectorImpl {};
87 :
88 : /// \ingroup DataStructuresGroup
89 : /// Default static size for vector impl
90 1 : constexpr size_t default_vector_impl_static_size = 0;
91 :
92 : /*!
93 : * \ingroup DataStructuresGroup
94 : * \brief Base class template for various DataVector and related types
95 : *
96 : * \details The `VectorImpl` class is the generic parent class for vectors
97 : * representing collections of related function values, such as `DataVector`s
98 : * for contiguous data over a computational domain.
99 : *
100 : * The `VectorImpl` does not itself define any particular mathematical
101 : * operations on the contained values. The `VectorImpl` template class and the
102 : * macros defined in `VectorImpl.hpp` assist in the construction of various
103 : * derived classes supporting a chosen set of mathematical operations.
104 : *
105 : * In addition, the equivalence operator `==` is inherited from the underlying
106 : * `blaze::CustomVector` type, and returns true if and only if the size and
107 : * contents of the two compared vectors are equivalent.
108 : *
109 : * Template parameters:
110 : * - `T` is the underlying stored type, e.g. `double`, `std::complex<double>`,
111 : * `float`, etc.
112 : * - `VectorType` is the type that should be associated with the VectorImpl
113 : * during mathematical computations. In most cases, inherited types should
114 : * have themselves as the second template argument, e.g.
115 : * ```
116 : * class DataVector : VectorImpl<double, DataVector> {
117 : * ```
118 : * - `StaticSize` is the size for the static part of the vector. If the vector
119 : * is constructed or resized with a size that is less than or equal to this
120 : * StaticSize, no heap allocations will be done. It will instead use the stack
121 : * allocation. Default is `default_vector_impl_static_size`.
122 : *
123 : * The second template parameter communicates arithmetic type restrictions to
124 : * the underlying Blaze framework. For example, if `VectorType` is
125 : * `DataVector`, then the underlying architecture will prevent addition with a
126 : * vector type whose `ResultType` (which is aliased to its `VectorType`) is
127 : * `ModalVector`. Since `DataVector`s and `ModalVector`s represent data in
128 : * different spaces, we wish to forbid several operations between them. This
129 : * vector-type-tracking through an expression prevents accidental mixing of
130 : * vector types in math expressions.
131 : *
132 : * \note
133 : * - If `SPECTRE_NAN_INIT` is defined, then the `VectorImpl` is default
134 : * initialized to `signaling_NaN()`. Otherwise, the vector is filled with
135 : * uninitialized memory for performance.
136 : */
137 : template <typename T, typename VectorType,
138 : size_t StaticSize = default_vector_impl_static_size>
139 1 : class VectorImpl
140 : : public blaze::CustomVector<
141 : T, blaze::AlignmentFlag::unaligned, blaze::PaddingFlag::unpadded,
142 : blaze::defaultTransposeFlag, blaze::GroupTag<0>, VectorType>,
143 : MarkAsVectorImpl {
144 : public:
145 0 : using value_type = T;
146 0 : using size_type = size_t;
147 0 : using difference_type = std::ptrdiff_t;
148 0 : using BaseType = blaze::CustomVector<
149 : T, blaze::AlignmentFlag::unaligned, blaze::PaddingFlag::unpadded,
150 : blaze::defaultTransposeFlag, blaze::GroupTag<0>, VectorType>;
151 0 : static constexpr bool transpose_flag = blaze::defaultTransposeFlag;
152 0 : static constexpr size_t static_size = StaticSize;
153 :
154 0 : using ElementType = T;
155 0 : using TransposeType = VectorImpl<T, VectorType, StaticSize>;
156 0 : using CompositeType = const VectorImpl<T, VectorType, StaticSize>&;
157 0 : using iterator = typename BaseType::Iterator;
158 0 : using const_iterator = typename BaseType::ConstIterator;
159 :
160 : using BaseType::operator[];
161 : using BaseType::begin;
162 : using BaseType::cbegin;
163 : using BaseType::cend;
164 : using BaseType::data;
165 : using BaseType::end;
166 : using BaseType::size;
167 :
168 : /// @{
169 : /// Upcast to `BaseType`
170 : /// \attention
171 : /// upcast should only be used when implementing a derived vector type, not in
172 : /// calling code
173 1 : const BaseType& operator*() const {
174 : return static_cast<const BaseType&>(*this);
175 : }
176 1 : BaseType& operator*() { return static_cast<BaseType&>(*this); }
177 : /// @}
178 :
179 : /// Create with the given size. In debug mode, the vector is initialized to
180 : /// 'NaN' by default. If not initialized to 'NaN', the memory is allocated but
181 : /// not initialized.
182 : ///
183 : /// - `set_size` number of values
184 1 : explicit VectorImpl(size_t set_size)
185 : : owned_data_(heap_alloc_if_necessary(set_size)) {
186 : reset_pointer_vector(set_size);
187 : #if defined(SPECTRE_NAN_INIT)
188 : std::fill(data(), data() + set_size,
189 : std::numeric_limits<value_type>::signaling_NaN());
190 : #endif // SPECTRE_NAN_INIT
191 : }
192 :
193 : /// Create with the given size and value.
194 : ///
195 : /// - `set_size` number of values
196 : /// - `value` the value to initialize each element
197 1 : VectorImpl(size_t set_size, T value)
198 : : owned_data_(heap_alloc_if_necessary(set_size)) {
199 : reset_pointer_vector(set_size);
200 : std::fill(data(), data() + set_size, value);
201 : }
202 :
203 : /// Create from a copy of the given container
204 : ///
205 : /// \param container A container with a `value_type` that is the same as `T`.
206 : /// Currently restricted to `std::vector<T>` and `std::array<T>`.
207 : template <
208 : typename Container,
209 : Requires<std::is_same_v<typename Container::value_type, T>> = nullptr>
210 1 : explicit VectorImpl(const Container& container)
211 : : owned_data_(heap_alloc_if_necessary(container.size())) {
212 : static_assert(std::is_same_v<Container, std::vector<T>> or
213 : tt::is_std_array_v<Container>,
214 : "This constructor is currently restricted to std::vector and "
215 : "std::array out of caution.");
216 : reset_pointer_vector(container.size());
217 : std::copy(container.begin(), container.end(), data());
218 : }
219 :
220 : /// Create a non-owning VectorImpl that points to `start`
221 1 : VectorImpl(T* start, size_t set_size)
222 : : BaseType(start, set_size), owning_(false) {}
223 :
224 : /// Create from an initializer list of `T`.
225 : template <class U, Requires<std::is_same_v<U, T>> = nullptr>
226 1 : VectorImpl(std::initializer_list<U> list)
227 : : owned_data_(heap_alloc_if_necessary(list.size())) {
228 : reset_pointer_vector(list.size());
229 : // Note: can't use memcpy with an initializer list.
230 : std::copy(list.begin(), list.end(), data());
231 : }
232 :
233 : /// Empty VectorImpl
234 1 : VectorImpl() = default;
235 : /// \cond HIDDEN_SYMBOLS
236 : ~VectorImpl() = default;
237 :
238 : VectorImpl(const VectorImpl<T, VectorType, StaticSize>& rhs);
239 : VectorImpl& operator=(const VectorImpl<T, VectorType, StaticSize>& rhs);
240 : VectorImpl(VectorImpl<T, VectorType, StaticSize>&& rhs);
241 : VectorImpl& operator=(VectorImpl<T, VectorType, StaticSize>&& rhs);
242 :
243 : // This is a converting constructor. clang-tidy complains that it's not
244 : // explicit, but we want it to allow conversion.
245 : // clang-tidy: mark as explicit (we want conversion to VectorImpl type)
246 : template <typename VT, bool VF,
247 : Requires<VectorImpl_detail::is_assignable_v<
248 : VectorType, typename VT::ResultType>> = nullptr>
249 : VectorImpl(const blaze::DenseVector<VT, VF>& expression); // NOLINT
250 :
251 : template <typename VT, bool VF>
252 : VectorImpl& operator=(const blaze::DenseVector<VT, VF>& expression);
253 : /// \endcond
254 :
255 0 : VectorImpl& operator=(const T& rhs);
256 :
257 0 : decltype(auto) SPECTRE_ALWAYS_INLINE operator[](const size_t index) {
258 : ASSERT(index < size(), "Out-of-range access to element "
259 : << index << " of a size " << size()
260 : << " Blaze vector.");
261 : return BaseType::operator[](index);
262 : }
263 :
264 0 : decltype(auto) SPECTRE_ALWAYS_INLINE operator[](const size_t index) const {
265 : ASSERT(index < size(), "Out-of-range access to element "
266 : << index << " of a size " << size()
267 : << " Blaze vector.");
268 : return BaseType::operator[](index);
269 : }
270 :
271 : /// @{
272 : /// Set the VectorImpl to be a reference to another VectorImpl object
273 1 : void set_data_ref(gsl::not_null<VectorType*> rhs) {
274 : set_data_ref(rhs->data(), rhs->size());
275 : }
276 :
277 1 : void set_data_ref(T* const start, const size_t set_size) {
278 : clear();
279 : if (start != nullptr) {
280 : (**this).reset(start, set_size);
281 : }
282 : owning_ = false;
283 : }
284 : /// @}
285 :
286 : /*!
287 : * \brief A common operation for checking the size and resizing a memory
288 : * buffer if needed to ensure that it has the desired size. This operation is
289 : * not permitted on a non-owning vector.
290 : *
291 : * \note This utility should NOT be used when it is anticipated that the
292 : * supplied buffer will typically be the wrong size (in that case, suggest
293 : * either manual checking or restructuring so that resizing is less common).
294 : * This uses `UNLIKELY` to perform the check most quickly when the buffer
295 : * needs no resizing, but will be slower when resizing is common.
296 : */
297 1 : void SPECTRE_ALWAYS_INLINE destructive_resize(const size_t new_size) {
298 : if (UNLIKELY(size() != new_size)) {
299 : ASSERT(owning_,
300 : MakeString{}
301 : << "Attempting to resize a non-owning vector from size: "
302 : << size() << " to size: " << new_size
303 : << " but we may not destructively resize a non-owning vector");
304 : owned_data_ = heap_alloc_if_necessary(new_size);
305 : reset_pointer_vector(new_size);
306 : }
307 : }
308 :
309 : /// Returns true if the class owns the data
310 1 : bool is_owning() const { return owning_; }
311 :
312 : /// Put the class in the default-constructed state.
313 1 : void clear();
314 :
315 : /// Serialization for Charm++
316 : // NOLINTNEXTLINE(google-runtime-references)
317 1 : void pup(PUP::er& p);
318 :
319 : protected:
320 0 : std::unique_ptr<value_type[]> owned_data_{};
321 0 : std::array<T, StaticSize> static_owned_data_{};
322 0 : bool owning_{true};
323 :
324 : // This should only be called if we are owning. If we are not owning, then
325 : // neither owned_data_ or static_owned_data_ actually has the data we want.
326 0 : SPECTRE_ALWAYS_INLINE void reset_pointer_vector(const size_t set_size) {
327 : if (set_size == 0) {
328 : BaseType::clear();
329 : return;
330 : }
331 : if (owned_data_ == nullptr and set_size > StaticSize) {
332 : ERROR(
333 : "VectorImpl::reset_pointer_vector cannot be called when owned_data_ "
334 : "is nullptr.");
335 : }
336 :
337 : if (set_size <= StaticSize) {
338 : this->reset(static_owned_data_.data(), set_size);
339 : // Free memory if downsizing
340 : owned_data_ = nullptr;
341 : } else {
342 : this->reset(owned_data_.get(), set_size);
343 : }
344 : }
345 :
346 0 : SPECTRE_ALWAYS_INLINE std::unique_ptr<value_type[]> heap_alloc_if_necessary(
347 : const size_t set_size) {
348 : return set_size > StaticSize
349 : ? cpp20::make_unique_for_overwrite<value_type[]>(set_size)
350 : : nullptr;
351 : }
352 : };
353 :
354 : /// \cond HIDDEN_SYMBOLS
355 : template <typename T, typename VectorType, size_t StaticSize>
356 : VectorImpl<T, VectorType, StaticSize>::VectorImpl(
357 : const VectorImpl<T, VectorType, StaticSize>& rhs)
358 : : BaseType{rhs}, owned_data_(heap_alloc_if_necessary(rhs.size())) {
359 : reset_pointer_vector(rhs.size());
360 : std::memcpy(data(), rhs.data(), size() * sizeof(value_type));
361 : }
362 :
363 : template <typename T, typename VectorType, size_t StaticSize>
364 : VectorImpl<T, VectorType, StaticSize>&
365 : VectorImpl<T, VectorType, StaticSize>::operator=(
366 : const VectorImpl<T, VectorType, StaticSize>& rhs) {
367 : if (this != &rhs) {
368 : if (owning_) {
369 : if (size() != rhs.size()) {
370 : owned_data_.reset();
371 : owned_data_ = heap_alloc_if_necessary(rhs.size());
372 : }
373 : reset_pointer_vector(rhs.size());
374 : } else {
375 : ASSERT(rhs.size() == size(), "Must copy into same size, not "
376 : << rhs.size() << " into " << size());
377 : }
378 : if (LIKELY(data() != rhs.data())) {
379 : std::memcpy(data(), rhs.data(), size() * sizeof(value_type));
380 : }
381 : }
382 : return *this;
383 : }
384 :
385 : template <typename T, typename VectorType, size_t StaticSize>
386 : VectorImpl<T, VectorType, StaticSize>::VectorImpl(
387 : VectorImpl<T, VectorType, StaticSize>&& rhs) {
388 : owned_data_ = std::move(rhs.owned_data_);
389 : static_owned_data_ = std::move(rhs.static_owned_data_);
390 : **this = std::move(*rhs);
391 : owning_ = rhs.owning_;
392 : if (owning_) {
393 : reset_pointer_vector(size());
394 : } else {
395 : this->reset(data(), size());
396 : }
397 : rhs.clear();
398 : }
399 :
400 : template <typename T, typename VectorType, size_t StaticSize>
401 : VectorImpl<T, VectorType, StaticSize>&
402 : VectorImpl<T, VectorType, StaticSize>::operator=(
403 : VectorImpl<T, VectorType, StaticSize>&& rhs) {
404 : ASSERT(rhs.is_owning(),
405 : "Cannot move assign from a non-owning vector, because the correct "
406 : "behavior is unclear.");
407 : if (this != &rhs) {
408 : if (owning_) {
409 : owned_data_ = std::move(rhs.owned_data_);
410 : static_owned_data_ = std::move(rhs.static_owned_data_);
411 : **this = std::move(*rhs);
412 : reset_pointer_vector(size());
413 : rhs.clear();
414 : } else {
415 : ASSERT(rhs.size() == size(), "Must move into same size, not "
416 : << rhs.size() << " into " << size());
417 : if (LIKELY(data() != rhs.data())) {
418 : std::memcpy(data(), rhs.data(), size() * sizeof(value_type));
419 : rhs.clear();
420 : }
421 : }
422 : }
423 : return *this;
424 : }
425 :
426 : // This is a converting constructor. clang-tidy complains that it's not
427 : // explicit, but we want it to allow conversion.
428 : // clang-tidy: mark as explicit (we want conversion to VectorImpl)
429 : template <typename T, typename VectorType, size_t StaticSize>
430 : template <typename VT, bool VF,
431 : Requires<VectorImpl_detail::is_assignable_v<VectorType,
432 : typename VT::ResultType>>>
433 : VectorImpl<T, VectorType, StaticSize>::VectorImpl(
434 : const blaze::DenseVector<VT, VF>& expression) // NOLINT
435 : : owned_data_(heap_alloc_if_necessary((*expression).size())) {
436 : static_assert(
437 : VectorImpl_detail::is_assignable_v<VectorType, typename VT::ResultType>,
438 : "Cannot construct the VectorImpl type from the given expression type.");
439 : reset_pointer_vector((*expression).size());
440 : **this = expression;
441 : }
442 :
443 : template <typename T, typename VectorType, size_t StaticSize>
444 : template <typename VT, bool VF>
445 : VectorImpl<T, VectorType, StaticSize>&
446 : VectorImpl<T, VectorType, StaticSize>::operator=(
447 : const blaze::DenseVector<VT, VF>& expression) {
448 : static_assert(
449 : VectorImpl_detail::is_assignable_v<VectorType, typename VT::ResultType>,
450 : "Cannot assign to the VectorImpl type from the given expression type.");
451 : if (owning_ and (*expression).size() != size()) {
452 : owned_data_ = heap_alloc_if_necessary((*expression).size());
453 : reset_pointer_vector((*expression).size());
454 : } else if (not owning_) {
455 : ASSERT((*expression).size() == size(), "Must assign into same size, not "
456 : << (*expression).size()
457 : << " into " << size());
458 : }
459 : **this = expression;
460 : return *this;
461 : }
462 : /// \endcond
463 :
464 : // The case of assigning a type apart from the same VectorImpl or a
465 : // `blaze::DenseVector` forwards the assignment to the `blaze::CustomVector`
466 : // base type. In the case of a single compatible value, this fills the vector
467 : // with that value.
468 : template <typename T, typename VectorType, size_t StaticSize>
469 : VectorImpl<T, VectorType, StaticSize>&
470 : VectorImpl<T, VectorType, StaticSize>::operator=(const T& rhs) {
471 : **this = rhs;
472 : return *this;
473 : }
474 :
475 : template <typename T, typename VectorType, size_t StaticSize>
476 : void VectorImpl<T, VectorType, StaticSize>::clear() {
477 : BaseType::clear();
478 : owning_ = true;
479 : owned_data_.reset();
480 : // The state of static_owned_data_ doesn't matter.
481 : }
482 :
483 : template <typename T, typename VectorType, size_t StaticSize>
484 : void VectorImpl<T, VectorType, StaticSize>::pup(PUP::er& p) { // NOLINT
485 : if (not owning_ and p.isSizing()) {
486 : return;
487 : }
488 : ASSERT(owning_, "Cannot pup a non-owning vector!");
489 : auto my_size = size();
490 : p | my_size;
491 : if (my_size > 0) {
492 : if (p.isUnpacking()) {
493 : owning_ = true;
494 : owned_data_ = heap_alloc_if_necessary(my_size);
495 : reset_pointer_vector(my_size);
496 : }
497 : PUParray(p, data(), size());
498 : }
499 : }
500 :
501 : /// Output operator for VectorImpl
502 : template <typename T, typename VectorType, size_t StaticSize>
503 1 : std::ostream& operator<<(std::ostream& os,
504 : const VectorImpl<T, VectorType, StaticSize>& d) {
505 : sequence_print_helper(os, d.begin(), d.end());
506 : return os;
507 : }
508 :
509 0 : #define DECLARE_GENERAL_VECTOR_BLAZE_TRAITS(VECTOR_TYPE) \
510 : template <> \
511 : struct IsDenseVector<VECTOR_TYPE> : public blaze::TrueType {}; \
512 : \
513 : template <> \
514 : struct IsVector<VECTOR_TYPE> : public blaze::TrueType {}; \
515 : \
516 : template <> \
517 : struct CustomTransposeType<VECTOR_TYPE> { \
518 : using Type = VECTOR_TYPE; \
519 : }
520 :
521 : /*!
522 : * \ingroup DataStructuresGroup
523 : * \brief Instructs Blaze to provide the appropriate vector result type after
524 : * math operations. This is accomplished by specializing Blaze's type traits
525 : * that are used for handling return type deduction and specifying the `using
526 : * Type =` nested type alias in the traits.
527 : *
528 : * \param VECTOR_TYPE The vector type, which matches the type of the operation
529 : * result (e.g. `DataVector`)
530 : *
531 : * \param BLAZE_MATH_TRAIT The blaze trait/expression for which you want to
532 : * specify the return type (e.g. `AddTrait`).
533 : */
534 1 : #define BLAZE_TRAIT_SPECIALIZE_BINARY_TRAIT(VECTOR_TYPE, BLAZE_MATH_TRAIT) \
535 : template <> \
536 : struct BLAZE_MATH_TRAIT<VECTOR_TYPE, VECTOR_TYPE> { \
537 : using Type = VECTOR_TYPE; \
538 : }; \
539 : template <> \
540 : struct BLAZE_MATH_TRAIT<VECTOR_TYPE, VECTOR_TYPE::value_type> { \
541 : using Type = VECTOR_TYPE; \
542 : }; \
543 : template <> \
544 : struct BLAZE_MATH_TRAIT<VECTOR_TYPE::value_type, VECTOR_TYPE> { \
545 : using Type = VECTOR_TYPE; \
546 : }
547 :
548 : /*!
549 : * \ingroup DataStructuresGroup
550 : * \brief Instructs Blaze to provide the appropriate vector result type of an
551 : * operator between `VECTOR_TYPE` and `COMPATIBLE`, where the operation is
552 : * represented by `BLAZE_MATH_TRAIT`
553 : *
554 : * \param VECTOR_TYPE The vector type, which matches the type of the operation
555 : * result (e.g. `ComplexDataVector`)
556 : *
557 : * \param COMPATIBLE the type for which you want math operations to work with
558 : * `VECTOR_TYPE` smoothly (e.g. `DataVector`)
559 : *
560 : * \param BLAZE_MATH_TRAIT The blaze trait for which you want declare the Type
561 : * field (e.g. `AddTrait`)
562 : *
563 : * \param RESULT_TYPE The type which should be used as the 'return' type for the
564 : * binary operation
565 : */
566 : #define BLAZE_TRAIT_SPECIALIZE_COMPATIBLE_BINARY_TRAIT( \
567 1 : VECTOR_TYPE, COMPATIBLE, BLAZE_MATH_TRAIT, RESULT_TYPE) \
568 : template <> \
569 : struct BLAZE_MATH_TRAIT<VECTOR_TYPE, COMPATIBLE> { \
570 : using Type = RESULT_TYPE; \
571 : }; \
572 : template <> \
573 : struct BLAZE_MATH_TRAIT<COMPATIBLE, VECTOR_TYPE> { \
574 : using Type = RESULT_TYPE; \
575 : }
576 :
577 : /*!
578 : * \ingroup DataStructuresGroup
579 : * \brief Instructs Blaze to provide the appropriate vector result type of
580 : * arithmetic operations for `VECTOR_TYPE`. This is accomplished by specializing
581 : * Blaze's type traits that are used for handling return type deduction.
582 : *
583 : * \details Type definitions here are suitable for contiguous data
584 : * (e.g. `DataVector`), but this macro might need to be tweaked for other types
585 : * of data, for instance Fourier coefficients.
586 : *
587 : * \param VECTOR_TYPE The vector type, which for the arithmetic operations is
588 : * the type of the operation result (e.g. `DataVector`)
589 : */
590 1 : #define VECTOR_BLAZE_TRAIT_SPECIALIZE_ARITHMETIC_TRAITS(VECTOR_TYPE) \
591 : template <> \
592 : struct TransposeFlag<VECTOR_TYPE> \
593 : : BoolConstant<VECTOR_TYPE::transpose_flag> {}; \
594 : BLAZE_TRAIT_SPECIALIZE_BINARY_TRAIT(VECTOR_TYPE, AddTrait); \
595 : BLAZE_TRAIT_SPECIALIZE_BINARY_TRAIT(VECTOR_TYPE, SubTrait); \
596 : BLAZE_TRAIT_SPECIALIZE_BINARY_TRAIT(VECTOR_TYPE, MultTrait); \
597 : BLAZE_TRAIT_SPECIALIZE_BINARY_TRAIT(VECTOR_TYPE, DivTrait)
598 :
599 : /*!
600 : * \ingroup DataStructuresGroup
601 : * \brief Instructs Blaze to provide the appropriate vector result type of `Map`
602 : * operations (unary and binary) acting on `VECTOR_TYPE`. This is accomplished
603 : * by specializing Blaze's type traits that are used for handling return type
604 : * deduction.
605 : *
606 : * \details Type declarations here are suitable for contiguous data (e.g.
607 : * `DataVector`), but this macro might need to be tweaked for other types of
608 : * data, for instance Fourier coefficients.
609 : *
610 : * \param VECTOR_TYPE The vector type, which for the `Map` operations is
611 : * the type of the operation result (e.g. `DataVector`)
612 : */
613 1 : #define VECTOR_BLAZE_TRAIT_SPECIALIZE_ALL_MAP_TRAITS(VECTOR_TYPE) \
614 : template <typename Operator> \
615 : struct MapTrait<VECTOR_TYPE, Operator> { \
616 : using Type = VECTOR_TYPE; \
617 : }; \
618 : template <typename Operator> \
619 : struct MapTrait<VECTOR_TYPE, VECTOR_TYPE, Operator> { \
620 : using Type = VECTOR_TYPE; \
621 : }
622 :
623 : /*!
624 : * \ingroup DataStructuresGroup
625 : * \brief Defines the set of binary operations often supported for
626 : * `std::array<VECTOR_TYPE, size>`, for arbitrary `size`.
627 : *
628 : * \param VECTOR_TYPE The vector type (e.g. `DataVector`)
629 : */
630 1 : #define MAKE_STD_ARRAY_VECTOR_BINOPS(VECTOR_TYPE) \
631 : DEFINE_STD_ARRAY_BINOP(VECTOR_TYPE, VECTOR_TYPE::value_type, \
632 : VECTOR_TYPE, operator+, std::plus<>()) \
633 : DEFINE_STD_ARRAY_BINOP(VECTOR_TYPE, VECTOR_TYPE, \
634 : VECTOR_TYPE::value_type, operator+, std::plus<>()) \
635 : DEFINE_STD_ARRAY_BINOP(VECTOR_TYPE, VECTOR_TYPE, VECTOR_TYPE, operator+, \
636 : std::plus<>()) \
637 : \
638 : DEFINE_STD_ARRAY_BINOP(VECTOR_TYPE, VECTOR_TYPE::value_type, \
639 : VECTOR_TYPE, operator-, std::minus<>()) \
640 : DEFINE_STD_ARRAY_BINOP(VECTOR_TYPE, VECTOR_TYPE, \
641 : VECTOR_TYPE::value_type, operator-, std::minus<>()) \
642 : DEFINE_STD_ARRAY_BINOP(VECTOR_TYPE, VECTOR_TYPE, VECTOR_TYPE, operator-, \
643 : std::minus<>()) \
644 : \
645 : DEFINE_STD_ARRAY_INPLACE_BINOP(VECTOR_TYPE, VECTOR_TYPE, operator-=, \
646 : std::minus<>()) \
647 : DEFINE_STD_ARRAY_INPLACE_BINOP( \
648 : VECTOR_TYPE, VECTOR_TYPE::value_type, operator-=, std::minus<>()) \
649 : DEFINE_STD_ARRAY_INPLACE_BINOP(VECTOR_TYPE, VECTOR_TYPE, operator+=, \
650 : std::plus<>()) \
651 : DEFINE_STD_ARRAY_INPLACE_BINOP( \
652 : VECTOR_TYPE, VECTOR_TYPE::value_type, operator+=, std::plus<>())
653 :
654 : /*!
655 : * \ingroup DataStructuresGroup
656 : * \brief Defines the `MakeWithValueImpl` `apply` specialization
657 : *
658 : * \details The `MakeWithValueImpl<VECTOR_TYPE, VECTOR_TYPE>` member
659 : * `apply(VECTOR_TYPE, VECTOR_TYPE::value_type)` specialization defined by this
660 : * macro produces an object with the same size as the `input` argument,
661 : * initialized with the `value` argument in every entry.
662 : *
663 : * \param VECTOR_TYPE The vector type (e.g. `DataVector`)
664 : */
665 1 : #define MAKE_WITH_VALUE_IMPL_DEFINITION_FOR(VECTOR_TYPE) \
666 : namespace MakeWithValueImpls { \
667 : template <> \
668 : struct NumberOfPoints<VECTOR_TYPE> { \
669 : static SPECTRE_ALWAYS_INLINE size_t apply(const VECTOR_TYPE& input) { \
670 : return input.size(); \
671 : } \
672 : }; \
673 : template <> \
674 : struct MakeWithSize<VECTOR_TYPE> { \
675 : static SPECTRE_ALWAYS_INLINE VECTOR_TYPE \
676 : apply(const size_t size, const VECTOR_TYPE::value_type value) { \
677 : return VECTOR_TYPE(size, value); \
678 : } \
679 : }; \
680 : } /* namespace MakeWithValueImpls */ \
681 : template <> \
682 : struct SetNumberOfGridPointsImpls::SetNumberOfGridPointsImpl<VECTOR_TYPE> { \
683 : static constexpr bool is_trivial = false; \
684 : static SPECTRE_ALWAYS_INLINE void apply( \
685 : const gsl::not_null<VECTOR_TYPE*> result, const size_t size) { \
686 : result->destructive_resize(size); \
687 : } \
688 : };
689 :
690 : /// @{
691 : /*!
692 : * \ingroup DataStructuresGroup
693 : * \ingroup TypeTraitsGroup
694 : * \brief Helper struct to determine the element type of a VectorImpl or
695 : * container of VectorImpl
696 : *
697 : * \details Extracts the element type of a `VectorImpl`, a std::array of
698 : * `VectorImpl`, or a reference or pointer to a `VectorImpl`. In any of these
699 : * cases, the `type` member is defined as the `ElementType` of the `VectorImpl`
700 : * in question. If, instead, `get_vector_element_type` is passed an arithmetic
701 : * or complex arithemetic type, the `type` member is defined as the passed type.
702 : *
703 : * \snippet DataStructures/Test_VectorImpl.cpp get_vector_element_type_example
704 : */
705 : // cast to bool needed to avoid the compiler mistaking the type to be determined
706 : // by T
707 : template <typename T,
708 : bool = static_cast<bool>(tt::is_complex_of_fundamental_v<T> or
709 : std::is_fundamental_v<T>)>
710 1 : struct get_vector_element_type;
711 : template <typename T>
712 0 : struct get_vector_element_type<T, true> {
713 0 : using type = T;
714 : };
715 : template <typename T>
716 0 : struct get_vector_element_type<const T, false> {
717 0 : using type = typename get_vector_element_type<T>::type;
718 : };
719 : template <typename T>
720 0 : struct get_vector_element_type<T, false> {
721 0 : using type = typename get_vector_element_type<
722 : typename T::ResultType::ElementType>::type;
723 : };
724 : template <typename T>
725 0 : struct get_vector_element_type<T*, false> {
726 0 : using type = typename get_vector_element_type<T>::type;
727 : };
728 : template <typename T>
729 0 : struct get_vector_element_type<T&, false> {
730 0 : using type = typename get_vector_element_type<T>::type;
731 : };
732 : template <typename T, size_t S>
733 : struct get_vector_element_type<std::array<T, S>, false> {
734 : using type = typename get_vector_element_type<T>::type;
735 : };
736 : /// @}
737 :
738 : template <typename T>
739 0 : using get_vector_element_type_t = typename get_vector_element_type<T>::type;
740 :
741 : namespace detail {
742 : template <typename T, typename VectorType, size_t StaticSize>
743 : std::true_type is_derived_of_vector_impl_impl(
744 : const VectorImpl<T, VectorType, StaticSize>*);
745 :
746 : std::false_type is_derived_of_vector_impl_impl(...);
747 : } // namespace detail
748 :
749 : /// \ingroup TypeTraitsGroup
750 : /// This is `std::true_type` if the provided type possesses an implicit
751 : /// conversion to any `VectorImpl`, which is the primary feature of SpECTRE
752 : /// vectors generally. Otherwise, it is `std::false_type`.
753 : template <typename T>
754 1 : using is_derived_of_vector_impl =
755 : decltype(detail::is_derived_of_vector_impl_impl(std::declval<T*>()));
756 :
757 : template <typename T>
758 0 : constexpr bool is_derived_of_vector_impl_v =
759 : is_derived_of_vector_impl<T>::value;
760 :
761 : // impose strict equality for derived classes of VectorImpl; note that this
762 : // overrides intended behavior in blaze for comparison operators to use
763 : // approximate equality in favor of equality between containers being
764 : // appropriately recursive. This form primarily works by using templates to
765 : // ensure that our comparison operator is resolved with higher priority than the
766 : // blaze form as of blaze 3.8
767 : template <
768 : typename Lhs, typename Rhs,
769 : Requires<(is_derived_of_vector_impl_v<Lhs> or
770 : is_derived_of_vector_impl_v<
771 : Rhs>)and not(std::is_base_of_v<blaze::Computation, Lhs> or
772 : std::is_base_of_v<blaze::Computation, Rhs>) and
773 : not(std::is_same_v<Rhs, typename Lhs::ElementType> or
774 : std::is_same_v<Lhs, typename Rhs::ElementType>)> = nullptr>
775 0 : bool operator==(const Lhs& lhs, const Rhs& rhs) {
776 : return blaze::equal<blaze::strict>(lhs, rhs);
777 : }
778 :
779 : template <
780 : typename Lhs, typename Rhs,
781 : Requires<(is_derived_of_vector_impl_v<Lhs> or
782 : is_derived_of_vector_impl_v<
783 : Rhs>)and not(std::is_base_of_v<blaze::Computation, Lhs> or
784 : std::is_base_of_v<blaze::Computation, Rhs>) and
785 : not(std::is_same_v<Rhs, typename Lhs::ElementType> or
786 : std::is_same_v<Lhs, typename Lhs::ElementType>)> = nullptr>
787 0 : bool operator!=(const Lhs& lhs, const Rhs& rhs) {
788 : return not(lhs == rhs);
789 : }
790 :
791 : // Impose strict equality for any expression templates; note that
792 : // this overrides intended behavior in blaze for comparison
793 : // operators to use approximate equality in favor of equality
794 : // between containers being appropriately recursive. This form
795 : // primarily works by using templates to ensure that our
796 : // comparison operator is resolved with higher priority than the
797 : // blaze form as of blaze 3.8
798 : template <typename Lhs, typename Rhs,
799 : Requires<std::is_base_of_v<blaze::Computation, Lhs> or
800 : std::is_base_of_v<blaze::Computation, Rhs>> = nullptr>
801 : bool operator==(const Lhs& lhs, const Rhs& rhs) {
802 : return blaze::equal<blaze::strict>(lhs, rhs);
803 : }
804 :
805 : template <typename Lhs, typename Rhs,
806 : Requires<std::is_base_of_v<blaze::Computation, Lhs> or
807 : std::is_base_of_v<blaze::Computation, Rhs>> = nullptr>
808 : bool operator!=(const Lhs& lhs, const Rhs& rhs) {
809 : return not(lhs == rhs);
810 : }
811 :
812 : template <typename Lhs, Requires<is_derived_of_vector_impl_v<Lhs>> = nullptr>
813 0 : bool operator==(const Lhs& lhs, const typename Lhs::ElementType& rhs) {
814 : for (const auto& element : lhs) {
815 : if (element != rhs) {
816 : return false;
817 : }
818 : }
819 : return true;
820 : }
821 :
822 : template <typename Lhs, Requires<is_derived_of_vector_impl_v<Lhs>> = nullptr>
823 0 : bool operator!=(const Lhs& lhs, const typename Lhs::ElementType& rhs) {
824 : return not(lhs == rhs);
825 : }
826 :
827 : template <typename Rhs, Requires<is_derived_of_vector_impl_v<Rhs>> = nullptr>
828 0 : bool operator==(const typename Rhs::ElementType& lhs, const Rhs& rhs) {
829 : return rhs == lhs;
830 : }
831 :
832 : template <typename Rhs, Requires<is_derived_of_vector_impl_v<Rhs>> = nullptr>
833 0 : bool operator!=(const typename Rhs::ElementType& lhs, const Rhs& rhs) {
834 : return not(lhs == rhs);
835 : }
836 :
837 : /// \ingroup DataStructuresGroup
838 : /// Make the input `view` a `const` view of the const data `vector`, at
839 : /// offset `offset` and length `extent`.
840 : ///
841 : /// \warning This DOES modify the (const) input `view`. The reason `view` is
842 : /// taken by const pointer is to try to insist that the object to be a `const`
843 : /// view is actually const. Of course, there are ways of subverting this
844 : /// intended functionality and editing the data pointed into by `view` after
845 : /// this function is called; doing so is highly discouraged and results in
846 : /// undefined behavior.
847 : template <typename VectorType,
848 : Requires<is_derived_of_vector_impl_v<VectorType>> = nullptr>
849 1 : void make_const_view(const gsl::not_null<const VectorType*> view,
850 : const VectorType& vector, const size_t offset,
851 : const size_t extent) {
852 : const_cast<VectorType*>(view.get()) // NOLINT
853 : ->set_data_ref(
854 : const_cast<typename VectorType::value_type*>(vector.data()) // NOLINT
855 : + offset, // NOLINT
856 : extent);
857 : }
858 :
859 : template <typename T, typename VectorType, size_t StaticSize>
860 0 : inline bool contains_allocations(
861 : const VectorImpl<T, VectorType, StaticSize>& value) {
862 : return value.size() > StaticSize and value.is_owning();
863 : }
|