Line data Source code
1 1 : // Distributed under the MIT License.
2 : // See LICENSE.txt for details.
3 :
4 : /// \file
5 : /// Defines classes for Tensor
6 :
7 : #pragma once
8 :
9 : #include <concepts>
10 : #include <cstddef>
11 : #include <pup.h>
12 : #include <pup_stl.h>
13 :
14 : #include "DataStructures/ComplexDataVector.hpp"
15 : #include "DataStructures/ComplexModalVector.hpp"
16 : #include "DataStructures/DataVector.hpp"
17 : #include "DataStructures/ModalVector.hpp"
18 : #include "DataStructures/SpinWeighted.hpp"
19 : #include "DataStructures/Tensor/Expressions/AddSubtract.hpp"
20 : #include "DataStructures/Tensor/Expressions/Contract.hpp"
21 : #include "DataStructures/Tensor/Expressions/DataTypeSupport.hpp"
22 : #include "DataStructures/Tensor/Expressions/Divide.hpp"
23 : #include "DataStructures/Tensor/Expressions/Evaluate.hpp"
24 : #include "DataStructures/Tensor/Expressions/IndexPropertyCheck.hpp"
25 : #include "DataStructures/Tensor/Expressions/LhsTensorSymmAndIndices.hpp"
26 : #include "DataStructures/Tensor/Expressions/Negate.hpp"
27 : #include "DataStructures/Tensor/Expressions/NumberAsExpression.hpp"
28 : #include "DataStructures/Tensor/Expressions/Product.hpp"
29 : #include "DataStructures/Tensor/Expressions/SpatialSpacetimeIndex.hpp"
30 : #include "DataStructures/Tensor/Expressions/SquareRoot.hpp"
31 : #include "DataStructures/Tensor/Expressions/TensorAsExpression.hpp"
32 : #include "DataStructures/Tensor/Expressions/TensorExpression.hpp"
33 : #include "DataStructures/Tensor/Expressions/TensorIndex.hpp"
34 : #include "DataStructures/Tensor/Expressions/TensorIndexTransformation.hpp"
35 : #include "DataStructures/Tensor/Expressions/TimeIndex.hpp"
36 : #include "DataStructures/Tensor/IndexType.hpp"
37 : #include "DataStructures/Tensor/Structure.hpp"
38 : #include "DataStructures/Tensor/TypeAliases.hpp"
39 : #include "Utilities/ErrorHandling/Error.hpp"
40 : #include "Utilities/ForceInline.hpp"
41 : #include "Utilities/Gsl.hpp"
42 : #include "Utilities/MakeArray.hpp"
43 : #include "Utilities/MakeWithValue.hpp"
44 : #include "Utilities/PrettyType.hpp"
45 : #include "Utilities/SetNumberOfGridPoints.hpp"
46 : #include "Utilities/Simd/Simd.hpp"
47 : #include "Utilities/StdHelpers.hpp"
48 : #include "Utilities/TMPL.hpp"
49 : #include "Utilities/TypeTraits.hpp"
50 : #include "Utilities/TypeTraits/IsStreamable.hpp"
51 : // Not all the TensorExpression includes are necessary, but we include them so
52 : // that Tensor.hpp provides a uniform interface to Tensors and TensorExpressions
53 :
54 : /// \cond
55 : template <typename X, typename Symm = Symmetry<>,
56 : typename IndexList = index_list<>>
57 : class Tensor;
58 : /// \endcond
59 :
60 : namespace Tensor_detail {
61 : template <typename T, typename = void>
62 : struct is_tensor : std::false_type {};
63 : template <typename X, typename Symm, typename IndexList>
64 : struct is_tensor<Tensor<X, Symm, IndexList>> : std::true_type {};
65 : } // namespace Tensor_detail
66 :
67 : /*!
68 : * \ingroup TensorGroup
69 : * \brief Represents an object with multiple components
70 : *
71 : * \details
72 : * Tensor is a container that represents indexable geometric objects. Each index
73 : * has its own dimension, valence, and frame and must be either spatial or
74 : * spacetime. Note that the dimension passed to `SpatialIndex` and
75 : * `SpacetimeIndex` is always the spatial dimension of the index. Tensors with
76 : * symmetric indices are stored only once and must be of the same
77 : * type. A list of available type aliases can be found in the ::tnsr namespace
78 : * where the adopted conventions are:
79 : *
80 : * 1. Upper case for contravariant or upper indices, lower case for covariant or
81 : * lower indices.
82 : *
83 : * 2. `a, b, c, d` are used for spacetime indices while `i, j, k, l` are used
84 : * for spatial indices.
85 : *
86 : * 3. `::Scalar` is not inside the `::tnsr` namespace but is used to represent
87 : * a scalar with no indices.
88 : *
89 : * \example
90 : * \snippet Test_Tensor.cpp scalar
91 : * \snippet Test_Tensor.cpp spatial_vector
92 : * \snippet Test_Tensor.cpp spacetime_vector
93 : * \snippet Test_Tensor.cpp rank_3_122
94 : *
95 : * \tparam X the type held
96 : * \tparam Symm the ::Symmetry of the indices
97 : * \tparam IndexList a typelist of \ref SpacetimeIndex "TensorIndexType"'s
98 : */
99 : template <typename X, typename Symm, template <typename...> class IndexList,
100 : typename... Indices>
101 1 : class Tensor<X, Symm, IndexList<Indices...>> {
102 : static_assert(sizeof...(Indices) < 5,
103 : "If you are sure you need rank 5 or higher Tensor's please "
104 : "file an issue on GitHub or discuss with a core developer of "
105 : "SpECTRE.");
106 : static_assert(
107 : std::is_same_v<X, std::complex<double>> or std::is_same_v<X, double> or
108 : std::is_same_v<X, ComplexDataVector> or
109 : std::is_same_v<X, ComplexModalVector> or
110 : std::is_same_v<X, DataVector> or std::is_same_v<X, ModalVector> or
111 : is_spin_weighted_of_v<ComplexDataVector, X> or
112 : is_spin_weighted_of_v<ComplexModalVector, X> or
113 : simd::is_batch<X>::value,
114 : "Unsupported type. While other types are technically possible it is not "
115 : "clear that Tensor is the correct container for them. Please seek advice "
116 : "on the topic by discussing with the SpECTRE developers.");
117 :
118 : public:
119 : /// The type of the sequence that holds the data
120 1 : using storage_type =
121 : std::array<X, Tensor_detail::Structure<Symm, Indices...>::size()>;
122 : /// The type that is stored by the Tensor
123 1 : using type = X;
124 : /// Typelist of the symmetry of the Tensor
125 : ///
126 : /// \details
127 : /// For a rank-3 tensor symmetric in the last two indices,
128 : /// \f$T_{a(bc)}\f$, the ::Symmetry is `<2, 1, 1>`. For a non-symmetric rank-2
129 : /// tensor the ::Symmetry is `<2, 1>`.
130 1 : using symmetry = Symm;
131 : /// Typelist of the \ref SpacetimeIndex "TensorIndexType"'s that the
132 : /// Tensor has
133 1 : using index_list = tmpl::list<Indices...>;
134 : /// The number of indices that the Tensor has
135 1 : static constexpr size_t num_tensor_indices = sizeof...(Indices);
136 : /// The Tensor_detail::Structure for the particular tensor index structure
137 : ///
138 : /// Each tensor index structure, e.g. \f$T_{ab}\f$, \f$T_a{}^b\f$ or
139 : /// \f$T^{ab}\f$ has its own Tensor_detail::TensorStructure that holds
140 : /// information about how the data is stored, what the multiplicity of the
141 : /// stored indices are, the number of (independent) components, etc.
142 1 : using structure = Tensor_detail::Structure<Symm, Indices...>;
143 : /// The type of the TensorExpression that would represent this Tensor in a
144 : /// tensor expression.
145 : template <typename ArgsList>
146 1 : using TE = tenex::TensorAsExpression<Tensor<X, Symm, IndexList<Indices...>>,
147 : ArgsList>;
148 :
149 0 : Tensor() = default;
150 0 : ~Tensor() = default;
151 0 : Tensor(const Tensor&) = default;
152 0 : Tensor(Tensor&&) = default;
153 0 : Tensor& operator=(const Tensor&) = default;
154 0 : Tensor& operator=(Tensor&&) = default;
155 :
156 : /// Initialize a vector or scalar from an array
157 : ///
158 : /// \example
159 : /// \snippet Test_Tensor.cpp init_vector
160 : /// \param data the values of the individual components of the Vector
161 1 : explicit Tensor(storage_type data)
162 : requires(sizeof...(Indices) <= 1);
163 :
164 : /// Constructor that passes "args" to constructor of X and initializes each
165 : /// component to be the same
166 : template <typename Arg0, typename... Args>
167 1 : explicit Tensor(Arg0&& arg0, Args&&... args)
168 : // NOLINTNEXTLINE(readability-simplify-boolean-expr)
169 : requires(not(sizeof...(Args) == 0 and
170 : std::same_as<std::decay_t<Arg0>, Tensor>) and
171 : std::constructible_from<X, Arg0, Args...>);
172 :
173 0 : using value_type = typename storage_type::value_type;
174 0 : using reference = typename storage_type::reference;
175 0 : using const_reference = typename storage_type::const_reference;
176 0 : using iterator = typename storage_type::iterator;
177 0 : using const_iterator = typename storage_type::const_iterator;
178 0 : using pointer = typename storage_type::pointer;
179 0 : using const_pointer = typename storage_type::const_pointer;
180 0 : using reverse_iterator = typename storage_type::reverse_iterator;
181 0 : using const_reverse_iterator = typename storage_type::const_reverse_iterator;
182 :
183 0 : iterator begin() { return data_.begin(); }
184 0 : const_iterator begin() const { return data_.begin(); }
185 0 : const_iterator cbegin() const { return data_.begin(); }
186 :
187 0 : iterator end() { return data_.end(); }
188 0 : const_iterator end() const { return data_.end(); }
189 0 : const_iterator cend() const { return data_.end(); }
190 :
191 0 : reverse_iterator rbegin() { return data_.rbegin(); }
192 0 : const_reverse_iterator rbegin() const { return data_.rbegin(); }
193 0 : const_reverse_iterator crbegin() const { return data_.rbegin(); }
194 :
195 0 : reverse_iterator rend() { return data_.rend(); }
196 0 : const_reverse_iterator rend() const { return data_.rend(); }
197 0 : const_reverse_iterator crend() const { return data_.rend(); }
198 :
199 : /// @{
200 : /// Get data entry using an array representing a tensor index
201 : ///
202 : /// \details
203 : /// Let \f$T_{abc}\f$ be a Tensor.
204 : /// Then `get({{0, 2, 1}})` returns the \f$T_{0 2 1}\f$ component.
205 : /// \param tensor_index the index at which to get the data
206 : template <typename T>
207 1 : SPECTRE_ALWAYS_INLINE constexpr reference get(
208 : const std::array<T, sizeof...(Indices)>& tensor_index) {
209 : return gsl::at(data_, structure::get_storage_index(tensor_index));
210 : }
211 : template <typename T>
212 1 : SPECTRE_ALWAYS_INLINE constexpr const_reference get(
213 : const std::array<T, sizeof...(Indices)>& tensor_index) const {
214 : return gsl::at(data_, structure::get_storage_index(tensor_index));
215 : }
216 : /// @}
217 : /// @{
218 : /// Get data entry using a list of integers representing a tensor index
219 : ///
220 : /// \details
221 : /// Let \f$T_{abc}\f$ be a Tensor.
222 : /// Then `get(0, 2, 1)` returns the \f$T_{0 2 1}\f$ component.
223 : /// \param n the index at which to get the data
224 : template <typename... N>
225 1 : constexpr reference get(N... n) {
226 : static_assert(
227 : sizeof...(Indices) == sizeof...(N),
228 : "the number of tensor indices specified must match the rank of "
229 : "the tensor");
230 : return gsl::at(data_, structure::get_storage_index(n...));
231 : }
232 : template <typename... N>
233 1 : constexpr const_reference get(N... n) const {
234 : static_assert(
235 : sizeof...(Indices) == sizeof...(N),
236 : "the number of tensor indices specified must match the rank of "
237 : "the tensor");
238 : return gsl::at(data_, structure::get_storage_index(n...));
239 : }
240 : /// @}
241 :
242 : /// @{
243 : /// Retrieve the index `N...` by computing the storage index at compile time
244 : // clang-tidy: redundant declaration (bug in clang-tidy)
245 : template <int... N, typename... Args>
246 : friend SPECTRE_ALWAYS_INLINE constexpr typename Tensor<Args...>::reference
247 1 : get(Tensor<Args...>& t); // NOLINT
248 : // clang-tidy: redundant declaration (bug in clang-tidy)
249 : template <int... N, typename... Args>
250 : friend SPECTRE_ALWAYS_INLINE constexpr
251 : typename Tensor<Args...>::const_reference
252 1 : get(const Tensor<Args...>& t); // NOLINT
253 : /// @}
254 :
255 : /// @{
256 : /// Retrieve a TensorExpression object with the index structure passed in
257 : template <typename... TensorIndices>
258 1 : SPECTRE_ALWAYS_INLINE constexpr auto operator()(
259 : TensorIndices... /*meta*/) const {
260 : static_assert((... and tt::is_tensor_index<TensorIndices>::value),
261 : "The tensor expression must be created using TensorIndex "
262 : "objects to represent generic indices, e.g. ti::a, ti::b, "
263 : "etc.");
264 : static_assert(
265 : tenex::tensorindex_list_is_valid<tmpl::list<TensorIndices...>>::value,
266 : "Cannot create a tensor expression with a repeated generic index. "
267 : "(Note that the concrete time indices, ti::T and ti::t, can be "
268 : "repeated.) If you intend to contract, ensure that the indices to "
269 : "contract have opposite valences.");
270 : static_assert(
271 : std::is_same_v<tmpl::integral_list<UpLo, TensorIndices::valence...>,
272 : tmpl::integral_list<UpLo, Indices::ul...>>,
273 : "The valences of the generic indices in the expression do "
274 : "not match the valences of the indices in the Tensor.");
275 : static_assert((... and (not(TensorIndices::is_spacetime and
276 : (Indices::index_type == IndexType::Spatial)))),
277 : "Cannot use a spacetime index for a spatial index. e.g. "
278 : "Cannot do R(ti::a), where R's index is spatial, because "
279 : "ti::a denotes a generic spacetime index.");
280 : return tenex::contract(TE<tmpl::list<TensorIndices...>>{*this});
281 : }
282 : /// @}
283 :
284 : /// @{
285 : /// Return i'th component of storage vector
286 1 : constexpr reference operator[](const size_t storage_index) {
287 : return gsl::at(data_, storage_index);
288 : }
289 1 : constexpr const_reference operator[](const size_t storage_index) const {
290 : return gsl::at(data_, storage_index);
291 : }
292 : /// @}
293 :
294 : /// Return the number of independent components of the Tensor
295 : ///
296 : /// \details
297 : /// Returns the number of independent components of the Tensor taking into
298 : /// account symmetries. For example, let \f$T_{ab}\f$ be a n-dimensional
299 : /// rank-2 symmetric tensor, then the number of independent components is
300 : /// \f$n(n+1)/2\f$.
301 1 : SPECTRE_ALWAYS_INLINE static constexpr size_t size() {
302 : return structure::size();
303 : }
304 :
305 : /// Returns the rank of the Tensor
306 : ///
307 : /// \details
308 : /// The rank of a tensor is the number of indices it has. For example, the
309 : /// tensor \f$v^a\f$ is rank-1, the tensor \f$\phi\f$ is rank-0, and the
310 : /// tensor \f$T_{abc}\f$ is rank-3.
311 1 : SPECTRE_ALWAYS_INLINE static constexpr size_t rank() {
312 : return sizeof...(Indices);
313 : }
314 :
315 : /// @{
316 : /// Given an iterator or storage index, get the canonical tensor index.
317 : /// For scalars this is defined to be std::array<int, 1>{{0}}
318 : SPECTRE_ALWAYS_INLINE constexpr std::array<size_t, sizeof...(Indices)>
319 1 : get_tensor_index(const const_iterator& iter) const {
320 : return structure::get_canonical_tensor_index(
321 : static_cast<size_t>(iter - begin()));
322 : }
323 : SPECTRE_ALWAYS_INLINE static constexpr std::array<size_t, sizeof...(Indices)>
324 1 : get_tensor_index(const size_t storage_index) {
325 : return structure::get_canonical_tensor_index(storage_index);
326 : }
327 : /// @}
328 :
329 : /// @{
330 : /// Get the storage index of the tensor index. Should only be used when
331 : /// optimizing code in which computing the storage index is a bottleneck.
332 : template <typename... N>
333 1 : SPECTRE_ALWAYS_INLINE static constexpr size_t get_storage_index(
334 : const N... args) {
335 : return structure::get_storage_index(args...);
336 : }
337 : template <typename I>
338 1 : SPECTRE_ALWAYS_INLINE static constexpr size_t get_storage_index(
339 : const std::array<I, sizeof...(Indices)>& tensor_index) {
340 : return structure::get_storage_index(tensor_index);
341 : }
342 : /// @}
343 :
344 : /// @{
345 : /// Given an iterator or storage index, get the multiplicity of an index
346 : ///
347 : /// \see TensorMetafunctions::compute_multiplicity
348 1 : SPECTRE_ALWAYS_INLINE constexpr size_t multiplicity(
349 : const const_iterator& iter) const {
350 : return structure::multiplicity(static_cast<size_t>(iter - begin()));
351 : }
352 1 : SPECTRE_ALWAYS_INLINE static constexpr size_t multiplicity(
353 : const size_t storage_index) {
354 : return structure::multiplicity(storage_index);
355 : }
356 : /// @}
357 :
358 : /// @{
359 : /// Get dimensionality of i'th tensor index
360 : ///
361 : /// \snippet Test_Tensor.cpp index_dim
362 : /// \see ::index_dim
363 1 : SPECTRE_ALWAYS_INLINE static constexpr size_t index_dim(const size_t i) {
364 : static_assert(sizeof...(Indices),
365 : "A scalar does not have any indices from which you can "
366 : "retrieve the dimensionality.");
367 : return gsl::at(structure::dims(), i);
368 : }
369 : /// @}
370 :
371 : /// @{
372 : /// Return an array corresponding to the ::Symmetry of the Tensor
373 : SPECTRE_ALWAYS_INLINE static constexpr std::array<int, sizeof...(Indices)>
374 1 : symmetries() {
375 : return structure::symmetries();
376 : }
377 : /// @}
378 :
379 : /// @{
380 : /// Return array of the ::IndexType's (spatial or spacetime)
381 : SPECTRE_ALWAYS_INLINE static constexpr std::array<IndexType,
382 : sizeof...(Indices)>
383 1 : index_types() {
384 : return structure::index_types();
385 : }
386 : /// @}
387 :
388 : /// @{
389 : /// Return array of dimensionality of each index
390 : ///
391 : /// \snippet Test_Tensor.cpp index_dim
392 : /// \see index_dim ::index_dim
393 : SPECTRE_ALWAYS_INLINE static constexpr std::array<size_t, sizeof...(Indices)>
394 1 : index_dims() {
395 : return structure::dims();
396 : }
397 : /// @}
398 :
399 : /// @{
400 : /// Return array of the valence of each index (::UpLo)
401 : SPECTRE_ALWAYS_INLINE static constexpr std::array<UpLo, sizeof...(Indices)>
402 1 : index_valences() {
403 : return structure::index_valences();
404 : }
405 : /// @}
406 :
407 : /// @{
408 : /// Returns std::tuple of the ::Frame of each index
409 1 : SPECTRE_ALWAYS_INLINE static constexpr auto index_frames() {
410 : return Tensor_detail::Structure<Symm, Indices...>::index_frames();
411 : }
412 : /// @}
413 :
414 : /// @{
415 : /// \brief Given a tensor index, get the canonical label associated with the
416 : /// canonical \ref SpacetimeIndex "TensorIndexType"
417 : ///
418 : /// \param tensor_index The index of the tensor component to label
419 : /// \param axis_labels The labels for the indices. Defaults to "t", "x", "y"
420 : /// and "z" for spacetime indices and "x", "y" and "z" for spatial indices.
421 : /// Note that a tensor can have indices of different types, so we specify
422 : /// labels for each index individually.
423 : template <typename T = int>
424 1 : static std::string component_name(
425 : const std::array<T, rank()>& tensor_index = std::array<T, rank()>{},
426 : const std::array<std::string, rank()>& axis_labels =
427 : make_array<rank()>(std::string(""))) {
428 : return structure::component_name(tensor_index, axis_labels);
429 : }
430 : /// @}
431 :
432 : /// @{
433 : /// \brief Suffix to append to the tensor name that indicates the component
434 : ///
435 : /// The suffix is empty for scalars, otherwise it is an underscore followed by
436 : /// the `Tensor::component_name` of either the `tensor_index` or the canonical
437 : /// tensor index obtained from the `storage_index`. Use `axis_labels` to
438 : /// overwrite the default labels for each component (see
439 : /// `Tensor::component_name`).
440 : ///
441 : /// An example use case for the suffix is to label tensor components in
442 : /// data files.
443 : ///
444 : /// \see Tensor::component_name
445 : template <typename IndexType = int>
446 1 : static std::string component_suffix(
447 : const std::array<IndexType, rank()>& tensor_index =
448 : std::array<IndexType, rank()>{},
449 : const std::array<std::string, rank()>& axis_labels =
450 : make_array<rank()>(std::string(""))) {
451 : return rank() == 0 ? "" : "_" + component_name(tensor_index, axis_labels);
452 : }
453 :
454 1 : static std::string component_suffix(
455 : const size_t storage_index,
456 : const std::array<std::string, rank()>& axis_labels) {
457 : return component_suffix(get_tensor_index(storage_index), axis_labels);
458 : }
459 :
460 1 : static std::string component_suffix(size_t storage_index);
461 : /// @}
462 :
463 : /// Copy tensor data into an `std::vector<X>` along with the
464 : /// component names into a `std::vector<std::string>`
465 : /// \requires `std::is_same<X, DataVector>::%value` is true
466 1 : std::pair<std::vector<std::string>, std::vector<X>> get_vector_of_data()
467 : const;
468 :
469 : /// \cond HIDDEN_SYMBOLS
470 : /// Serialization function used by Charm++
471 : // NOLINTNEXTLINE(google-runtime-references)
472 : void pup(PUP::er& p) { p | data_; }
473 : /// \endcond
474 :
475 : private:
476 : // clang-tidy: redundant declaration
477 : /// \cond
478 : template <int I, class... Ts>
479 : friend SPECTRE_ALWAYS_INLINE constexpr size_t index_dim( // NOLINT
480 : const Tensor<Ts...>& /*t*/);
481 : /// \endcond
482 :
483 0 : storage_type data_;
484 : };
485 :
486 : // ================================================================
487 : // Template Definitions - Variadic templates must be in header
488 : // ================================================================
489 :
490 : template <typename X, typename Symm, template <typename...> class IndexList,
491 : typename... Indices>
492 : // Implementation note: we explicitly prevent inlining and mark the function
493 : // as used so that GDB's pretty printing facilities have access to this
494 : // function.
495 0 : __attribute__((noinline)) __attribute__((used)) std::string
496 : Tensor<X, Symm, IndexList<Indices...>>::component_suffix(
497 : const size_t storage_index) {
498 : return component_suffix(get_tensor_index(storage_index),
499 : make_array<rank()>(std::string("")));
500 : }
501 :
502 : template <typename X, typename Symm, template <typename...> class IndexList,
503 : typename... Indices>
504 : Tensor<X, Symm, IndexList<Indices...>>::Tensor(storage_type data)
505 : requires(sizeof...(Indices) <= 1)
506 : : data_(std::move(data)) {}
507 :
508 : template <typename X, typename Symm, template <typename...> class IndexList,
509 : typename... Indices>
510 : template <typename Arg0, typename... Args>
511 : Tensor<X, Symm, IndexList<Indices...>>::Tensor(Arg0&& arg0, Args&&... args)
512 : // NOLINTNEXTLINE(readability-simplify-boolean-expr)
513 : requires(not(sizeof...(Args) == 0 and
514 : std::same_as<std::decay_t<Arg0>, Tensor>) and
515 : std::constructible_from<X, Arg0, Args...>)
516 : : data_(make_array<size(), X>(std::forward<Arg0>(arg0),
517 : std::forward<Args>(args)...)) {}
518 :
519 : template <typename X, typename Symm, template <typename...> class IndexList,
520 : typename... Indices>
521 : std::pair<std::vector<std::string>, std::vector<X>>
522 : Tensor<X, Symm, IndexList<Indices...>>::get_vector_of_data() const {
523 : std::vector<value_type> serialized_tensor(size());
524 : std::vector<std::string> component_names(size());
525 : for (size_t i = 0; i < data_.size(); ++i) {
526 : component_names[i] = component_name(get_tensor_index(i));
527 : serialized_tensor[i] = gsl::at(data_, i);
528 : }
529 : return std::make_pair(component_names, serialized_tensor);
530 : }
531 :
532 : template <int... N, typename... Args>
533 0 : SPECTRE_ALWAYS_INLINE constexpr typename Tensor<Args...>::reference get(
534 : Tensor<Args...>& t) {
535 : static_assert(Tensor<Args...>::rank() == sizeof...(N),
536 : "the number of tensor indices specified must match the rank "
537 : "of the tensor");
538 : return gsl::at(
539 : t.data_, Tensor<Args...>::structure::template get_storage_index<N...>());
540 : }
541 :
542 : template <int... N, typename... Args>
543 0 : SPECTRE_ALWAYS_INLINE constexpr typename Tensor<Args...>::const_reference get(
544 : const Tensor<Args...>& t) {
545 : static_assert(Tensor<Args...>::rank() == sizeof...(N),
546 : "the number of tensor indices specified must match the rank "
547 : "of the tensor");
548 : return gsl::at(
549 : t.data_, Tensor<Args...>::structure::template get_storage_index<N...>());
550 : }
551 :
552 : template <typename X, typename Symm, template <typename...> class IndexList,
553 : typename... Indices>
554 0 : bool operator==(const Tensor<X, Symm, IndexList<Indices...>>& lhs,
555 : const Tensor<X, Symm, IndexList<Indices...>>& rhs) {
556 : return std::equal(lhs.begin(), lhs.end(), rhs.begin());
557 : }
558 : template <typename X, typename Symm, template <typename...> class IndexList,
559 : typename... Indices>
560 0 : bool operator!=(const Tensor<X, Symm, IndexList<Indices...>>& lhs,
561 : const Tensor<X, Symm, IndexList<Indices...>>& rhs) {
562 : return not(lhs == rhs);
563 : }
564 :
565 : /// \ingroup TensorGroup
566 : /// Get dimensionality of i'th tensor index
567 : ///
568 : /// \snippet Test_Tensor.cpp index_dim
569 : template <int I, class... Ts>
570 1 : SPECTRE_ALWAYS_INLINE constexpr size_t index_dim(const Tensor<Ts...>& /*t*/) {
571 : return Tensor<Ts...>::structure::template dim<I>();
572 : }
573 :
574 : // We place the stream operators in the header file so they do not need to be
575 : // explicitly instantiated.
576 : template <typename X, typename Symm, template <typename...> class IndexList,
577 : typename... Indices>
578 0 : std::ostream& operator<<(std::ostream& os,
579 : const Tensor<X, Symm, IndexList<Indices...>>& x) {
580 : static_assert(tt::is_streamable_v<decltype(os), X>,
581 : "operator<< is not defined for the type you are trying to "
582 : "stream in Tensor");
583 : for (size_t i = 0; i < x.size() - 1; ++i) {
584 : os << "T" << x.get_tensor_index(i) << "=" << x[i] << "\n";
585 : }
586 : size_t i = x.size() - 1;
587 : os << "T" << x.get_tensor_index(i) << "=" << x[i];
588 : return os;
589 : }
590 :
591 : namespace MakeWithValueImpls {
592 : template <typename T, typename... Structure>
593 0 : struct NumberOfPoints<Tensor<T, Structure...>> {
594 : static SPECTRE_ALWAYS_INLINE size_t
595 0 : apply(const Tensor<T, Structure...>& input) {
596 : return number_of_points(*input.begin());
597 : }
598 : };
599 :
600 : template <typename T, typename... Structure>
601 0 : struct MakeWithSize<Tensor<T, Structure...>> {
602 : template <typename U>
603 0 : static SPECTRE_ALWAYS_INLINE Tensor<T, Structure...> apply(const size_t size,
604 : const U value) {
605 : return Tensor<T, Structure...>(make_with_value<T>(size, value));
606 : }
607 : };
608 :
609 : template <typename... Structure, typename T>
610 0 : struct MakeWithValueImpl<Tensor<double, Structure...>, T> {
611 0 : static SPECTRE_ALWAYS_INLINE Tensor<double, Structure...> apply(
612 : const T& /*input*/, const double value) {
613 : return Tensor<double, Structure...>(value);
614 : }
615 : };
616 :
617 : template <typename... Structure, typename T>
618 : struct MakeWithValueImpl<Tensor<std::complex<double>, Structure...>, T> {
619 : static SPECTRE_ALWAYS_INLINE Tensor<std::complex<double>, Structure...> apply(
620 : const T& /*input*/, const std::complex<double> value) {
621 : return Tensor<std::complex<double>, Structure...>(value);
622 : }
623 : };
624 : } // namespace MakeWithValueImpls
625 :
626 : template <typename T, typename... Structure>
627 0 : struct SetNumberOfGridPointsImpls::SetNumberOfGridPointsImpl<
628 : Tensor<T, Structure...>> {
629 0 : static constexpr bool is_trivial = SetNumberOfGridPointsImpl<T>::is_trivial;
630 0 : static SPECTRE_ALWAYS_INLINE void apply(
631 : const gsl::not_null<Tensor<T, Structure...>*> result, const size_t size) {
632 : for (auto& component : *result) {
633 : set_number_of_grid_points(make_not_null(&component), size);
634 : }
635 : }
636 : };
|
