Line data Source code
1 0 : // Distributed under the MIT License. 2 : // See LICENSE.txt for details. 3 : 4 : #pragma once 5 : 6 : #include <complex> 7 : #include <cstddef> 8 : 9 : #include "DataStructures/ComplexDataVector.hpp" 10 : #include "DataStructures/DataVector.hpp" 11 : #include "Utilities/Gsl.hpp" 12 : 13 : namespace Spectral { 14 : /// \ingroup SpectralGroup 15 : /// Namespace for spin-weighted spherical harmonic utilities. 16 1 : namespace Swsh { 17 : 18 : /// \brief A set of labels for the possible representations of complex numbers 19 : /// for the `ComplexDataView` and the computations performed in the 20 : /// spin-weighted spherical harmonic transform library. 21 : /// 22 : /// \details The representation describes one of two behaviors: 23 : /// - `Interleaved`: The vectors of complex numbers will be represented by 24 : /// alternating doubles in memory. This causes both the real and imaginary part 25 : /// at a given gridpoint to be near one another, but successive real values 26 : /// farther. This is the native representation of complex data in the C++ 27 : /// standard, and is the representation needed for Blaze math 28 : /// operations. Therefore, using this representation type in libsharp 29 : /// computations will cause operations which access only the real or imaginary 30 : /// parts individually to trace over larger memory regions. However, this 31 : /// representation will give rise to fewer copying operations to perform the 32 : /// libsharp operations. 33 : /// 34 : /// - `RealsThenImags`: The vectors of complex numbers will primarily be 35 : /// represented by a pair of vectors of doubles, one for the real values and 36 : /// one for the imaginary values (the full computation cannot be performed 37 : /// exclusively in this representation, as it must return to a vector of 38 : /// `std::complex<double>` for Blaze math operations). This causes the 39 : /// successive real values for different gridpoints to be closer in memory, but 40 : /// the real and imaginary parts for a given gridpoint to be farther in 41 : /// memory. This is not the native representation for complex data in C++, so 42 : /// the data must be transformed between operations which use Blaze and the 43 : /// transform operations which use `RealsThenImags`. Therefore, using this 44 : /// representation in libsharp computations will cause operations which act on 45 : /// real or imaginary parts individually to have better memory locality (so 46 : /// likely improved cache performance, but such statements are highly 47 : /// hardware-dependent). However, this representation will give rise to more 48 : /// copying operations to perform the libsharp operations. 49 : /// 50 : /// \note The pair of representations is provided as a means to 'turn a dial' in 51 : /// optimizations. It is unclear which of these representations will be 52 : /// preferable, and it may well be the case that different representations are 53 : /// better for different calculation types or different hardware. Therefore, 54 : /// when optimizing code which uses libsharp, it is encouraged to profile the 55 : /// cost of each representation for a computation and choose the one which 56 : /// performs best. 57 0 : enum class ComplexRepresentation { Interleaved, RealsThenImags }; 58 : 59 : namespace detail { 60 : // A storage container for storing sub-vector references ("views") of a 61 : // ComplexDataVector. 62 : // 63 : // This class takes as a template argument a ComplexRepresentation to use when 64 : // returning pointers to the components of the complex data. The representation 65 : // is either: 66 : // - `ComplexRepresentation::Interleaved`, which indicates that the complex data 67 : // is represented as a vector of `std::complex<double>`, and then 68 : // ComplexDataView acts as a view (a pure reference to a subvector). 69 : // - `ComplexRepresentation::RealsThenImags`, which indicates that the complex 70 : // data is represented as a pair of vectors of `double`. This is no longer a 71 : // strict view in the typical sense, as the memory layout is different from 72 : // the data which it references. In order to minimize copies, the user must 73 : // specify when edits to the ComplexDataView are complete by calling the 74 : // `copy_back_to_source()` member function. 75 : // 76 : // Warning: For optimization reasons, mutations applied via member functions or 77 : // manipulations of data pointed to by pointers returned by member functions 78 : // *may or may not* be immediately applied to the source vector used to 79 : // construct the ComplexDataView. Correct use of this class will first perform 80 : // (potentially several) manipulations of the data via the ComplexDataView, then 81 : // as a final step flush the data to the source vector via the member function 82 : // `copy_back_to_source()`. At that point, the vector is guaranteed to be 83 : // placed in the same state as the ComplexDataView. **The process of 84 : // constructing a ComplexDataView, then mutating data in both the 85 : // ComplexDataView and the source vector, then calling `copy_back_to_source()` 86 : // is considered undefined behavior, and depends on the details of the memory 87 : // layout chosen.** 88 : template <ComplexRepresentation Representation> 89 : class ComplexDataView { 90 : public: 91 : // The Representation determines the internal data representation 92 : static const ComplexRepresentation complex_representation = Representation; 93 : 94 : // The internal storage type 95 : using value_type = std::complex<double>; 96 : 97 : // Construct a view which starts at index `offset` of the supplied 98 : // vector, and extends for `size` elements 99 : ComplexDataView(gsl::not_null<ComplexDataVector*> vector, size_t size, 100 : size_t offset = 0); 101 : 102 : // Construct a view which starts at pointer `start` and extends for 103 : // `size` elements. Need not be a part of a ComplexDataVector. 104 : ComplexDataView(std::complex<double>* start, size_t size); 105 : 106 : // For the lifetime of the data view, it points to the same portion of a 107 : // single vector. We disallow default move assignment, as it would exhibit 108 : // behavior that would contradict that. All assignment operations permitted by 109 : // the ComplexDataView are copying operations which act only on the data, not 110 : // the reference. 111 : ComplexDataView() = delete; 112 : ComplexDataView(const ComplexDataView&) = default; 113 : ComplexDataView(ComplexDataView&&) = default; 114 : ComplexDataView operator=(ComplexDataView&&) = delete; 115 : ~ComplexDataView() = default; 116 : 117 : // assign into the view the values of a same-sized ComplexDataVector 118 : ComplexDataView<Representation>& operator=(const ComplexDataVector& vector); 119 : 120 : // Assign into the view the values from a same-sized view 121 : ComplexDataView<Representation>& operator=( 122 : const ComplexDataView<Representation>& view); 123 : 124 : // Conjugate the data. 125 : void conjugate(); 126 : 127 : // Assign into the real components of a view the values from a 128 : // provided `DataVector` 129 : ComplexDataView<Representation>& assign_real(const DataVector& vector); 130 : 131 : // Assign into the imaginary components of a view the values from a 132 : // provided `DataVector` 133 : ComplexDataView<Representation>& assign_imag(const DataVector& vector); 134 : 135 : // Gets the size of the view (the number of complex entries). 136 : size_t size() const { return size_; } 137 : 138 : // Gets the stride between successive real or successive imaginary components 139 : // in memory. 140 : static constexpr size_t stride() { return stride_; } 141 : 142 : // Gets the raw pointer to the start of the real data, which are separated 143 : // from one another by `stride()`. 144 : double* real_data(); 145 : const double* real_data() const; 146 : 147 : // Gets the raw pointer to the start of the imaginary data, which are 148 : // separated from one another by `stride()`. 149 : double* imag_data(); 150 : const double* imag_data() const; 151 : 152 : // A function for manually flushing the data back to the source. This class 153 : // minimizes copies wherever possible, so the manual control of bringing the 154 : // data back to the source vector only when necessary is important for 155 : // optimization. 156 : // Warning: Until this function is called, mutations applied via other member 157 : // functions or mutations applied via the pointers obtained from other member 158 : // functions are not guaranteed to be applied to the source vector. 159 : void copy_back_to_source(); 160 : 161 : private: 162 : size_t size_; 163 : 164 : static const size_t stride_ = 165 : Representation == ComplexRepresentation::RealsThenImags ? 1 : 2; 166 : 167 : // In either case, we want to avoid unnecessary copies, so we track whether 168 : // the data has been copied from one representation to the other. 169 : bool real_slices_up_to_date_; 170 : // These two DataVectors are unused in the case of `Interleaved` 171 : // representation 172 : DataVector real_slice_; 173 : DataVector imag_slice_; 174 : 175 : double* data_real_; 176 : double* data_imag_; 177 : 178 : ComplexDataVector complex_view_; 179 : }; 180 : 181 : } // namespace detail 182 : } // namespace Swsh 183 : } // namespace Spectral