Line data    Source code

       1           0 : // Distributed under the MIT License.
       2             : // See LICENSE.txt for details.
       3             : 
       4             : #pragma once
       5             : 
       6             : #include <memory>
       7             : 
       8             : #include "DataStructures/Tensor/TypeAliases.hpp"
       9             : #include "Evolution/Systems/ForceFree/Tags.hpp"
      10             : #include "Options/Context.hpp"
      11             : #include "Options/String.hpp"
      12             : #include "PointwiseFunctions/AnalyticSolutions/AnalyticSolution.hpp"
      13             : #include "PointwiseFunctions/AnalyticSolutions/GeneralRelativity/Minkowski.hpp"
      14             : #include "PointwiseFunctions/InitialDataUtilities/InitialData.hpp"
      15             : #include "Utilities/Serialization/CharmPupable.hpp"
      16             : #include "Utilities/TMPL.hpp"
      17             : #include "Utilities/TaggedTuple.hpp"
      18             : 
      19             : /// \cond
      20             : class DataVector;
      21             : // IWYU pragma: no_forward_declare Tensor
      22             : namespace PUP {
      23             : class er;
      24             : }  // namespace PUP
      25             : /// \endcond
      26             : 
      27           1 : namespace ForceFree::Solutions {
      28             : 
      29             : /*!
      30             :  * \brief Alfven wave propagating along \f$x\f$ direction in flat spacetime with
      31             :  * the wave speed \f$\mu\f$.
      32             :  *
      33             :  * This test problem was introduced in \cite Komissarov2004 with $\mu=0$.
      34             :  *
      35             :  * In the wave frame (with prime superscript), the stationary solution is given
      36             :  * by
      37             :  *
      38             :  * \f{align*}
      39             :  *  B'_x & = B'_y = 1.0 ,\\
      40             :  *  B'_z(x') & = \left\{\begin{array}{ll}
      41             :  *      1.0         & \text{if } x' < -0.1   \\
      42             :  *      1.15 + 0.15 \sin (5\pi x') & \text{if } -0.1 < x' < 0.1 \\
      43             :  *      1.3         & \text{if } x' > 0.1 \end{array}\right\} ,\\
      44             :  *      E'_x & = -B'_z ,\\
      45             :  *      E'_y & = 0 \\
      46             :  *      E'_z & = 1.0
      47             :  * \f}
      48             :  *
      49             :  * and
      50             :  *
      51             :  * \f{align*}
      52             :  *  q'   & = \left\{\begin{array}{ll}
      53             :  *      - 0.75 \pi \cos (5 \pi x')  & \text{if } -0.1 < x' < 0.1 \\
      54             :  *      0         & \text{otherwise}    \end{array}\right\} , \\
      55             :  *  J_x' & = 0 , \\
      56             :  *  J_y' & = \left\{\begin{array}{ll}
      57             :  *      - 0.75 \pi \cos (5 \pi x')  & \text{if } -0.1 < x' < 0.1 \\
      58             :  *      0         & \text{otherwise} \end{array}\right\} , \\
      59             :  *  J_z' & = 0 .
      60             :  * \f}
      61             :  *
      62             :  * Applying the Lorentz transformation, electromagnetic fields and 4-current in
      63             :  * the grid frame at \f$t=0\f$ are given by
      64             :  *
      65             :  * \f{align*}
      66             :  *  E_x(x) & = E'_x(\gamma x) , \\
      67             :  *  E_y(x) & = \gamma[E'_y(\gamma x) + \mu B'_z(\gamma x)] , \\
      68             :  *  E_z(x) & = \gamma[E'_z(\gamma x) - \mu B'_y(\gamma x)] , \\
      69             :  *  B_x(x) & = B'_x(\gamma x), \\
      70             :  *  B_y(x) & = \gamma[ B'_y(\gamma x) - \mu E'_z(\gamma x) ] , \\
      71             :  *  B_z(x) & = \gamma[ B'_z(\gamma x) + \mu E'_y(\gamma x) ] .
      72             :  * \f}
      73             :  *
      74             :  * and
      75             :  *
      76             :  * \f{align*}
      77             :  *  q(x)    & = \gamma q'(\gamma x) , \\
      78             :  *  J_x(x)  & = \gamma \mu q'(\gamma x) , \\
      79             :  *  J_y(x)  & = J_y'(\gamma x) , \\
      80             :  *  J_z(x)  & = 0 .
      81             :  * \f}
      82             :  *
      83             :  * The wave speed can be chosen any value $-1 < \mu < 1$, and the solution at
      84             :  * time $t$ is $f(x,t) = f(x-\mu t, 0)$ for any physical quantities.
      85             :  *
      86             :  */
      87           1 : class AlfvenWave : public evolution::initial_data::InitialData,
      88             :                    public MarkAsAnalyticSolution {
      89             :  public:
      90             :   /// The wave speed
      91           1 :   struct WaveSpeed {
      92           0 :     using type = double;
      93           0 :     static constexpr Options::String help = {
      94             :         "The wave speed along x direction"};
      95           0 :     static type lower_bound() { return -1.0; }
      96           0 :     static type upper_bound() { return 1.0; }
      97             :   };
      98             : 
      99           0 :   using options = tmpl::list<WaveSpeed>;
     100             : 
     101           0 :   static constexpr Options::String help{
     102             :       "Alfven wave propagating along x direction in flat spacetime with the "
     103             :       "wave speed mu"};
     104             : 
     105           0 :   AlfvenWave() = default;
     106           0 :   AlfvenWave(const AlfvenWave&) = default;
     107           0 :   AlfvenWave& operator=(const AlfvenWave&) = default;
     108           0 :   AlfvenWave(AlfvenWave&&) = default;
     109           0 :   AlfvenWave& operator=(AlfvenWave&&) = default;
     110           0 :   ~AlfvenWave() override = default;
     111             : 
     112           0 :   AlfvenWave(double wave_speed, const Options::Context& context = {});
     113             : 
     114           0 :   auto get_clone() const
     115             :       -> std::unique_ptr<evolution::initial_data::InitialData> override;
     116             : 
     117             :   /// \cond
     118             :   explicit AlfvenWave(CkMigrateMessage* msg);
     119             :   using PUP::able::register_constructor;
     120             :   WRAPPED_PUPable_decl_template(AlfvenWave);
     121             :   /// \endcond
     122             : 
     123             :   // NOLINTNEXTLINE(google-runtime-references)
     124           0 :   void pup(PUP::er& p) override;
     125             : 
     126             :   /// @{
     127             :   /// Retrieve the EM variables at (x,t).
     128           1 :   auto variables(const tnsr::I<DataVector, 3>& x, double t,
     129             :                  tmpl::list<Tags::TildeE> /*meta*/) const
     130             :       -> tuples::TaggedTuple<Tags::TildeE>;
     131             : 
     132           1 :   auto variables(const tnsr::I<DataVector, 3>& x, double t,
     133             :                  tmpl::list<Tags::TildeB> /*meta*/) const
     134             :       -> tuples::TaggedTuple<Tags::TildeB>;
     135             : 
     136           1 :   static auto variables(const tnsr::I<DataVector, 3>& x, double t,
     137             :                         tmpl::list<Tags::TildePsi> /*meta*/)
     138             :       -> tuples::TaggedTuple<Tags::TildePsi>;
     139             : 
     140           1 :   static auto variables(const tnsr::I<DataVector, 3>& x, double t,
     141             :                         tmpl::list<Tags::TildePhi> /*meta*/)
     142             :       -> tuples::TaggedTuple<Tags::TildePhi>;
     143             : 
     144           1 :   auto variables(const tnsr::I<DataVector, 3>& x, double t,
     145             :                  tmpl::list<Tags::TildeQ> /*meta*/) const
     146             :       -> tuples::TaggedTuple<Tags::TildeQ>;
     147             :   /// @}
     148             : 
     149             :   /// Retrieve a collection of EM variables at `(x, t)`
     150             :   template <typename... Tags>
     151           1 :   tuples::TaggedTuple<Tags...> variables(const tnsr::I<DataVector, 3>& x,
     152             :                                          const double t,
     153             :                                          tmpl::list<Tags...> /*meta*/) const {
     154             :     static_assert(sizeof...(Tags) > 1,
     155             :                   "The generic template will recurse infinitely if only one "
     156             :                   "tag is being retrieved.");
     157             :     return {get<Tags>(variables(x, t, tmpl::list<Tags>{}))...};
     158             :   }
     159             : 
     160             :   /// Retrieve the metric variables
     161             :   template <typename Tag>
     162           1 :   tuples::TaggedTuple<Tag> variables(const tnsr::I<DataVector, 3>& x, double t,
     163             :                                      tmpl::list<Tag> /*meta*/) const {
     164             :     return background_spacetime_.variables(x, t, tmpl::list<Tag>{});
     165             :   }
     166             : 
     167             :  private:
     168           0 :   static DataVector wave_profile(const DataVector& x_prime);
     169           0 :   static DataVector charge_density(const DataVector& x_prime);
     170           0 :   double wave_speed_ = std::numeric_limits<double>::signaling_NaN();
     171           0 :   double lorentz_factor_ = std::numeric_limits<double>::signaling_NaN();
     172           0 :   gr::Solutions::Minkowski<3> background_spacetime_{};
     173             : 
     174           0 :   friend bool operator==(const AlfvenWave& lhs, const AlfvenWave& rhs);
     175             : };
     176             : 
     177           0 : bool operator!=(const AlfvenWave& lhs, const AlfvenWave& rhs);
     178             : 
     179             : }  // namespace ForceFree::Solutions