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