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 <string>
       8             : 
       9             : #include "DataStructures/Tensor/Tensor.hpp"
      10             : #include "Domain/CoordinateMaps/Distribution.hpp"
      11             : #include "Domain/Tags.hpp"
      12             : #include "Options/Auto.hpp"
      13             : #include "Options/String.hpp"
      14             : #include "Utilities/Gsl.hpp"
      15             : #include "Utilities/TMPL.hpp"
      16             : 
      17             : /// \cond
      18             : class DataVector;
      19             : /// \endcond
      20             : 
      21             : namespace domain {
      22             : 
      23             : /// @{
      24             : /*!
      25             :  * \ingroup ComputationalDomainGroup
      26             :  * \brief Coordinates suitable for visualizing large radii by compressing them
      27             :  * logarithmically or inversely.
      28             :  *
      29             :  * Rescales the coordinates $\boldsymbol{x}$ as
      30             :  *
      31             :  * \begin{equation}
      32             :  * \hat{\boldsymbol{x}} = \frac{\hat{r}}{r} \boldsymbol{x}
      33             :  * \text{,}
      34             :  * \end{equation}
      35             :  *
      36             :  * for $r > r_0$, where $r=\sqrt{x^2+y^2+z^2}$ is the Euclidean coordinate
      37             :  * radius and $r_0$ is the `inner_radius`.
      38             :  * The coordinates are compressed from $r \in [r_0, r_1]$ to $\hat{r} \in [r_0,
      39             :  * \hat{r}_1]$, where the `outer_radius` $r_1$ can be incomprehensibly large
      40             :  * like $10^9$ and the compressed outer radius $\hat{r}_1$ is reasonably small
      41             :  * so it can be visualized well. We choose
      42             :  *
      43             :  * \begin{equation}
      44             :  * \hat{r}_1 = r_0 \log_{10}(r_1)
      45             :  * \end{equation}
      46             :  *
      47             :  * so the compressed outer radius is a multiple of the inner radius and
      48             :  * increases with the outer radius as well, but exponentials are tamed.
      49             :  *
      50             :  * The radial compression map $\hat{r}(r)$ is just the inverse of the
      51             :  * `domain::CoordinateMaps::Interval` map, which is also used to distribute grid
      52             :  * points radially. Therefore, radial grid points will be distributed linearly
      53             :  * in the radially compressed coordinates if you use the same `compression`
      54             :  * distribution that you used to distribute radial grid points in the
      55             :  * `CoordsFrame`.
      56             :  *
      57             :  * \see domain::CoordinateMaps::Interval
      58             :  */
      59             : template <typename DataType, size_t Dim, typename CoordsFrame>
      60           1 : void radially_compressed_coordinates(
      61             :     gsl::not_null<tnsr::I<DataType, Dim, CoordsFrame>*> result,
      62             :     const tnsr::I<DataType, Dim, CoordsFrame>& coordinates, double inner_radius,
      63             :     double outer_radius, CoordinateMaps::Distribution compression);
      64             : template <typename DataType, size_t Dim, typename CoordsFrame>
      65           1 : tnsr::I<DataType, Dim, CoordsFrame> radially_compressed_coordinates(
      66             :     const tnsr::I<DataType, Dim, CoordsFrame>& coordinates, double inner_radius,
      67             :     double outer_radius, CoordinateMaps::Distribution compression);
      68             : /// @}
      69             : 
      70             : /// Options for radially compressed coordinates
      71             : ///
      72             : /// \see radially_compressed_coordinates
      73           1 : struct RadiallyCompressedCoordinatesOptions {
      74           0 :   static constexpr Options::String help =
      75             :       "Define radially compressed coordinates for visualizing large outer "
      76             :       "radii.";
      77           0 :   struct InnerRadius {
      78           0 :     using type = double;
      79           0 :     static constexpr Options::String help =
      80             :         "Radially compressed coordinates begin at this radius, and coincide "
      81             :         "with the original coordinates for smaller radii.";
      82             :   };
      83           0 :   struct OuterRadius {
      84           0 :     using type = double;
      85           0 :     static constexpr Options::String help =
      86             :         "Outer radius of the domain which will be compressed down to a "
      87             :         "comprehensible radius, namely to r_inner * log10(r_outer).";
      88             :   };
      89           0 :   struct Compression {
      90           0 :     using type = CoordinateMaps::Distribution;
      91           0 :     static constexpr Options::String help =
      92             :         "Compression mode: 'Logarithmic' or 'Inverse'. If you use the same "
      93             :         "mode that you used to distribute radial grid points then the "
      94             :         "grid points will be distributed linearly in the radially compressed "
      95             :         "coordinates.";
      96             :   };
      97           0 :   using options = tmpl::list<InnerRadius, OuterRadius, Compression>;
      98           0 :   void pup(PUP::er& p);
      99           0 :   double inner_radius;
     100           0 :   double outer_radius;
     101           0 :   CoordinateMaps::Distribution compression;
     102             : };
     103             : 
     104             : namespace OptionTags {
     105             : 
     106           0 : struct RadiallyCompressedCoordinates {
     107           0 :   using type = Options::Auto<domain::RadiallyCompressedCoordinatesOptions,
     108             :                              Options::AutoLabel::None>;
     109           0 :   static constexpr Options::String help =
     110             :       "Define radially compressed coordinates for visualizing large outer "
     111             :       "radii.";
     112             : };
     113             : 
     114             : }  // namespace OptionTags
     115             : 
     116             : namespace Tags {
     117             : 
     118             : /// Options for radially compressed coordinates, or `std::nullopt` if none
     119             : /// were specified in the input file
     120             : ///
     121             : /// \see radially_compressed_coordinates
     122           1 : struct RadiallyCompressedCoordinatesOptions : db::SimpleTag {
     123           0 :   using type = std::optional<domain::RadiallyCompressedCoordinatesOptions>;
     124           0 :   static constexpr bool pass_metavariables = false;
     125           0 :   using option_tags = tmpl::list<OptionTags::RadiallyCompressedCoordinates>;
     126           0 :   static type create_from_options(type value) { return value; };
     127             : };
     128             : 
     129             : /*!
     130             :  * \brief Coordinates suitable for visualizing large radii by compressing them
     131             :  * logarithmically or inversely.
     132             :  *
     133             :  * The coordinate map reduces to the identity if no options for radially
     134             :  * compressed coordinates were specified in the input file.
     135             :  *
     136             :  * The template parameter `RadialDim` on the compute tag determines which
     137             :  * dimensions to compress. If `RadialDim` is `void`, all dimensions are
     138             :  * compressed radially (Cartesian coordinates). If `RadialDim` is a
     139             :  * `tmpl::size_t`, only that dimension is compressed.
     140             :  *
     141             :  * \see radially_compressed_coordinates
     142             :  */
     143             : template <size_t Dim, typename CoordsFrame>
     144           1 : struct RadiallyCompressedCoordinates : db::SimpleTag {
     145           0 :   using type = tnsr::I<DataVector, Dim, CoordsFrame>;
     146             : };
     147             : 
     148             : /// \copydoc RadiallyCompressedCoordinates
     149             : template <size_t Dim, typename CoordsFrame, typename RadialDim = void>
     150           1 : struct RadiallyCompressedCoordinatesCompute
     151             :     : RadiallyCompressedCoordinates<Dim, CoordsFrame>,
     152             :       db::ComputeTag {
     153           0 :   using base = RadiallyCompressedCoordinates<Dim, CoordsFrame>;
     154           0 :   using return_type = tnsr::I<DataVector, Dim, CoordsFrame>;
     155           0 :   using argument_tags = tmpl::list<Tags::Coordinates<Dim, CoordsFrame>,
     156             :                                    Tags::RadiallyCompressedCoordinatesOptions>;
     157           0 :   static void function(
     158             :       const gsl::not_null<tnsr::I<DataVector, Dim, CoordsFrame>*> result,
     159             :       const tnsr::I<DataVector, Dim, CoordsFrame>& coords,
     160             :       const std::optional<domain::RadiallyCompressedCoordinatesOptions>&
     161             :           options) {
     162             :     if (options.has_value()) {
     163             :       if constexpr (std::is_same_v<RadialDim, void>) {
     164             :         radially_compressed_coordinates(result, coords, options->inner_radius,
     165             :                                         options->outer_radius,
     166             :                                         options->compression);
     167             :       } else {
     168             :         *result = coords;
     169             :         // Only compress the specified radial dimension
     170             :         tnsr::I<DataVector, 1, CoordsFrame> result_radial{};
     171             :         get<0>(result_radial).set_data_ref(&get<RadialDim::value>(*result));
     172             :         radially_compressed_coordinates(
     173             :             make_not_null(&result_radial),
     174             :             tnsr::I<DataVector, 1, CoordsFrame>{get<RadialDim::value>(coords)},
     175             :             options->inner_radius, options->outer_radius, options->compression);
     176             :       }
     177             :     } else {
     178             :       *result = coords;
     179             :     }
     180             :   }
     181             : };
     182             : 
     183             : }  // namespace Tags
     184             : }  // namespace domain