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             : }