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/TaggedTuple.hpp"
9 : #include "DataStructures/Tensor/TypeAliases.hpp"
10 : #include "Evolution/Systems/ForceFree/Tags.hpp"
11 : #include "Options/Context.hpp"
12 : #include "Options/String.hpp"
13 : #include "PointwiseFunctions/AnalyticSolutions/AnalyticSolution.hpp"
14 : #include "PointwiseFunctions/AnalyticSolutions/GeneralRelativity/Minkowski.hpp"
15 : #include "PointwiseFunctions/InitialDataUtilities/InitialData.hpp"
16 : #include "Utilities/Serialization/CharmPupable.hpp"
17 : #include "Utilities/TMPL.hpp"
18 :
19 : /// \cond
20 : class DataVector;
21 :
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
|
