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              : #include <ostream>
       8              : 
       9              : #include "DataStructures/Index.hpp"
      10              : #include "DataStructures/SliceIterator.hpp"
      11              : #include "DataStructures/Variables.hpp"
      12              : #include "Utilities/ErrorHandling/Assert.hpp"
      13              : #include "Utilities/Gsl.hpp"
      14              : #include "Utilities/TMPL.hpp"
      15              : 
      16              : /// @{
      17              : /*!
      18              :  * \ingroup DataStructuresGroup
      19              :  * \brief Slices the data within `vars` to a codimension 1 slice. The
      20              :  * slice has a constant logical coordinate in direction `sliced_dim`,
      21              :  * slicing the volume at `fixed_index` in that dimension.  For
      22              :  * example, to get the lower boundary of `sliced_dim`, pass `0` for
      23              :  * `fixed_index`; to get the upper boundary, pass
      24              :  * `extents[sliced_dim] - 1`.
      25              :  *
      26              :  * \see add_slice_to_data
      27              :  *
      28              :  * Returns Variables class sliced to a hypersurface.
      29              :  */
      30              : template <std::size_t VolumeDim, typename TagsList>
      31            1 : void data_on_slice(const gsl::not_null<Variables<TagsList>*> interface_vars,
      32              :                    const Variables<TagsList>& vars,
      33              :                    const Index<VolumeDim>& element_extents,
      34              :                    const size_t sliced_dim, const size_t fixed_index) {
      35              :   const size_t interface_grid_points =
      36              :       element_extents.slice_away(sliced_dim).product();
      37              :   const size_t volume_grid_points = vars.number_of_grid_points();
      38              :   constexpr const size_t number_of_independent_components =
      39              :       Variables<TagsList>::number_of_independent_components;
      40              : 
      41              :   if (interface_vars->number_of_grid_points() != interface_grid_points) {
      42              :     *interface_vars = Variables<TagsList>(interface_grid_points);
      43              :   }
      44              :   using value_type = typename Variables<TagsList>::value_type;
      45              :   const value_type* vars_data = vars.data();
      46              :   value_type* interface_vars_data = interface_vars->data();
      47              :   for (SliceIterator si(element_extents, sliced_dim, fixed_index); si; ++si) {
      48              :     for (size_t i = 0; i < number_of_independent_components; ++i) {
      49              :       // clang-tidy: do not use pointer arithmetic
      50              :       interface_vars_data[si.slice_offset() +                      // NOLINT
      51              :                           i * interface_grid_points] =             // NOLINT
      52              :           vars_data[si.volume_offset() + i * volume_grid_points];  // NOLINT
      53              :     }
      54              :   }
      55              : }
      56              : 
      57              : template <std::size_t VolumeDim, typename TagsList>
      58            1 : Variables<TagsList> data_on_slice(const Variables<TagsList>& vars,
      59              :                                   const Index<VolumeDim>& element_extents,
      60              :                                   const size_t sliced_dim,
      61              :                                   const size_t fixed_index) {
      62              :   Variables<TagsList> interface_vars(
      63              :       element_extents.slice_away(sliced_dim).product());
      64              :   data_on_slice(make_not_null(&interface_vars), vars, element_extents,
      65              :                 sliced_dim, fixed_index);
      66              :   return interface_vars;
      67              : }
      68              : /// @}
      69              : 
      70              : /*!
      71              :  * \ingroup DataStructuresGroup
      72              :  * \brief Adds data on a codimension 1 slice to a volume quantity. The
      73              :  * slice has a constant logical coordinate in direction `sliced_dim`,
      74              :  * slicing the volume at `fixed_index` in that dimension.  For
      75              :  * example, to add to the lower boundary of `sliced_dim`, pass `0` for
      76              :  * `fixed_index`; to add to the upper boundary, pass
      77              :  * `extents[sliced_dim] - 1`.
      78              :  *
      79              :  * \see data_on_slice
      80              :  */
      81              : template <std::size_t VolumeDim, typename... VolumeTags, typename... SliceTags>
      82            1 : void add_slice_to_data(
      83              :     const gsl::not_null<Variables<tmpl::list<VolumeTags...>>*> volume_vars,
      84              :     const Variables<tmpl::list<SliceTags...>>& vars_on_slice,
      85              :     const Index<VolumeDim>& extents, const size_t sliced_dim,
      86              :     const size_t fixed_index) {
      87              :   static_assert((std::is_same_v<db::remove_all_prefixes<VolumeTags>,
      88              :                                 db::remove_all_prefixes<SliceTags>> and
      89              :                  ...));
      90              :   static_assert(
      91              :       (std::is_same_v<typename VolumeTags::type, typename SliceTags::type> and
      92              :        ...),
      93              :       "Tensor types do not match.");
      94              :   constexpr const size_t number_of_independent_components =
      95              :       Variables<tmpl::list<VolumeTags...>>::number_of_independent_components;
      96              :   const size_t volume_grid_points = extents.product();
      97              :   const size_t slice_grid_points = extents.slice_away(sliced_dim).product();
      98              :   ASSERT(volume_vars->number_of_grid_points() == volume_grid_points,
      99              :          "volume_vars has wrong number of grid points.  Expected "
     100              :              << volume_grid_points << ", got "
     101              :              << volume_vars->number_of_grid_points());
     102              :   ASSERT(vars_on_slice.number_of_grid_points() == slice_grid_points,
     103              :          "vars_on_slice has wrong number of grid points.  Expected "
     104              :              << slice_grid_points << ", got "
     105              :              << vars_on_slice.number_of_grid_points());
     106              :   using value_type = typename Variables<tmpl::list<VolumeTags...>>::value_type;
     107              :   value_type* const volume_data = volume_vars->data();
     108              :   const value_type* const slice_data = vars_on_slice.data();
     109              :   for (SliceIterator si(extents, sliced_dim, fixed_index); si; ++si) {
     110              :     for (size_t i = 0; i < number_of_independent_components; ++i) {
     111              :       // clang-tidy: do not use pointer arithmetic
     112              :       volume_data[si.volume_offset() + i * volume_grid_points] +=  // NOLINT
     113              :           slice_data[si.slice_offset() + i * slice_grid_points];   // NOLINT
     114              :     }
     115              :   }
     116              : }