Line data Source code
1 1 : // Distributed under the MIT License.
2 : // See LICENSE.txt for details.
3 :
4 : /// \file
5 : /// Defines MathFunctions::Gaussian.
6 :
7 : #pragma once
8 :
9 : #include <array>
10 : #include <memory>
11 : #include <pup.h>
12 :
13 : #include "Options/String.hpp"
14 : #include "PointwiseFunctions/MathFunctions/MathFunction.hpp"
15 : #include "Utilities/Serialization/CharmPupable.hpp"
16 : #include "Utilities/TMPL.hpp"
17 :
18 : /// \cond
19 : class DataVector;
20 : /// \endcond
21 :
22 : namespace MathFunctions {
23 : template <size_t VolumeDim, typename Fr>
24 : class Gaussian;
25 :
26 : /*!
27 : * \ingroup MathFunctionsGroup
28 : * \brief 1D Gaussian \f$f = A \exp\left(-\frac{(x-x_0)^2}{w^2}\right)\f$
29 : *
30 : * \details Input file options are: Amplitude, Width, and Center. The function
31 : * takes input of type `double` or `DataVector` and returns
32 : * the same type as the input type.
33 : */
34 : template <typename Fr>
35 1 : class Gaussian<1, Fr> : public MathFunction<1, Fr> {
36 : public:
37 0 : struct Amplitude {
38 0 : using type = double;
39 0 : static constexpr Options::String help = {"The amplitude."};
40 : };
41 :
42 0 : struct Width {
43 0 : using type = double;
44 0 : static constexpr Options::String help = {"The width."};
45 0 : static type lower_bound() { return 0.; }
46 : };
47 :
48 0 : struct Center {
49 0 : using type = double;
50 0 : static constexpr Options::String help = {"The center."};
51 : };
52 0 : using options = tmpl::list<Amplitude, Width, Center>;
53 :
54 0 : static constexpr Options::String help = {
55 : "Computes a Gaussian about an arbitrary coordinate center with given "
56 : "width and amplitude"};
57 :
58 0 : WRAPPED_PUPable_decl_base_template(SINGLE_ARG(MathFunction<1, Fr>),
59 : Gaussian); // NOLINT
60 :
61 0 : explicit Gaussian(CkMigrateMessage* /*unused*/) {}
62 :
63 0 : Gaussian(double amplitude, double width, double center);
64 0 : Gaussian(double amplitude, double width, const std::array<double, 1>& center);
65 :
66 0 : Gaussian() = default;
67 0 : std::unique_ptr<MathFunction<1, Fr>> get_clone() const override;
68 :
69 1 : double operator()(const double& x) const override;
70 0 : DataVector operator()(const DataVector& x) const override;
71 : using MathFunction<1, Fr>::operator();
72 :
73 1 : double first_deriv(const double& x) const override;
74 0 : DataVector first_deriv(const DataVector& x) const override;
75 : using MathFunction<1, Fr>::first_deriv;
76 :
77 1 : double second_deriv(const double& x) const override;
78 0 : DataVector second_deriv(const DataVector& x) const override;
79 : using MathFunction<1, Fr>::second_deriv;
80 :
81 1 : double third_deriv(const double& x) const override;
82 0 : DataVector third_deriv(const DataVector& x) const override;
83 : using MathFunction<1, Fr>::third_deriv;
84 :
85 0 : bool operator==(const MathFunction<1, Fr>& other) const override;
86 0 : bool operator!=(const MathFunction<1, Fr>& other) const override;
87 :
88 0 : bool operator==(const Gaussian<1, Fr>& other) const;
89 0 : bool operator!=(const Gaussian<1, Fr>& other) const;
90 : // NOLINTNEXTLINE(google-runtime-references)
91 0 : void pup(PUP::er& p) override;
92 :
93 : private:
94 0 : double amplitude_{};
95 0 : double inverse_width_{};
96 0 : double center_{};
97 :
98 : template <typename T>
99 0 : T apply_call_operator(const T& x) const;
100 : template <typename T>
101 0 : T apply_first_deriv(const T& x) const;
102 : template <typename T>
103 0 : T apply_second_deriv(const T& x) const;
104 : template <typename T>
105 0 : T apply_third_deriv(const T& x) const;
106 : };
107 :
108 : /*!
109 : * \ingroup MathFunctionsGroup
110 : * \brief Gaussian \f$f = A \exp\left(-\frac{(x-x_0)^2}{w^2}\right)\f$
111 : *
112 : * \details Input file options are: Amplitude, Width, and Center. The function
113 : * takes input coordinates of type `tnsr::I<T, VolumeDim, Fr>`, where `T` is
114 : * e.g. `double` or `DataVector`, `Fr` is a frame (e.g. `Frame::Inertial`), and
115 : * `VolumeDim` is the dimension of the spatial volume, i.e. 2 or 3. (The case of
116 : * VolumeDim == 1 is handled specially by Gaussian<1, T, Frame::Inertial>.)
117 : */
118 : template <size_t VolumeDim, typename Fr>
119 1 : class Gaussian : public MathFunction<VolumeDim, Fr> {
120 : public:
121 0 : struct Amplitude {
122 0 : using type = double;
123 0 : static constexpr Options::String help = {"The amplitude."};
124 : };
125 :
126 0 : struct Width {
127 0 : using type = double;
128 0 : static constexpr Options::String help = {"The width."};
129 0 : static type lower_bound() { return 0.; }
130 : };
131 :
132 0 : struct Center {
133 0 : using type = std::array<double, VolumeDim>;
134 0 : static constexpr Options::String help = {"The center."};
135 : };
136 0 : using options = tmpl::list<Amplitude, Width, Center>;
137 :
138 0 : static constexpr Options::String help = {
139 : "Computes a Gaussian about an arbitrary coordinate center with given "
140 : "width and amplitude"};
141 :
142 0 : WRAPPED_PUPable_decl_base_template(SINGLE_ARG(MathFunction<VolumeDim, Fr>),
143 : Gaussian); // NOLINT
144 :
145 0 : explicit Gaussian(CkMigrateMessage* /*unused*/) {}
146 :
147 0 : Gaussian(double amplitude, double width,
148 : const std::array<double, VolumeDim>& center);
149 0 : Gaussian() = default;
150 0 : std::unique_ptr<MathFunction<VolumeDim, Fr>> get_clone() const override;
151 :
152 1 : Scalar<double> operator()(
153 : const tnsr::I<double, VolumeDim, Fr>& x) const override;
154 1 : Scalar<DataVector> operator()(
155 : const tnsr::I<DataVector, VolumeDim, Fr>& x) const override;
156 :
157 1 : tnsr::i<double, VolumeDim, Fr> first_deriv(
158 : const tnsr::I<double, VolumeDim, Fr>& x) const override;
159 1 : tnsr::i<DataVector, VolumeDim, Fr> first_deriv(
160 : const tnsr::I<DataVector, VolumeDim, Fr>& x) const override;
161 :
162 1 : tnsr::ii<double, VolumeDim, Fr> second_deriv(
163 : const tnsr::I<double, VolumeDim, Fr>& x) const override;
164 1 : tnsr::ii<DataVector, VolumeDim, Fr> second_deriv(
165 : const tnsr::I<DataVector, VolumeDim, Fr>& x) const override;
166 :
167 1 : tnsr::iii<double, VolumeDim, Fr> third_deriv(
168 : const tnsr::I<double, VolumeDim, Fr>& x) const override;
169 1 : tnsr::iii<DataVector, VolumeDim, Fr> third_deriv(
170 : const tnsr::I<DataVector, VolumeDim, Fr>& x) const override;
171 :
172 0 : bool operator==(const MathFunction<VolumeDim, Fr>& other) const override;
173 0 : bool operator!=(const MathFunction<VolumeDim, Fr>& other) const override;
174 :
175 0 : bool operator==(const Gaussian<VolumeDim, Fr>& other) const;
176 0 : bool operator!=(const Gaussian<VolumeDim, Fr>& other) const;
177 : // NOLINTNEXTLINE(google-runtime-references)
178 0 : void pup(PUP::er& p) override;
179 :
180 : private:
181 0 : double amplitude_{};
182 0 : double inverse_width_{};
183 0 : std::array<double, VolumeDim> center_{};
184 :
185 : template <typename T>
186 0 : tnsr::I<T, VolumeDim, Fr> centered_coordinates(
187 : const tnsr::I<T, VolumeDim, Fr>& x) const;
188 : template <typename T>
189 0 : Scalar<T> apply_call_operator(
190 : const tnsr::I<T, VolumeDim, Fr>& centered_coords) const;
191 : template <typename T>
192 0 : tnsr::i<T, VolumeDim, Fr> apply_first_deriv(
193 : const tnsr::I<T, VolumeDim, Fr>& centered_coords,
194 : const Scalar<T>& gaussian) const;
195 : template <typename T>
196 0 : tnsr::ii<T, VolumeDim, Fr> apply_second_deriv(
197 : const tnsr::I<T, VolumeDim, Fr>& centered_coords,
198 : const Scalar<T>& gaussian,
199 : const tnsr::i<T, VolumeDim, Fr>& d_gaussian) const;
200 : template <typename T>
201 0 : tnsr::iii<T, VolumeDim, Fr> apply_third_deriv(
202 : const tnsr::I<T, VolumeDim, Fr>& centered_coords,
203 : const Scalar<T>& gaussian, const tnsr::i<T, VolumeDim, Fr>& d_gaussian,
204 : const tnsr::ii<T, VolumeDim, Fr>& d2_gaussian) const;
205 : };
206 :
207 : template <size_t VolumeDim, typename Fr>
208 0 : bool operator!=(const Gaussian<VolumeDim, Fr>& lhs,
209 : const Gaussian<VolumeDim, Fr>& rhs) {
210 : return not(lhs == rhs);
211 : }
212 : } // namespace MathFunctions
213 :
214 : /// \cond
215 : template <size_t VolumeDim, typename Fr>
216 : PUP::able::PUP_ID MathFunctions::Gaussian<VolumeDim, Fr>::my_PUP_ID =
217 : 0; // NOLINT
218 :
219 : template <typename Fr>
220 : PUP::able::PUP_ID MathFunctions::Gaussian<1, Fr>::my_PUP_ID = 0; // NOLINT
221 : /// \endcond
|
