Line data Source code
1 0 : // Distributed under the MIT License.
2 : // See LICENSE.txt for details.
3 :
4 : #pragma once
5 :
6 : #include <array>
7 : #include <cstddef>
8 : #include <limits>
9 : #include <memory>
10 : #include <string>
11 : #include <unordered_map>
12 :
13 : #include "DataStructures/Tensor/TypeAliases.hpp"
14 : #include "Evolution/Systems/GeneralizedHarmonic/ConstraintDamping/DampingFunction.hpp"
15 : #include "Options/String.hpp"
16 : #include "Utilities/Serialization/CharmPupable.hpp"
17 : #include "Utilities/TMPL.hpp"
18 :
19 : /// \cond
20 : class DataVector;
21 : namespace PUP {
22 : class er;
23 : } // namespace PUP
24 : namespace domain::FunctionsOfTime {
25 : class FunctionOfTime;
26 : } // namespace domain::FunctionsOfTime
27 : /// \endcond
28 :
29 : namespace gh::ConstraintDamping {
30 : /*!
31 : * \brief A sum of three Gaussians plus a constant, where the Gaussian widths
32 : * are scaled by a domain::FunctionsOfTime::FunctionOfTime.
33 : *
34 : * \details The function \f$f\f$ is given by
35 : * \f{align}{
36 : * f = C + \sum_{\alpha=1}^3
37 : * A_\alpha \exp\left(-\frac{(x-(x_0)_\alpha)^2}{w_\alpha^2(t)}\right).
38 : * \f}
39 : * Input file options are: `Constant` \f$C\f$, `Amplitude[1-3]`
40 : * \f$A_\alpha\f$, `Width[1-3]` \f$w_\alpha\f$, and `Center[1-3]
41 : * `\f$(x_0)_\alpha\f$. The function takes input
42 : * coordinates \f$x\f$ of type `tnsr::I<T, 3, Frame::Grid>`, where `T` is e.g.
43 : * `double` or `DataVector`; note that this DampingFunction is only defined
44 : * for three spatial dimensions and for the grid frame. The Gaussian widths
45 : * \f$w_\alpha\f$ are scaled by the inverse of the value of a scalar
46 : * domain::FunctionsOfTime::FunctionOfTime \f$f(t)\f$: \f$w_\alpha(t) = w_\alpha
47 : * / f(t)\f$.
48 : *
49 : * The name of the domain::FunctionsOfTime::FunctionOfTime is hardcoded to be
50 : * `Expansion` as this class will typically be used with a domain that has an
51 : * expansion map and potentially an expansion control system, the names of which
52 : * will be `Expansion`.
53 : */
54 1 : class TimeDependentTripleGaussian : public DampingFunction<3, Frame::Grid> {
55 : public:
56 : template <size_t GaussianNumber>
57 0 : struct Gaussian {
58 0 : static constexpr Options::String help = {
59 : "Parameters for one of the Gaussians."};
60 0 : static std::string name() {
61 : return "Gaussian" + std::to_string(GaussianNumber);
62 : };
63 : };
64 0 : struct Constant {
65 0 : using type = double;
66 0 : static constexpr Options::String help = {"The constant."};
67 : };
68 :
69 : template <typename Group>
70 0 : struct Amplitude {
71 0 : using group = Group;
72 0 : using type = double;
73 0 : static constexpr Options::String help = {"The amplitude of the Gaussian."};
74 : };
75 :
76 : template <typename Group>
77 0 : struct Width {
78 0 : using group = Group;
79 0 : using type = double;
80 0 : static constexpr Options::String help = {
81 : "The unscaled width of the Gaussian."};
82 0 : static type lower_bound() { return 0.; }
83 : };
84 :
85 : template <typename Group>
86 0 : struct Center {
87 0 : using group = Group;
88 0 : using type = std::array<double, 3>;
89 0 : static constexpr Options::String help = {"The center of the Gaussian."};
90 : };
91 :
92 0 : using options = tmpl::list<
93 : Constant, Amplitude<Gaussian<1>>, Width<Gaussian<1>>, Center<Gaussian<1>>,
94 : Amplitude<Gaussian<2>>, Width<Gaussian<2>>, Center<Gaussian<2>>,
95 : Amplitude<Gaussian<3>>, Width<Gaussian<3>>, Center<Gaussian<3>>>;
96 :
97 0 : static constexpr Options::String help = {
98 : "Computes a sum of a constant and 3 Gaussians (each with its own "
99 : "amplitude, width, and coordinate center), with the Gaussian widths "
100 : "scaled by the inverse of a FunctionOfTime."};
101 :
102 : /// \cond
103 : WRAPPED_PUPable_decl_base_template(
104 : SINGLE_ARG(DampingFunction<3, Frame::Grid>),
105 : TimeDependentTripleGaussian); // NOLINT
106 : /// \endcond
107 :
108 0 : explicit TimeDependentTripleGaussian(CkMigrateMessage* msg);
109 :
110 0 : TimeDependentTripleGaussian(double constant, double amplitude_1,
111 : double width_1,
112 : const std::array<double, 3>& center_1,
113 : double amplitude_2, double width_2,
114 : const std::array<double, 3>& center_2,
115 : double amplitude_3, double width_3,
116 : const std::array<double, 3>& center_3);
117 :
118 0 : TimeDependentTripleGaussian() = default;
119 0 : ~TimeDependentTripleGaussian() override = default;
120 0 : TimeDependentTripleGaussian(const TimeDependentTripleGaussian& /*rhs*/) =
121 : default;
122 0 : TimeDependentTripleGaussian& operator=(
123 : const TimeDependentTripleGaussian& /*rhs*/) = default;
124 0 : TimeDependentTripleGaussian(TimeDependentTripleGaussian&& /*rhs*/) = default;
125 0 : TimeDependentTripleGaussian& operator=(
126 : TimeDependentTripleGaussian&& /*rhs*/) = default;
127 :
128 0 : void operator()(const gsl::not_null<Scalar<double>*> value_at_x,
129 : const tnsr::I<double, 3, Frame::Grid>& x, double time,
130 : const std::unordered_map<
131 : std::string,
132 : std::unique_ptr<domain::FunctionsOfTime::FunctionOfTime>>&
133 : functions_of_time) const override;
134 0 : void operator()(const gsl::not_null<Scalar<DataVector>*> value_at_x,
135 : const tnsr::I<DataVector, 3, Frame::Grid>& x, double time,
136 : const std::unordered_map<
137 : std::string,
138 : std::unique_ptr<domain::FunctionsOfTime::FunctionOfTime>>&
139 : functions_of_time) const override;
140 :
141 0 : auto get_clone() const
142 : -> std::unique_ptr<DampingFunction<3, Frame::Grid>> override;
143 :
144 : // NOLINTNEXTLINE(google-runtime-references)
145 0 : void pup(PUP::er& p) override;
146 :
147 : private:
148 0 : friend bool operator==(const TimeDependentTripleGaussian& lhs,
149 : const TimeDependentTripleGaussian& rhs) {
150 : return lhs.constant_ == rhs.constant_ and
151 : lhs.amplitude_1_ == rhs.amplitude_1_ and
152 : lhs.inverse_width_1_ == rhs.inverse_width_1_ and
153 : lhs.center_1_ == rhs.center_1_ and
154 : lhs.amplitude_2_ == rhs.amplitude_2_ and
155 : lhs.inverse_width_2_ == rhs.inverse_width_2_ and
156 : lhs.center_2_ == rhs.center_2_ and
157 : lhs.amplitude_3_ == rhs.amplitude_3_ and
158 : lhs.inverse_width_3_ == rhs.inverse_width_3_ and
159 : lhs.center_3_ == rhs.center_3_;
160 : }
161 :
162 0 : double constant_ = std::numeric_limits<double>::signaling_NaN();
163 0 : double amplitude_1_ = std::numeric_limits<double>::signaling_NaN();
164 0 : double inverse_width_1_ = std::numeric_limits<double>::signaling_NaN();
165 0 : std::array<double, 3> center_1_{};
166 0 : double amplitude_2_ = std::numeric_limits<double>::signaling_NaN();
167 0 : double inverse_width_2_ = std::numeric_limits<double>::signaling_NaN();
168 0 : std::array<double, 3> center_2_{};
169 0 : double amplitude_3_ = std::numeric_limits<double>::signaling_NaN();
170 0 : double inverse_width_3_ = std::numeric_limits<double>::signaling_NaN();
171 0 : std::array<double, 3> center_3_{};
172 0 : inline static const std::string function_of_time_for_scaling_{"Expansion"};
173 :
174 : template <typename T>
175 0 : void apply_call_operator(
176 : const gsl::not_null<Scalar<T>*> value_at_x,
177 : const tnsr::I<T, 3, Frame::Grid>& x, double time,
178 : const std::unordered_map<
179 : std::string,
180 : std::unique_ptr<domain::FunctionsOfTime::FunctionOfTime>>&
181 : functions_of_time) const;
182 : };
183 :
184 0 : bool operator!=(const TimeDependentTripleGaussian& lhs,
185 : const TimeDependentTripleGaussian& rhs);
186 : } // namespace gh::ConstraintDamping
|
