Line data    Source code

       1           1 : // Distributed under the MIT License.
       2             : // See LICENSE.txt for details.
       3             : 
       4             : /// \file
       5             : /// Defines tags related to domain quantities
       6             : 
       7             : #pragma once
       8             : 
       9             : #include <array>
      10             : #include <cstddef>
      11             : #include <memory>
      12             : #include <string>
      13             : #include <unordered_map>
      14             : #include <unordered_set>
      15             : #include <vector>
      16             : 
      17             : #include "DataStructures/DataBox/Subitems.hpp"
      18             : #include "DataStructures/DataBox/Tag.hpp"
      19             : #include "DataStructures/DataBox/TagName.hpp"
      20             : #include "DataStructures/Tensor/EagerMath/Determinant.hpp"
      21             : #include "DataStructures/Tensor/EagerMath/DeterminantAndInverse.hpp"
      22             : #include "DataStructures/Tensor/TypeAliases.hpp"
      23             : #include "DataStructures/Variables.hpp"
      24             : #include "Domain/Structure/Direction.hpp"
      25             : #include "Domain/Structure/Element.hpp"
      26             : #include "Utilities/GetOutput.hpp"
      27             : #include "Utilities/Gsl.hpp"
      28             : #include "Utilities/NoSuchType.hpp"
      29             : #include "Utilities/Requires.hpp"
      30             : #include "Utilities/TMPL.hpp"
      31             : #include "Utilities/TypeTraits.hpp"
      32             : #include "Utilities/TypeTraits/IsA.hpp"
      33             : 
      34             : /// \cond
      35             : class DataVector;
      36             : template <size_t VolumeDim>
      37             : class Domain;
      38             : template <size_t VolumeDim, typename Frame>
      39             : class ElementMap;
      40             : template <size_t VolumeDim>
      41             : class Mesh;
      42             : /// \endcond
      43             : 
      44             : namespace domain {
      45             : /// \ingroup ComputationalDomainGroup
      46             : /// \brief %Tags for the domain.
      47             : namespace Tags {
      48             : /// \ingroup DataBoxTagsGroup
      49             : /// \ingroup ComputationalDomainGroup
      50             : /// The ::Element associated with the DataBox
      51             : template <size_t VolumeDim>
      52           1 : struct Element : db::SimpleTag {
      53           0 :   using type = ::Element<VolumeDim>;
      54             : };
      55             : 
      56             : /// \ingroup DataBoxTagsGroup
      57             : /// \ingroup ComputationalDomainGroup
      58             : /// \brief The computational grid of the Element in the DataBox
      59             : /// \details The corresponding interface tag uses Mesh::slice_through to compute
      60             : /// the mesh on the face of the element.
      61             : template <size_t VolumeDim>
      62           1 : struct Mesh : db::SimpleTag {
      63           0 :   using type = ::Mesh<VolumeDim>;
      64             : };
      65             : 
      66             : /// \ingroup DataBoxTagsGroup
      67             : /// \ingroup ComputationalDomainGroup
      68             : /// The coordinate map from the ElementLogical frame to the TargetFrame
      69             : template <size_t VolumeDim, typename TargetFrame = Frame::Inertial>
      70           1 : struct ElementMap : db::SimpleTag {
      71           0 :   static constexpr size_t dim = VolumeDim;
      72           0 :   using target_frame = TargetFrame;
      73           0 :   using source_frame = Frame::ElementLogical;
      74             : 
      75           0 :   static std::string name() {
      76             :     return "ElementMap(" + get_output(TargetFrame{}) + ")";
      77             :   }
      78           0 :   using type = ::ElementMap<VolumeDim, TargetFrame>;
      79             : };
      80             : 
      81             : /// \ingroup DataBoxTagsGroup
      82             : /// \ingroup ComputationalDomainGroup
      83             : /// The coordinates in a given frame.
      84             : template <size_t Dim, typename Frame>
      85           1 : struct Coordinates : db::SimpleTag {
      86           0 :   static std::string name() { return get_output(Frame{}) + "Coordinates"; }
      87           0 :   using type = tnsr::I<DataVector, Dim, Frame>;
      88             : };
      89             : 
      90             : /// \ingroup DataBoxTagsGroup
      91             : /// \ingroup ComputationalDomainGroup
      92             : /// The coordinates in the target frame of `MapTag`. The `SourceCoordsTag`'s
      93             : /// frame must be the source frame of `MapTag`
      94             : template <class MapTag, class SourceCoordsTag,
      95             :           template <size_t, class> class CoordinatesTag = Coordinates>
      96           1 : struct MappedCoordinates
      97             :     : CoordinatesTag<MapTag::dim, typename MapTag::target_frame>,
      98             :       db::ComputeTag {
      99           0 :   using base = CoordinatesTag<MapTag::dim, typename MapTag::target_frame>;
     100           0 :   using return_type = typename base::type;
     101           0 :   using argument_tags = tmpl::list<MapTag, SourceCoordsTag>;
     102           0 :   static constexpr auto function(
     103             :       const gsl::not_null<return_type*> target_coords,
     104             :       const typename MapTag::type& element_map,
     105             :       const tnsr::I<DataVector, MapTag::dim, typename MapTag::source_frame>&
     106             :           source_coords) {
     107             :     *target_coords = element_map(source_coords);
     108             :   }
     109             : };
     110             : 
     111             : /// \ingroup DataBoxTagsGroup
     112             : /// \ingroup ComputationalDomainGroup
     113             : /// \brief The inverse Jacobian from the source frame to the target frame.
     114             : ///
     115             : /// Specifically, \f$\partial x^{\bar{i}} / \partial x^i\f$, where \f$\bar{i}\f$
     116             : /// denotes the source frame and \f$i\f$ denotes the target frame.
     117             : template <size_t Dim, typename SourceFrame, typename TargetFrame>
     118           1 : struct InverseJacobian : db::SimpleTag {
     119           0 :   static std::string name() {
     120             :     return "InverseJacobian(" + get_output(SourceFrame{}) + "," +
     121             :            get_output(TargetFrame{}) + ")";
     122             :   }
     123           0 :   using type = ::InverseJacobian<DataVector, Dim, SourceFrame, TargetFrame>;
     124             : };
     125             : 
     126             : /// \ingroup DataBoxTagsGroup
     127             : /// \ingroup ComputationalDomainGroup
     128             : /// Computes the inverse Jacobian of the map held by `MapTag` at the coordinates
     129             : /// held by `SourceCoordsTag`. The coordinates must be in the source frame of
     130             : /// the map.
     131             : template <typename MapTag, typename SourceCoordsTag>
     132           1 : struct InverseJacobianCompute
     133             :     : InverseJacobian<MapTag::dim, typename MapTag::source_frame,
     134             :                       typename MapTag::target_frame>,
     135             :       db::ComputeTag {
     136           0 :   using base = InverseJacobian<MapTag::dim, typename MapTag::source_frame,
     137             :                                typename MapTag::target_frame>;
     138           0 :   using return_type = typename base::type;
     139           0 :   using argument_tags = tmpl::list<MapTag, SourceCoordsTag>;
     140           0 :   static constexpr auto function(
     141             :       const gsl::not_null<return_type*> inv_jacobian,
     142             :       const typename MapTag::type& element_map,
     143             :       const tnsr::I<DataVector, MapTag::dim, typename MapTag::source_frame>&
     144             :           source_coords) {
     145             :     *inv_jacobian = element_map.inv_jacobian(source_coords);
     146             :   }
     147             : };
     148             : 
     149             : /// \ingroup DataBoxTagsGroup
     150             : /// \ingroup ComputationalDomainGroup
     151             : /// \brief The Jacobian from the source frame to the target frame.
     152             : ///
     153             : /// Specifically, \f$\partial x^{i} / \partial \xi^{\bar{i}}\f$, where
     154             : /// \f$\xi^\bar{i}\f$ denotes the source frame and \f$x^i\f$ denotes the target
     155             : /// frame.
     156             : template <size_t Dim, typename SourceFrame, typename TargetFrame>
     157           1 : struct Jacobian : db::SimpleTag {
     158           0 :   static std::string name() {
     159             :     return "Jacobian(" + get_output(SourceFrame{}) + "," +
     160             :            get_output(TargetFrame{}) + ")";
     161             :   }
     162           0 :   using type = ::Jacobian<DataVector, Dim, SourceFrame, TargetFrame>;
     163             : };
     164             : 
     165             : /// \ingroup DataBoxTagsGroup
     166             : /// \ingroup ComputationalDomainGroup
     167             : /// Computes the Jacobian of the map from the `InverseJacobian<Dim, SourceFrame,
     168             : /// TargetFrame>` tag.
     169             : template <size_t Dim, typename SourceFrame, typename TargetFrame>
     170           1 : struct JacobianCompute : Jacobian<Dim, SourceFrame, TargetFrame>,
     171             :                          db::ComputeTag {
     172           0 :   using base = Jacobian<Dim, SourceFrame, TargetFrame>;
     173           0 :   using return_type = typename base::type;
     174           0 :   using argument_tags =
     175             :       tmpl::list<InverseJacobian<Dim, SourceFrame, TargetFrame>>;
     176           0 :   static constexpr auto function(
     177             :       const gsl::not_null<return_type*> jacobian,
     178             :       const ::InverseJacobian<DataVector, Dim, SourceFrame, TargetFrame>&
     179             :           inv_jac) {
     180             :     *jacobian = determinant_and_inverse(inv_jac).second;
     181             :   }
     182             : };
     183             : 
     184             : /// \ingroup DataBoxTagsGroup
     185             : /// \ingroup ComputationalDomainGroup
     186             : /// \brief The determinant of the inverse Jacobian from the source frame to the
     187             : /// target frame.
     188             : template <typename SourceFrame, typename TargetFrame>
     189           1 : struct DetInvJacobian : db::SimpleTag {
     190           0 :   using type = Scalar<DataVector>;
     191           0 :   static std::string name() {
     192             :     return "DetInvJacobian(" + get_output(SourceFrame{}) + "," +
     193             :            get_output(TargetFrame{}) + ")";
     194             :   }
     195             : };
     196             : 
     197             : /// \ingroup DataBoxTagsGroup
     198             : /// \ingroup ComputationalDomainGroup
     199             : /// Computes the determinant of the inverse Jacobian.
     200             : template <size_t Dim, typename SourceFrame, typename TargetFrame>
     201           1 : struct DetInvJacobianCompute : db::ComputeTag,
     202             :                                DetInvJacobian<SourceFrame, TargetFrame> {
     203           0 :   using base = DetInvJacobian<SourceFrame, TargetFrame>;
     204           0 :   using return_type = typename base::type;
     205           0 :   using argument_tags =
     206             :       tmpl::list<InverseJacobian<Dim, SourceFrame, TargetFrame>>;
     207           0 :   static void function(const gsl::not_null<return_type*> det_inv_jac,
     208             :                        const ::InverseJacobian<DataVector, Dim, SourceFrame,
     209             :                                                TargetFrame>& inv_jac) {
     210             :     determinant(det_inv_jac, inv_jac);
     211             :   }
     212             : };
     213             : 
     214             : /// \ingroup DataBoxTagsGroup
     215             : /// \ingroup ComputationalDomainGroup
     216             : /// \brief The determinant of the Jacobian from the source frame to the target
     217             : /// frame.
     218             : template <typename SourceFrame, typename TargetFrame>
     219           1 : struct DetJacobian : db::SimpleTag {
     220           0 :   using type = Scalar<DataVector>;
     221           0 :   static std::string name() {
     222             :     return "DetJacobian(" + get_output(SourceFrame{}) + "," +
     223             :            get_output(TargetFrame{}) + ")";
     224             :   }
     225             : };
     226             : 
     227             : /// \ingroup DataBoxTagsGroup
     228             : /// \ingroup ComputationalDomainGroup
     229             : /// \brief The inverse Jacobian times the determinant of the Jacobian.
     230             : ///
     231             : /// This quantity is divergence-free analytically. See
     232             : /// `::dg::metric_identity_det_jac_times_inv_jac` for more information.
     233             : template <size_t Dim, typename SourceFrame, typename TargetFrame>
     234           1 : struct DetTimesInvJacobian : db::SimpleTag {
     235           0 :   using type = ::InverseJacobian<DataVector, Dim, SourceFrame, TargetFrame>;
     236           0 :   static std::string name() {
     237             :     return "DetTimesInvJacobian(" + get_output(SourceFrame{}) + "," +
     238             :            get_output(TargetFrame{}) + ")";
     239             :   }
     240             : };
     241             : 
     242             : /// \ingroup DataBoxTagsGroup
     243             : /// \ingroup ComputationalDomainGroup
     244             : /// Base tag for boundary data needed for updating the variables.
     245           1 : struct VariablesBoundaryData : db::BaseTag {};
     246             : 
     247             : /// @{
     248             : /// \ingroup DataBoxTagsGroup
     249             : /// \ingroup ComputationalDomainGroup
     250             : /// The set of directions to neighboring Elements
     251             : template <size_t VolumeDim>
     252           1 : struct InternalDirections : db::SimpleTag {
     253           0 :   static constexpr size_t volume_dim = VolumeDim;
     254           0 :   using type = std::unordered_set<::Direction<VolumeDim>>;
     255             : };
     256             : 
     257             : template <size_t VolumeDim>
     258           0 : struct InternalDirectionsCompute : InternalDirections<VolumeDim>,
     259             :                                    db::ComputeTag {
     260           0 :   static constexpr size_t volume_dim = VolumeDim;
     261           0 :   using base = InternalDirections<VolumeDim>;
     262           0 :   using return_type = std::unordered_set<::Direction<VolumeDim>>;
     263           0 :   using argument_tags = tmpl::list<Element<VolumeDim>>;
     264           0 :   static void function(const gsl::not_null<return_type*> directions,
     265             :                        const ::Element<VolumeDim>& element) {
     266             :     for (const auto& direction_neighbors : element.neighbors()) {
     267             :       directions->insert(direction_neighbors.first);
     268             :     }
     269             :   }
     270             : };
     271             : /// @}
     272             : 
     273             : /// @{
     274             : /// \ingroup DataBoxTagsGroup
     275             : /// \ingroup ComputationalDomainGroup
     276             : /// The set of directions which correspond to external boundaries.
     277             : /// Used for representing data on the interior side of the external boundary
     278             : /// faces.
     279             : template <size_t VolumeDim>
     280           1 : struct BoundaryDirectionsInterior : db::SimpleTag {
     281           0 :   static constexpr size_t volume_dim = VolumeDim;
     282           0 :   using type = std::unordered_set<::Direction<VolumeDim>>;
     283             : };
     284             : 
     285             : template <size_t VolumeDim>
     286           0 : struct BoundaryDirectionsInteriorCompute
     287             :     : BoundaryDirectionsInterior<VolumeDim>,
     288             :       db::ComputeTag {
     289           0 :   static constexpr size_t volume_dim = VolumeDim;
     290           0 :   using base = BoundaryDirectionsInterior<VolumeDim>;
     291           0 :   using return_type = std::unordered_set<::Direction<VolumeDim>>;
     292           0 :   using argument_tags = tmpl::list<Element<VolumeDim>>;
     293           0 :   static void function(const gsl::not_null<return_type*> directions,
     294             :                        const ::Element<VolumeDim>& element) {
     295             :     *directions = element.external_boundaries();
     296             :   }
     297             : };
     298             : /// @}
     299             : 
     300             : /// @{
     301             : /// \ingroup DataBoxTagsGroup
     302             : /// \ingroup ComputationalDomainGroup
     303             : /// The set of directions which correspond to external boundaries. To be used
     304             : /// to represent data which exists on the exterior side of the external boundary
     305             : /// faces.
     306             : template <size_t VolumeDim>
     307           1 : struct BoundaryDirectionsExterior : db::SimpleTag {
     308           0 :   static constexpr size_t volume_dim = VolumeDim;
     309           0 :   using type = std::unordered_set<::Direction<VolumeDim>>;
     310             : };
     311             : 
     312             : template <size_t VolumeDim>
     313           0 : struct BoundaryDirectionsExteriorCompute
     314             :     : BoundaryDirectionsExterior<VolumeDim>,
     315             :       db::ComputeTag {
     316           0 :   static constexpr size_t volume_dim = VolumeDim;
     317           0 :   using base = BoundaryDirectionsExterior<VolumeDim>;
     318           0 :   using return_type = std::unordered_set<::Direction<VolumeDim>>;
     319           0 :   using argument_tags = tmpl::list<Element<VolumeDim>>;
     320           0 :   static constexpr auto function(const gsl::not_null<return_type*> directions,
     321             :                                  const ::Element<VolumeDim>& element) {
     322             :     *directions = element.external_boundaries();
     323             :   }
     324             : };
     325             : /// @}
     326             : 
     327             : /// \ingroup DataBoxTagsGroup
     328             : /// \ingroup ComputationalDomainGroup
     329             : /// \brief Tag which is either a SimpleTag for quantities on an
     330             : /// interface, base tag to a compute item which acts on tags on an interface, or
     331             : /// base tag to a compute item which slices a tag from the volume to an
     332             : /// interface.
     333             : ///
     334             : /// The contained object will be a map from ::Direction to the item type of
     335             : /// `Tag`, with the set of directions being those produced by `DirectionsTag`.
     336             : ///
     337             : /// If a SimpleTag is desired on the interface, then this tag can be added to
     338             : /// the DataBox directly. If a ComputeTag which acts on Tags on the interface is
     339             : /// desired, then the tag should be added using `InterfaceCompute`. If a
     340             : /// ComputeTag which slices a TensorTag or a VariablesTag in the volume to an
     341             : /// Interface is desired, then it should be added using `Slice`. In all cases,
     342             : /// the tag can then be retrieved using `Tags::Interface<DirectionsTag, Tag>`.
     343             : ///
     344             : /// \tparam DirectionsTag the item of directions
     345             : /// \tparam Tag the tag labeling the item
     346             : ///
     347             : /// \see InterfaceCompute, Slice
     348             : template <typename DirectionsTag, typename Tag>
     349           1 : struct Interface : db::SimpleTag {
     350             :   static_assert(db::is_simple_tag_v<DirectionsTag>);
     351             :   static_assert(db::is_simple_tag_v<Tag>);
     352           0 :   static std::string name() {
     353             :     return "Interface<" + db::tag_name<DirectionsTag>() + ", " +
     354             :            db::tag_name<Tag>() + ">";
     355             :   };
     356           0 :   using tag = Tag;
     357           0 :   using type = std::unordered_map<::Direction<DirectionsTag::volume_dim>,
     358             :                                   typename Tag::type>;
     359             : };
     360             : 
     361             : /// \ingroup DataBoxTagsGroup
     362             : /// \ingroup ComputationalDomainGroup
     363             : /// ::Direction to an interface
     364             : template <size_t VolumeDim>
     365           1 : struct Direction : db::SimpleTag {
     366           0 :   using type = ::Direction<VolumeDim>;
     367             : };
     368             : 
     369             : }  // namespace Tags
     370             : }  // namespace domain
     371             : 
     372             : namespace db {
     373             : namespace detail {
     374             : template <typename DirectionsTag, typename VariablesTag>
     375             : struct InterfaceSubitemsImpl {
     376             :   using type = tmpl::transform<
     377             :       typename VariablesTag::type::tags_list,
     378             :       tmpl::bind<domain::Tags::Interface, tmpl::pin<DirectionsTag>, tmpl::_1>>;
     379             : 
     380             :   using tag = domain::Tags::Interface<DirectionsTag, VariablesTag>;
     381             : 
     382             :   template <typename Subtag>
     383             :   static void create_item(
     384             :       const gsl::not_null<typename tag::type*> parent_value,
     385             :       const gsl::not_null<typename Subtag::type*> sub_value) {
     386             :     sub_value->clear();
     387             :     for (auto& direction_vars : *parent_value) {
     388             :       const auto& direction = direction_vars.first;
     389             :       auto& parent_vars = get<typename Subtag::tag>(direction_vars.second);
     390             :       auto& sub_var = (*sub_value)[direction];
     391             :       for (auto vars_it = parent_vars.begin(), sub_var_it = sub_var.begin();
     392             :            vars_it != parent_vars.end(); ++vars_it, ++sub_var_it) {
     393             :         sub_var_it->set_data_ref(&*vars_it);
     394             :       }
     395             :     }
     396             :   }
     397             : 
     398             :   // The `return_type` can be anything for Subitems because the DataBox figures
     399             :   // out the correct return type, we just use the `return_type` type alias to
     400             :   // signal to the DataBox we want mutating behavior.
     401             :   using return_type = NoSuchType;
     402             : 
     403             :   template <typename Subtag>
     404             :   static void create_compute_item(
     405             :       const gsl::not_null<typename Subtag::type*> sub_value,
     406             :       const typename tag::type& parent_value) {
     407             :     for (const auto& direction_vars : parent_value) {
     408             :       const auto& direction = direction_vars.first;
     409             :       const auto& parent_vars =
     410             :           get<typename Subtag::tag>(direction_vars.second);
     411             :       auto& sub_var = (*sub_value)[direction];
     412             :       auto sub_var_it = sub_var.begin();
     413             :       for (auto vars_it = parent_vars.begin();
     414             :            vars_it != parent_vars.end(); ++vars_it, ++sub_var_it) {
     415             :         // clang-tidy: do not use const_cast
     416             :         // The DataBox will only give out a const reference to the
     417             :         // result of a compute item.  Here, that is a reference to a
     418             :         // const map to Tensors of DataVectors.  There is no (publicly
     419             :         // visible) indirection there, so having the map const will
     420             :         // allow only allow const access to the contained DataVectors,
     421             :         // so no modification through the pointer cast here is
     422             :         // possible.
     423             :         sub_var_it->set_data_ref(const_cast<DataVector*>(&*vars_it));  // NOLINT
     424             :       }
     425             :     }
     426             :   }
     427             : };
     428             : }  // namespace detail
     429             : 
     430             : template <typename DirectionsTag, typename VariablesTag>
     431           0 : struct Subitems<domain::Tags::Interface<DirectionsTag, VariablesTag>,
     432             :                 Requires<tt::is_a_v<Variables, typename VariablesTag::type>>>
     433             :     : detail::InterfaceSubitemsImpl<DirectionsTag, VariablesTag> {};
     434             : 
     435             : }  // namespace db