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 either `SPECTRE_DEBUG` or `SPECTRE_NAN_INIT` are defined, then the
134 : * `VectorImpl` is default initialized to `signaling_NaN()`. Otherwise, the
135 : * vector is filled with 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_DEBUG) || defined(SPECTRE_NAN_INIT)
188 : std::fill(data(), data() + set_size,
189 : std::numeric_limits<value_type>::signaling_NaN());
190 : #endif // SPECTRE_DEBUG
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 : return;
329 : }
330 : if (owned_data_ == nullptr and set_size > StaticSize) {
331 : ERROR(
332 : "VectorImpl::reset_pointer_vector cannot be called when owned_data_ "
333 : "is nullptr.");
334 : }
335 :
336 : if (set_size <= StaticSize) {
337 : this->reset(static_owned_data_.data(), set_size);
338 : // Free memory if downsizing
339 : owned_data_ = nullptr;
340 : } else {
341 : this->reset(owned_data_.get(), set_size);
342 : }
343 : }
344 :
345 0 : SPECTRE_ALWAYS_INLINE std::unique_ptr<value_type[]> heap_alloc_if_necessary(
346 : const size_t set_size) {
347 : return set_size > StaticSize
348 : ? cpp20::make_unique_for_overwrite<value_type[]>(set_size)
349 : : nullptr;
350 : }
351 : };
352 :
353 : /// \cond HIDDEN_SYMBOLS
354 : template <typename T, typename VectorType, size_t StaticSize>
355 : VectorImpl<T, VectorType, StaticSize>::VectorImpl(
356 : const VectorImpl<T, VectorType, StaticSize>& rhs)
357 : : BaseType{rhs}, owned_data_(heap_alloc_if_necessary(rhs.size())) {
358 : reset_pointer_vector(rhs.size());
359 : std::memcpy(data(), rhs.data(), size() * sizeof(value_type));
360 : }
361 :
362 : template <typename T, typename VectorType, size_t StaticSize>
363 : VectorImpl<T, VectorType, StaticSize>&
364 : VectorImpl<T, VectorType, StaticSize>::operator=(
365 : const VectorImpl<T, VectorType, StaticSize>& rhs) {
366 : if (this != &rhs) {
367 : if (owning_) {
368 : if (size() != rhs.size()) {
369 : owned_data_.reset();
370 : owned_data_ = heap_alloc_if_necessary(rhs.size());
371 : }
372 : reset_pointer_vector(rhs.size());
373 : } else {
374 : ASSERT(rhs.size() == size(), "Must copy into same size, not "
375 : << rhs.size() << " into " << size());
376 : }
377 : if (LIKELY(data() != rhs.data())) {
378 : std::memcpy(data(), rhs.data(), size() * sizeof(value_type));
379 : }
380 : }
381 : return *this;
382 : }
383 :
384 : template <typename T, typename VectorType, size_t StaticSize>
385 : VectorImpl<T, VectorType, StaticSize>::VectorImpl(
386 : VectorImpl<T, VectorType, StaticSize>&& rhs) {
387 : owned_data_ = std::move(rhs.owned_data_);
388 : static_owned_data_ = std::move(rhs.static_owned_data_);
389 : **this = std::move(*rhs);
390 : owning_ = rhs.owning_;
391 : if (owning_) {
392 : reset_pointer_vector(size());
393 : } else {
394 : this->reset(data(), size());
395 : }
396 : rhs.clear();
397 : }
398 :
399 : template <typename T, typename VectorType, size_t StaticSize>
400 : VectorImpl<T, VectorType, StaticSize>&
401 : VectorImpl<T, VectorType, StaticSize>::operator=(
402 : VectorImpl<T, VectorType, StaticSize>&& rhs) {
403 : ASSERT(rhs.is_owning(),
404 : "Cannot move assign from a non-owning vector, because the correct "
405 : "behavior is unclear.");
406 : if (this != &rhs) {
407 : if (owning_) {
408 : owned_data_ = std::move(rhs.owned_data_);
409 : static_owned_data_ = std::move(rhs.static_owned_data_);
410 : **this = std::move(*rhs);
411 : reset_pointer_vector(size());
412 : rhs.clear();
413 : } else {
414 : ASSERT(rhs.size() == size(), "Must move into same size, not "
415 : << rhs.size() << " into " << size());
416 : if (LIKELY(data() != rhs.data())) {
417 : std::memcpy(data(), rhs.data(), size() * sizeof(value_type));
418 : rhs.clear();
419 : }
420 : }
421 : }
422 : return *this;
423 : }
424 :
425 : // This is a converting constructor. clang-tidy complains that it's not
426 : // explicit, but we want it to allow conversion.
427 : // clang-tidy: mark as explicit (we want conversion to VectorImpl)
428 : template <typename T, typename VectorType, size_t StaticSize>
429 : template <typename VT, bool VF,
430 : Requires<VectorImpl_detail::is_assignable_v<VectorType,
431 : typename VT::ResultType>>>
432 : VectorImpl<T, VectorType, StaticSize>::VectorImpl(
433 : const blaze::DenseVector<VT, VF>& expression) // NOLINT
434 : : owned_data_(heap_alloc_if_necessary((*expression).size())) {
435 : static_assert(
436 : VectorImpl_detail::is_assignable_v<VectorType, typename VT::ResultType>,
437 : "Cannot construct the VectorImpl type from the given expression type.");
438 : reset_pointer_vector((*expression).size());
439 : **this = expression;
440 : }
441 :
442 : template <typename T, typename VectorType, size_t StaticSize>
443 : template <typename VT, bool VF>
444 : VectorImpl<T, VectorType, StaticSize>&
445 : VectorImpl<T, VectorType, StaticSize>::operator=(
446 : const blaze::DenseVector<VT, VF>& expression) {
447 : static_assert(
448 : VectorImpl_detail::is_assignable_v<VectorType, typename VT::ResultType>,
449 : "Cannot assign to the VectorImpl type from the given expression type.");
450 : if (owning_ and (*expression).size() != size()) {
451 : owned_data_ = heap_alloc_if_necessary((*expression).size());
452 : reset_pointer_vector((*expression).size());
453 : } else if (not owning_) {
454 : ASSERT((*expression).size() == size(), "Must assign into same size, not "
455 : << (*expression).size()
456 : << " into " << size());
457 : }
458 : **this = expression;
459 : return *this;
460 : }
461 : /// \endcond
462 :
463 : // The case of assigning a type apart from the same VectorImpl or a
464 : // `blaze::DenseVector` forwards the assignment to the `blaze::CustomVector`
465 : // base type. In the case of a single compatible value, this fills the vector
466 : // with that value.
467 : template <typename T, typename VectorType, size_t StaticSize>
468 : VectorImpl<T, VectorType, StaticSize>&
469 : VectorImpl<T, VectorType, StaticSize>::operator=(const T& rhs) {
470 : **this = rhs;
471 : return *this;
472 : }
473 :
474 : template <typename T, typename VectorType, size_t StaticSize>
475 : void VectorImpl<T, VectorType, StaticSize>::clear() {
476 : BaseType::clear();
477 : owning_ = true;
478 : owned_data_.reset();
479 : // The state of static_owned_data_ doesn't matter.
480 : }
481 :
482 : template <typename T, typename VectorType, size_t StaticSize>
483 : void VectorImpl<T, VectorType, StaticSize>::pup(PUP::er& p) { // NOLINT
484 : if (not owning_ and p.isSizing()) {
485 : return;
486 : }
487 : ASSERT(owning_, "Cannot pup a non-owning vector!");
488 : auto my_size = size();
489 : p | my_size;
490 : if (my_size > 0) {
491 : if (p.isUnpacking()) {
492 : owning_ = true;
493 : owned_data_ = heap_alloc_if_necessary(my_size);
494 : reset_pointer_vector(my_size);
495 : }
496 : PUParray(p, data(), size());
497 : }
498 : }
499 :
500 : /// Output operator for VectorImpl
501 : template <typename T, typename VectorType, size_t StaticSize>
502 1 : std::ostream& operator<<(std::ostream& os,
503 : const VectorImpl<T, VectorType, StaticSize>& d) {
504 : sequence_print_helper(os, d.begin(), d.end());
505 : return os;
506 : }
507 :
508 0 : #define DECLARE_GENERAL_VECTOR_BLAZE_TRAITS(VECTOR_TYPE) \
509 : template <> \
510 : struct IsDenseVector<VECTOR_TYPE> : public blaze::TrueType {}; \
511 : \
512 : template <> \
513 : struct IsVector<VECTOR_TYPE> : public blaze::TrueType {}; \
514 : \
515 : template <> \
516 : struct CustomTransposeType<VECTOR_TYPE> { \
517 : using Type = VECTOR_TYPE; \
518 : }
519 :
520 : /*!
521 : * \ingroup DataStructuresGroup
522 : * \brief Instructs Blaze to provide the appropriate vector result type after
523 : * math operations. This is accomplished by specializing Blaze's type traits
524 : * that are used for handling return type deduction and specifying the `using
525 : * Type =` nested type alias in the traits.
526 : *
527 : * \param VECTOR_TYPE The vector type, which matches the type of the operation
528 : * result (e.g. `DataVector`)
529 : *
530 : * \param BLAZE_MATH_TRAIT The blaze trait/expression for which you want to
531 : * specify the return type (e.g. `AddTrait`).
532 : */
533 1 : #define BLAZE_TRAIT_SPECIALIZE_BINARY_TRAIT(VECTOR_TYPE, BLAZE_MATH_TRAIT) \
534 : template <> \
535 : struct BLAZE_MATH_TRAIT<VECTOR_TYPE, VECTOR_TYPE> { \
536 : using Type = VECTOR_TYPE; \
537 : }; \
538 : template <> \
539 : struct BLAZE_MATH_TRAIT<VECTOR_TYPE, VECTOR_TYPE::value_type> { \
540 : using Type = VECTOR_TYPE; \
541 : }; \
542 : template <> \
543 : struct BLAZE_MATH_TRAIT<VECTOR_TYPE::value_type, VECTOR_TYPE> { \
544 : using Type = VECTOR_TYPE; \
545 : }
546 :
547 : /*!
548 : * \ingroup DataStructuresGroup
549 : * \brief Instructs Blaze to provide the appropriate vector result type of an
550 : * operator between `VECTOR_TYPE` and `COMPATIBLE`, where the operation is
551 : * represented by `BLAZE_MATH_TRAIT`
552 : *
553 : * \param VECTOR_TYPE The vector type, which matches the type of the operation
554 : * result (e.g. `ComplexDataVector`)
555 : *
556 : * \param COMPATIBLE the type for which you want math operations to work with
557 : * `VECTOR_TYPE` smoothly (e.g. `DataVector`)
558 : *
559 : * \param BLAZE_MATH_TRAIT The blaze trait for which you want declare the Type
560 : * field (e.g. `AddTrait`)
561 : *
562 : * \param RESULT_TYPE The type which should be used as the 'return' type for the
563 : * binary operation
564 : */
565 : #define BLAZE_TRAIT_SPECIALIZE_COMPATIBLE_BINARY_TRAIT( \
566 1 : VECTOR_TYPE, COMPATIBLE, BLAZE_MATH_TRAIT, RESULT_TYPE) \
567 : template <> \
568 : struct BLAZE_MATH_TRAIT<VECTOR_TYPE, COMPATIBLE> { \
569 : using Type = RESULT_TYPE; \
570 : }; \
571 : template <> \
572 : struct BLAZE_MATH_TRAIT<COMPATIBLE, VECTOR_TYPE> { \
573 : using Type = RESULT_TYPE; \
574 : }
575 :
576 : /*!
577 : * \ingroup DataStructuresGroup
578 : * \brief Instructs Blaze to provide the appropriate vector result type of
579 : * arithmetic operations for `VECTOR_TYPE`. This is accomplished by specializing
580 : * Blaze's type traits that are used for handling return type deduction.
581 : *
582 : * \details Type definitions here are suitable for contiguous data
583 : * (e.g. `DataVector`), but this macro might need to be tweaked for other types
584 : * of data, for instance Fourier coefficients.
585 : *
586 : * \param VECTOR_TYPE The vector type, which for the arithmetic operations is
587 : * the type of the operation result (e.g. `DataVector`)
588 : */
589 1 : #define VECTOR_BLAZE_TRAIT_SPECIALIZE_ARITHMETIC_TRAITS(VECTOR_TYPE) \
590 : template <> \
591 : struct TransposeFlag<VECTOR_TYPE> \
592 : : BoolConstant<VECTOR_TYPE::transpose_flag> {}; \
593 : BLAZE_TRAIT_SPECIALIZE_BINARY_TRAIT(VECTOR_TYPE, AddTrait); \
594 : BLAZE_TRAIT_SPECIALIZE_BINARY_TRAIT(VECTOR_TYPE, SubTrait); \
595 : BLAZE_TRAIT_SPECIALIZE_BINARY_TRAIT(VECTOR_TYPE, MultTrait); \
596 : BLAZE_TRAIT_SPECIALIZE_BINARY_TRAIT(VECTOR_TYPE, DivTrait)
597 :
598 : /*!
599 : * \ingroup DataStructuresGroup
600 : * \brief Instructs Blaze to provide the appropriate vector result type of `Map`
601 : * operations (unary and binary) acting on `VECTOR_TYPE`. This is accomplished
602 : * by specializing Blaze's type traits that are used for handling return type
603 : * deduction.
604 : *
605 : * \details Type declarations here are suitable for contiguous data (e.g.
606 : * `DataVector`), but this macro might need to be tweaked for other types of
607 : * data, for instance Fourier coefficients.
608 : *
609 : * \param VECTOR_TYPE The vector type, which for the `Map` operations is
610 : * the type of the operation result (e.g. `DataVector`)
611 : */
612 1 : #define VECTOR_BLAZE_TRAIT_SPECIALIZE_ALL_MAP_TRAITS(VECTOR_TYPE) \
613 : template <typename Operator> \
614 : struct MapTrait<VECTOR_TYPE, Operator> { \
615 : using Type = VECTOR_TYPE; \
616 : }; \
617 : template <typename Operator> \
618 : struct MapTrait<VECTOR_TYPE, VECTOR_TYPE, Operator> { \
619 : using Type = VECTOR_TYPE; \
620 : }
621 :
622 : /*!
623 : * \ingroup DataStructuresGroup
624 : * \brief Defines the set of binary operations often supported for
625 : * `std::array<VECTOR_TYPE, size>`, for arbitrary `size`.
626 : *
627 : * \param VECTOR_TYPE The vector type (e.g. `DataVector`)
628 : */
629 1 : #define MAKE_STD_ARRAY_VECTOR_BINOPS(VECTOR_TYPE) \
630 : DEFINE_STD_ARRAY_BINOP(VECTOR_TYPE, VECTOR_TYPE::value_type, \
631 : VECTOR_TYPE, operator+, std::plus<>()) \
632 : DEFINE_STD_ARRAY_BINOP(VECTOR_TYPE, VECTOR_TYPE, \
633 : VECTOR_TYPE::value_type, operator+, std::plus<>()) \
634 : DEFINE_STD_ARRAY_BINOP(VECTOR_TYPE, VECTOR_TYPE, VECTOR_TYPE, operator+, \
635 : std::plus<>()) \
636 : \
637 : DEFINE_STD_ARRAY_BINOP(VECTOR_TYPE, VECTOR_TYPE::value_type, \
638 : VECTOR_TYPE, operator-, std::minus<>()) \
639 : DEFINE_STD_ARRAY_BINOP(VECTOR_TYPE, VECTOR_TYPE, \
640 : VECTOR_TYPE::value_type, operator-, std::minus<>()) \
641 : DEFINE_STD_ARRAY_BINOP(VECTOR_TYPE, VECTOR_TYPE, VECTOR_TYPE, operator-, \
642 : std::minus<>()) \
643 : \
644 : DEFINE_STD_ARRAY_INPLACE_BINOP(VECTOR_TYPE, VECTOR_TYPE, operator-=, \
645 : std::minus<>()) \
646 : DEFINE_STD_ARRAY_INPLACE_BINOP( \
647 : VECTOR_TYPE, VECTOR_TYPE::value_type, operator-=, std::minus<>()) \
648 : DEFINE_STD_ARRAY_INPLACE_BINOP(VECTOR_TYPE, VECTOR_TYPE, operator+=, \
649 : std::plus<>()) \
650 : DEFINE_STD_ARRAY_INPLACE_BINOP( \
651 : VECTOR_TYPE, VECTOR_TYPE::value_type, operator+=, std::plus<>())
652 :
653 : /*!
654 : * \ingroup DataStructuresGroup
655 : * \brief Defines the `MakeWithValueImpl` `apply` specialization
656 : *
657 : * \details The `MakeWithValueImpl<VECTOR_TYPE, VECTOR_TYPE>` member
658 : * `apply(VECTOR_TYPE, VECTOR_TYPE::value_type)` specialization defined by this
659 : * macro produces an object with the same size as the `input` argument,
660 : * initialized with the `value` argument in every entry.
661 : *
662 : * \param VECTOR_TYPE The vector type (e.g. `DataVector`)
663 : */
664 1 : #define MAKE_WITH_VALUE_IMPL_DEFINITION_FOR(VECTOR_TYPE) \
665 : namespace MakeWithValueImpls { \
666 : template <> \
667 : struct NumberOfPoints<VECTOR_TYPE> { \
668 : static SPECTRE_ALWAYS_INLINE size_t apply(const VECTOR_TYPE& input) { \
669 : return input.size(); \
670 : } \
671 : }; \
672 : template <> \
673 : struct MakeWithSize<VECTOR_TYPE> { \
674 : static SPECTRE_ALWAYS_INLINE VECTOR_TYPE \
675 : apply(const size_t size, const VECTOR_TYPE::value_type value) { \
676 : return VECTOR_TYPE(size, value); \
677 : } \
678 : }; \
679 : } /* namespace MakeWithValueImpls */ \
680 : template <> \
681 : struct SetNumberOfGridPointsImpls::SetNumberOfGridPointsImpl<VECTOR_TYPE> { \
682 : static constexpr bool is_trivial = false; \
683 : static SPECTRE_ALWAYS_INLINE void apply( \
684 : const gsl::not_null<VECTOR_TYPE*> result, const size_t size) { \
685 : result->destructive_resize(size); \
686 : } \
687 : };
688 :
689 : /// @{
690 : /*!
691 : * \ingroup DataStructuresGroup
692 : * \ingroup TypeTraitsGroup
693 : * \brief Helper struct to determine the element type of a VectorImpl or
694 : * container of VectorImpl
695 : *
696 : * \details Extracts the element type of a `VectorImpl`, a std::array of
697 : * `VectorImpl`, or a reference or pointer to a `VectorImpl`. In any of these
698 : * cases, the `type` member is defined as the `ElementType` of the `VectorImpl`
699 : * in question. If, instead, `get_vector_element_type` is passed an arithmetic
700 : * or complex arithemetic type, the `type` member is defined as the passed type.
701 : *
702 : * \snippet DataStructures/Test_VectorImpl.cpp get_vector_element_type_example
703 : */
704 : // cast to bool needed to avoid the compiler mistaking the type to be determined
705 : // by T
706 : template <typename T,
707 : bool = static_cast<bool>(tt::is_complex_of_fundamental_v<T> or
708 : std::is_fundamental_v<T>)>
709 1 : struct get_vector_element_type;
710 : template <typename T>
711 0 : struct get_vector_element_type<T, true> {
712 0 : using type = T;
713 : };
714 : template <typename T>
715 0 : struct get_vector_element_type<const T, false> {
716 0 : using type = typename get_vector_element_type<T>::type;
717 : };
718 : template <typename T>
719 0 : struct get_vector_element_type<T, false> {
720 0 : using type = typename get_vector_element_type<
721 : typename T::ResultType::ElementType>::type;
722 : };
723 : template <typename T>
724 0 : struct get_vector_element_type<T*, false> {
725 0 : using type = typename get_vector_element_type<T>::type;
726 : };
727 : template <typename T>
728 0 : struct get_vector_element_type<T&, false> {
729 0 : using type = typename get_vector_element_type<T>::type;
730 : };
731 : template <typename T, size_t S>
732 : struct get_vector_element_type<std::array<T, S>, false> {
733 : using type = typename get_vector_element_type<T>::type;
734 : };
735 : /// @}
736 :
737 : template <typename T>
738 0 : using get_vector_element_type_t = typename get_vector_element_type<T>::type;
739 :
740 : namespace detail {
741 : template <typename T, typename VectorType, size_t StaticSize>
742 : std::true_type is_derived_of_vector_impl_impl(
743 : const VectorImpl<T, VectorType, StaticSize>*);
744 :
745 : std::false_type is_derived_of_vector_impl_impl(...);
746 : } // namespace detail
747 :
748 : /// \ingroup TypeTraitsGroup
749 : /// This is `std::true_type` if the provided type possesses an implicit
750 : /// conversion to any `VectorImpl`, which is the primary feature of SpECTRE
751 : /// vectors generally. Otherwise, it is `std::false_type`.
752 : template <typename T>
753 1 : using is_derived_of_vector_impl =
754 : decltype(detail::is_derived_of_vector_impl_impl(std::declval<T*>()));
755 :
756 : template <typename T>
757 0 : constexpr bool is_derived_of_vector_impl_v =
758 : is_derived_of_vector_impl<T>::value;
759 :
760 : // impose strict equality for derived classes of VectorImpl; note that this
761 : // overrides intended behavior in blaze for comparison operators to use
762 : // approximate equality in favor of equality between containers being
763 : // appropriately recursive. This form primarily works by using templates to
764 : // ensure that our comparison operator is resolved with higher priority than the
765 : // blaze form as of blaze 3.8
766 : template <
767 : typename Lhs, typename Rhs,
768 : Requires<(is_derived_of_vector_impl_v<Lhs> or
769 : is_derived_of_vector_impl_v<
770 : Rhs>)and not(std::is_base_of_v<blaze::Computation, Lhs> or
771 : std::is_base_of_v<blaze::Computation, Rhs>) and
772 : not(std::is_same_v<Rhs, typename Lhs::ElementType> or
773 : std::is_same_v<Lhs, typename Rhs::ElementType>)> = nullptr>
774 0 : bool operator==(const Lhs& lhs, const Rhs& rhs) {
775 : return blaze::equal<blaze::strict>(lhs, rhs);
776 : }
777 :
778 : template <
779 : typename Lhs, typename Rhs,
780 : Requires<(is_derived_of_vector_impl_v<Lhs> or
781 : is_derived_of_vector_impl_v<
782 : Rhs>)and not(std::is_base_of_v<blaze::Computation, Lhs> or
783 : std::is_base_of_v<blaze::Computation, Rhs>) and
784 : not(std::is_same_v<Rhs, typename Lhs::ElementType> or
785 : std::is_same_v<Lhs, typename Lhs::ElementType>)> = nullptr>
786 0 : bool operator!=(const Lhs& lhs, const Rhs& rhs) {
787 : return not(lhs == rhs);
788 : }
789 :
790 : // Impose strict equality for any expression templates; note that
791 : // this overrides intended behavior in blaze for comparison
792 : // operators to use approximate equality in favor of equality
793 : // between containers being appropriately recursive. This form
794 : // primarily works by using templates to ensure that our
795 : // comparison operator is resolved with higher priority than the
796 : // blaze form as of blaze 3.8
797 : template <typename Lhs, typename Rhs,
798 : Requires<std::is_base_of_v<blaze::Computation, Lhs> or
799 : std::is_base_of_v<blaze::Computation, Rhs>> = nullptr>
800 : bool operator==(const Lhs& lhs, const Rhs& rhs) {
801 : return blaze::equal<blaze::strict>(lhs, rhs);
802 : }
803 :
804 : template <typename Lhs, typename Rhs,
805 : Requires<std::is_base_of_v<blaze::Computation, Lhs> or
806 : std::is_base_of_v<blaze::Computation, Rhs>> = nullptr>
807 : bool operator!=(const Lhs& lhs, const Rhs& rhs) {
808 : return not(lhs == rhs);
809 : }
810 :
811 : template <typename Lhs, Requires<is_derived_of_vector_impl_v<Lhs>> = nullptr>
812 0 : bool operator==(const Lhs& lhs, const typename Lhs::ElementType& rhs) {
813 : for (const auto& element : lhs) {
814 : if (element != rhs) {
815 : return false;
816 : }
817 : }
818 : return true;
819 : }
820 :
821 : template <typename Lhs, Requires<is_derived_of_vector_impl_v<Lhs>> = nullptr>
822 0 : bool operator!=(const Lhs& lhs, const typename Lhs::ElementType& rhs) {
823 : return not(lhs == rhs);
824 : }
825 :
826 : template <typename Rhs, Requires<is_derived_of_vector_impl_v<Rhs>> = nullptr>
827 0 : bool operator==(const typename Rhs::ElementType& lhs, const Rhs& rhs) {
828 : return rhs == lhs;
829 : }
830 :
831 : template <typename Rhs, Requires<is_derived_of_vector_impl_v<Rhs>> = nullptr>
832 0 : bool operator!=(const typename Rhs::ElementType& lhs, const Rhs& rhs) {
833 : return not(lhs == rhs);
834 : }
835 :
836 : /// \ingroup DataStructuresGroup
837 : /// Make the input `view` a `const` view of the const data `vector`, at
838 : /// offset `offset` and length `extent`.
839 : ///
840 : /// \warning This DOES modify the (const) input `view`. The reason `view` is
841 : /// taken by const pointer is to try to insist that the object to be a `const`
842 : /// view is actually const. Of course, there are ways of subverting this
843 : /// intended functionality and editing the data pointed into by `view` after
844 : /// this function is called; doing so is highly discouraged and results in
845 : /// undefined behavior.
846 : template <typename VectorType,
847 : Requires<is_derived_of_vector_impl_v<VectorType>> = nullptr>
848 1 : void make_const_view(const gsl::not_null<const VectorType*> view,
849 : const VectorType& vector, const size_t offset,
850 : const size_t extent) {
851 : const_cast<VectorType*>(view.get()) // NOLINT
852 : ->set_data_ref(
853 : const_cast<typename VectorType::value_type*>(vector.data()) // NOLINT
854 : + offset, // NOLINT
855 : extent);
856 : }
857 :
858 : template <typename T, typename VectorType, size_t StaticSize>
859 0 : inline bool contains_allocations(
860 : const VectorImpl<T, VectorType, StaticSize>& value) {
861 : return value.size() > StaticSize and value.is_owning();
862 : }
|