Line data    Source code

       1           0 : // Distributed under the MIT License.
       2             : // See LICENSE.txt for details.
       3             : 
       4             : #pragma once
       5             : 
       6             : #include <cstddef>
       7             : 
       8             : #include "DataStructures/DataVector.hpp"
       9             : #include "DataStructures/Tensor/Tensor.hpp"
      10             : #include "DataStructures/Variables.hpp"
      11             : #include "Utilities/ErrorHandling/Assert.hpp"
      12             : #include "Utilities/Gsl.hpp"
      13             : #include "Utilities/TMPL.hpp"
      14             : 
      15             : /// \ingroup DataStructuresGroup
      16             : /// Copy a given index of each component of a `Tensor<DataVector>` or
      17             : /// `Variables<DataVector>` into a `Tensor<double>`, single point
      18             : /// `Tensor<DataVector>`, or single-point `Variables<DataVector>`.
      19             : ///
      20             : /// \note There is no by-value overload extracting to a
      21             : /// `Tensor<DataVector>`.  This is both for the practical reason that
      22             : /// it would be ambiguous with the `Tensor<double>` overload and
      23             : /// because allocating multiple `DataVector`s for the return type
      24             : /// would usually be very inefficient.
      25             : ///
      26             : /// \see overwrite_point
      27             : /// @{
      28             : template <typename... Structure>
      29           1 : void extract_point(
      30             :     const gsl::not_null<Tensor<double, Structure...>*> destination,
      31             :     const Tensor<DataVector, Structure...>& source, const size_t index) {
      32             :   for (size_t i = 0; i < destination->size(); ++i) {
      33             :     (*destination)[i] = source[i][index];
      34             :   }
      35             : }
      36             : 
      37             : template <typename... Structure>
      38           1 : Tensor<double, Structure...> extract_point(
      39             :     const Tensor<DataVector, Structure...>& tensor, const size_t index) {
      40             :   Tensor<double, Structure...> result;
      41             :   extract_point(make_not_null(&result), tensor, index);
      42             :   return result;
      43             : }
      44             : 
      45             : template <typename... Structure>
      46           1 : void extract_point(
      47             :     const gsl::not_null<Tensor<DataVector, Structure...>*> destination,
      48             :     const Tensor<DataVector, Structure...>& source, const size_t index) {
      49             :   ASSERT(destination->begin()->size() == 1,
      50             :          "Output tensor components have wrong size: "
      51             :          << destination->begin()->size());
      52             :   for (size_t i = 0; i < destination->size(); ++i) {
      53             :     (*destination)[i][0] = source[i][index];
      54             :   }
      55             : }
      56             : 
      57             : template <typename... Tags>
      58           1 : void extract_point(
      59             :     const gsl::not_null<Variables<tmpl::list<Tags...>>*> result,
      60             :     const Variables<tmpl::list<Tags...>>& variables, const size_t index) {
      61             :   result->initialize(1);
      62             :   expand_pack((extract_point(
      63             :       make_not_null(&get<Tags>(*result)), get<Tags>(variables), index), 0)...);
      64             : }
      65             : 
      66             : template <typename... Tags>
      67           1 : Variables<tmpl::list<Tags...>> extract_point(
      68             :     const Variables<tmpl::list<Tags...>>& variables, const size_t index) {
      69             :   Variables<tmpl::list<Tags...>> result(1);
      70             :   extract_point(make_not_null(&result), variables, index);
      71             :   return result;
      72             : }
      73             : /// @}
      74             : 
      75             : /// \ingroup DataStructuresGroup
      76             : /// Copy a `Tensor<double>`, single point `Tensor<DataVector>`, or
      77             : /// single-point `Variables<DataVector>` into the given index of each
      78             : /// component of a `Tensor<DataVector>` or `Variables<DataVector>`.
      79             : ///
      80             : /// \see extract_point
      81             : /// @{
      82             : template <typename... Structure>
      83           1 : void overwrite_point(
      84             :     const gsl::not_null<Tensor<DataVector, Structure...>*> destination,
      85             :     const Tensor<double, Structure...>& source, const size_t index) {
      86             :   for (size_t i = 0; i < destination->size(); ++i) {
      87             :     (*destination)[i][index] = source[i];
      88             :   }
      89             : }
      90             : 
      91             : template <typename... Structure>
      92           1 : void overwrite_point(
      93             :     const gsl::not_null<Tensor<DataVector, Structure...>*> destination,
      94             :     const Tensor<DataVector, Structure...>& source, const size_t index) {
      95             :   ASSERT(source.begin()->size() == 1,
      96             :          "Cannot overwrite with " << source.begin()->size() << " points.");
      97             :   for (size_t i = 0; i < destination->size(); ++i) {
      98             :     (*destination)[i][index] = source[i][0];
      99             :   }
     100             : }
     101             : 
     102             : template <typename... Tags>
     103           1 : void overwrite_point(
     104             :     const gsl::not_null<Variables<tmpl::list<Tags...>>*> destination,
     105             :     const Variables<tmpl::list<Tags...>>& source, const size_t index) {
     106             :   ASSERT(source.number_of_grid_points() == 1,
     107             :          "Must overwrite with a single point.");
     108             :   expand_pack((overwrite_point(make_not_null(&get<Tags>(*destination)),
     109             :                                extract_point(get<Tags>(source), 0), index),
     110             :                0)...);
     111             : }
     112             : /// @}