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/DataBox/Tag.hpp"
      10             : #include "DataStructures/Tensor/EagerMath/Determinant.hpp"
      11             : #include "DataStructures/Tensor/TypeAliases.hpp"
      12             : #include "Domain/ElementMap.hpp"
      13             : #include "Evolution/DgSubcell/Tags/Coordinates.hpp"
      14             : #include "Utilities/SetNumberOfGridPoints.hpp"
      15             : 
      16             : /// \cond
      17             : namespace Tags {
      18             : struct Time;
      19             : }  // namespace Tags
      20             : /// \endcond
      21             : 
      22           1 : namespace evolution::dg::subcell::fd::Tags {
      23             : /// \brief The inverse Jacobian from the element logical frame to the grid frame
      24             : /// at the cell centers.
      25             : ///
      26             : /// Specifically, \f$\partial x^{\bar{i}} / \partial x^i\f$, where \f$\bar{i}\f$
      27             : /// denotes the element logical frame and \f$i\f$ denotes the grid frame.
      28             : template <size_t Dim>
      29           1 : struct InverseJacobianLogicalToGrid : db::SimpleTag {
      30           0 :   static std::string name() { return "InverseJacobian(Logical,Grid)"; }
      31           0 :   using type =
      32             :       ::InverseJacobian<DataVector, Dim, Frame::ElementLogical, Frame::Grid>;
      33             : };
      34             : 
      35             : /// \brief The determinant of the inverse Jacobian from the element logical
      36             : /// frame to the grid frame at the cell centers.
      37           1 : struct DetInverseJacobianLogicalToGrid : db::SimpleTag {
      38           0 :   static std::string name() { return "Det(InverseJacobian(Logical,Grid))"; }
      39           0 :   using type = Scalar<DataVector>;
      40             : };
      41             : 
      42             : /// \brief The determinant of the inverse Jacobian from the element logical
      43             : /// frame to the inertial frame at the cell centers.
      44           1 : struct DetInverseJacobianLogicalToInertial : db::SimpleTag {
      45           0 :   static std::string name() { return "Det(InverseJacobian(Logical,Inertial))"; }
      46           0 :   using type = Scalar<DataVector>;
      47             : };
      48             : 
      49             : /// \brief The inverse Jacobian from the element logical frame to the inertial
      50             : /// frame at the cell centers.
      51             : ///
      52             : /// Specifically, \f$\partial x^{\bar{i}} / \partial x^i\f$, where \f$\bar{i}\f$
      53             : /// denotes the element logical frame and \f$i\f$ denotes the inertial frame.
      54             : template <size_t Dim>
      55           1 : struct InverseJacobianLogicalToInertial : db::SimpleTag {
      56           0 :   static std::string name() { return "InverseJacobian(Logical,Inertial)"; }
      57           0 :   using type = ::InverseJacobian<DataVector, Dim, Frame::ElementLogical,
      58             :                                  Frame::Inertial>;
      59             : };
      60             : 
      61             : /// Compute item for the inverse jacobian matrix from logical to
      62             : /// grid coordinates
      63             : template <typename MapTagLogicalToGrid, size_t Dim>
      64           1 : struct InverseJacobianLogicalToGridCompute : InverseJacobianLogicalToGrid<Dim>,
      65             :                                              db::ComputeTag {
      66           0 :   static constexpr size_t dim = Dim;
      67           0 :   using base = InverseJacobianLogicalToGrid<Dim>;
      68           0 :   using return_type = typename base::type;
      69           0 :   using argument_tags = tmpl::list<
      70             :       MapTagLogicalToGrid,
      71             :       evolution::dg::subcell::Tags::Coordinates<dim, Frame::ElementLogical>>;
      72           0 :   static void function(
      73             :       const gsl::not_null<return_type*> inverse_jacobian_logical_to_grid,
      74             :       const ElementMap<Dim, Frame::Grid>& logical_to_grid_map,
      75             :       const tnsr::I<DataVector, dim, Frame::ElementLogical>& logical_coords) {
      76             :     *inverse_jacobian_logical_to_grid =
      77             :         logical_to_grid_map.inv_jacobian(logical_coords);
      78             :   }
      79             : };
      80             : 
      81             : /// Compute item for the determinant of the inverse jacobian matrix
      82             : /// from logical to grid coordinates
      83             : template <size_t Dim>
      84           1 : struct DetInverseJacobianLogicalToGridCompute
      85             :     : DetInverseJacobianLogicalToGrid,
      86             :       db::ComputeTag {
      87           0 :   static constexpr size_t dim = Dim;
      88           0 :   using base = DetInverseJacobianLogicalToGrid;
      89           0 :   using return_type = typename base::type;
      90           0 :   using argument_tags = tmpl::list<InverseJacobianLogicalToGrid<Dim>>;
      91           0 :   static void function(
      92             :       const gsl::not_null<return_type*> det_inverse_jacobian_logical_to_grid,
      93             :       const ::InverseJacobian<DataVector, Dim, Frame::ElementLogical,
      94             :                               Frame::Grid>& inverse_jacobian_logical_to_grid) {
      95             :     *det_inverse_jacobian_logical_to_grid =
      96             :         determinant(inverse_jacobian_logical_to_grid);
      97             :   }
      98             : };
      99             : 
     100             : /// Compute item for the inverse jacobian matrix from logical to
     101             : /// inertial coordinates
     102             : template <typename MapTagGridToInertial, size_t Dim>
     103           1 : struct InverseJacobianLogicalToInertialCompute
     104             :     : InverseJacobianLogicalToInertial<Dim>,
     105             :       db::ComputeTag {
     106           0 :   static constexpr size_t dim = Dim;
     107           0 :   using base = InverseJacobianLogicalToInertial<Dim>;
     108           0 :   using return_type = typename base::type;
     109           0 :   using argument_tags =
     110             :       tmpl::list<MapTagGridToInertial,
     111             :                  evolution::dg::subcell::Tags::Coordinates<dim, Frame::Grid>,
     112             :                  ::Tags::Time, ::domain::Tags::FunctionsOfTime,
     113             :                  InverseJacobianLogicalToGrid<Dim>>;
     114           0 :   static void function(
     115             :       const gsl::not_null<return_type*> inverse_jacobian_logical_to_inertial,
     116             :       const ::domain::CoordinateMapBase<Frame::Grid, Frame::Inertial, dim>&
     117             :           grid_to_inertial_map,
     118             :       const tnsr::I<DataVector, dim, Frame::Grid>& grid_coords,
     119             :       const double time,
     120             :       const std::unordered_map<
     121             :           std::string,
     122             :           std::unique_ptr<::domain::FunctionsOfTime::FunctionOfTime>>&
     123             :           functions_of_time,
     124             :       const ::InverseJacobian<DataVector, Dim, Frame::ElementLogical,
     125             :                               Frame::Grid>& inverse_jacobian_logical_to_grid) {
     126             :     if (grid_to_inertial_map.is_identity()) {
     127             :       // Optimization for time-independent maps; we just point to the
     128             :       // logical-to-grid inverse jacobian.
     129             :       const size_t num_pts = inverse_jacobian_logical_to_grid[0].size();
     130             :       for (size_t storage_index = 0;
     131             :            storage_index < inverse_jacobian_logical_to_grid.size();
     132             :            ++storage_index) {
     133             :         make_const_view(
     134             :             make_not_null(&std::as_const(
     135             :                 (*inverse_jacobian_logical_to_inertial)[storage_index])),
     136             :             inverse_jacobian_logical_to_grid[storage_index], 0, num_pts);
     137             :       }
     138             :     } else {
     139             :       set_number_of_grid_points(inverse_jacobian_logical_to_inertial,
     140             :                                 get<0>(grid_coords).size());
     141             :       // Get grid to inertial inverse jacobian
     142             :       const auto& inverse_jacobian_grid_to_inertial =
     143             :           grid_to_inertial_map.inv_jacobian(grid_coords, time,
     144             :                                             functions_of_time);
     145             :       for (size_t i = 0; i < Dim; i++) {
     146             :         for (size_t j = 0; j < Dim; j++) {
     147             :           auto& inv_jacobian_component =
     148             :               inverse_jacobian_logical_to_inertial->get(i, j);
     149             :           inv_jacobian_component = 0.;
     150             :           for (size_t k = 0; k < Dim; k++) {
     151             :             inv_jacobian_component +=
     152             :                 inverse_jacobian_logical_to_grid.get(i, k) *
     153             :                 inverse_jacobian_grid_to_inertial.get(k, j);
     154             :           }
     155             :         }
     156             :       }
     157             :     }
     158             :   }
     159             : };
     160             : 
     161             : /// Compute item for the determinant of the inverse jacobian matrix
     162             : /// from logical to inertial coordinates
     163             : template <typename MapTagGridToInertial, size_t Dim>
     164           1 : struct DetInverseJacobianLogicalToInertialCompute
     165             :     : DetInverseJacobianLogicalToInertial,
     166             :       db::ComputeTag {
     167           0 :   static constexpr size_t dim = Dim;
     168           0 :   using base = DetInverseJacobianLogicalToInertial;
     169           0 :   using return_type = typename base::type;
     170           0 :   using argument_tags =
     171             :       tmpl::list<MapTagGridToInertial, InverseJacobianLogicalToInertial<Dim>,
     172             :                  DetInverseJacobianLogicalToGrid>;
     173           0 :   static void function(
     174             :       const gsl::not_null<return_type*>
     175             :           det_inverse_jacobian_logical_to_inertial,
     176             :       const ::domain::CoordinateMapBase<Frame::Grid, Frame::Inertial, dim>&
     177             :           grid_to_inertial_map,
     178             :       const ::InverseJacobian<DataVector, Dim, Frame::ElementLogical,
     179             :                               Frame::Inertial>&
     180             :           inverse_jacobian_logical_to_inertial,
     181             :       const Scalar<DataVector>& det_inverse_jacobian_logical_to_grid) {
     182             :     if (grid_to_inertial_map.is_identity()) {
     183             :       const size_t num_pts = get(det_inverse_jacobian_logical_to_grid).size();
     184             :       make_const_view(make_not_null(&std::as_const(
     185             :                           get(*det_inverse_jacobian_logical_to_inertial))),
     186             :                       get(det_inverse_jacobian_logical_to_grid), 0, num_pts);
     187             :     } else {
     188             :       *det_inverse_jacobian_logical_to_inertial =
     189             :           determinant(inverse_jacobian_logical_to_inertial);
     190             :     }
     191             :   }
     192             : };
     193             : }  // namespace evolution::dg::subcell::fd::Tags