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