Line data Source code
1 0 : // Distributed under the MIT License.
2 : // See LICENSE.txt for details.
3 :
4 : #pragma once
5 :
6 : #include <cstddef>
7 : #include <limits>
8 : #include <string>
9 :
10 : #include "Evolution/DgSubcell/ReconstructionMethod.hpp"
11 : #include "NumericalAlgorithms/FiniteDifference/DerivativeOrder.hpp"
12 : #include "Options/String.hpp"
13 : #include "Utilities/Literals.hpp"
14 : #include "Utilities/TMPL.hpp"
15 :
16 : /// \cond
17 : namespace PUP {
18 : class er;
19 : } // namespace PUP
20 : /// \endcond
21 :
22 : namespace evolution::dg::subcell {
23 : /*!
24 : * \brief Holds the system-agnostic subcell parameters, such as numbers
25 : * controlling when to switch between DG and subcell.
26 : */
27 1 : class SubcellOptions {
28 : public:
29 : /// Parameters related to the troubled cell indicator (TCI) that determines
30 : /// when to switch between DG and FD.
31 1 : struct TroubledCellIndicator {
32 0 : static constexpr Options::String help =
33 : "Parameters related to the troubled cell indicator (TCI) that "
34 : "determines when to switch between DG and FD.";
35 : };
36 :
37 0 : struct PerssonTci {
38 0 : static constexpr Options::String help =
39 : "Parameters related to the Persson TCI";
40 0 : using group = TroubledCellIndicator;
41 : };
42 : /// The exponent \f$\alpha\f$ passed to the Persson troubled-cell indicator
43 1 : struct PerssonExponent {
44 0 : static std::string name() { return "Exponent"; }
45 0 : static constexpr Options::String help{
46 : "The exponent at which the error should decrease with (N+1-M)"};
47 0 : using type = double;
48 0 : static constexpr type lower_bound() { return 1.0; }
49 0 : static constexpr type upper_bound() { return 10.0; }
50 0 : using group = PerssonTci;
51 : };
52 : /// The number of highest modes the Persson troubled-cell indicator monitors
53 1 : struct PerssonNumHighestModes {
54 0 : static std::string name() { return "NumHighestModes"; }
55 0 : static constexpr Options::String help{
56 : "The number of highest modes M the Persson TCI monitors."};
57 0 : using type = size_t;
58 0 : static constexpr type lower_bound() { return 1_st; }
59 0 : static constexpr type upper_bound() { return 10_st; }
60 0 : using group = PerssonTci;
61 : };
62 :
63 0 : struct RdmpTci {
64 0 : static constexpr Options::String help =
65 : "Parameters related to the relaxed discrete maximum principle TCI";
66 0 : using group = TroubledCellIndicator;
67 : };
68 : /// The \f$\delta_0\f$ parameter in the relaxed discrete maximum principle
69 : /// troubled-cell indicator
70 1 : struct RdmpDelta0 {
71 0 : static std::string name() { return "Delta0"; }
72 0 : static constexpr Options::String help{"Absolute jump tolerance parameter."};
73 0 : using type = double;
74 0 : static type lower_bound() { return 0.0; }
75 0 : using group = RdmpTci;
76 : };
77 : /// The \f$\epsilon\f$ parameter in the relaxed discrete maximum principle
78 : /// troubled-cell indicator
79 1 : struct RdmpEpsilon {
80 0 : static std::string name() { return "Epsilon"; }
81 0 : static constexpr Options::String help{
82 : "The jump-dependent relaxation constant."};
83 0 : using type = double;
84 0 : static type lower_bound() { return 0.0; }
85 0 : static type upper_bound() { return 1.0; }
86 0 : using group = RdmpTci;
87 : };
88 : /// If true, then we always use the subcell method, not DG.
89 1 : struct AlwaysUseSubcells {
90 0 : static constexpr Options::String help{
91 : "If true, then always use the subcell method (e.g. finite-difference) "
92 : "instead of DG."};
93 0 : using type = bool;
94 0 : using group = TroubledCellIndicator;
95 : };
96 : /// If true, then we allow extension directions to be used for
97 : /// computing the ghost points in the FD subcell method to avoid
98 : /// extrapolation.
99 1 : struct EnableExtensionDirections {
100 0 : static constexpr Options::String help{
101 : "If true, then we allow extension directions to be used for "
102 : "computing the ghost points in the FD subcell method to avoid "
103 : "extrapolation."};
104 0 : using type = bool;
105 0 : using group = TroubledCellIndicator;
106 : };
107 : /// Method to use for reconstructing the DG solution from the subcell
108 : /// solution.
109 1 : struct SubcellToDgReconstructionMethod {
110 0 : static constexpr Options::String help{
111 : "Method to use for reconstructing the DG solution from the subcell "
112 : "solution."};
113 0 : using type = fd::ReconstructionMethod;
114 : };
115 : /// \brief Use a width-one halo of FD elements around any troubled element.
116 : ///
117 : /// This provides a buffer of FD subcells so that as a discontinuity moves
118 : /// from one element to another we do not get any Gibbs phenomenon. In the
119 : /// case where we evolve the spacetime metric (e.g. GH+GRMHD) a halo region
120 : /// provides a buffer in case the stellar surface is near an element
121 : /// boundary. Since the GH variables are interpolated using high-order
122 : /// unlimited reconstruction, they can run into issues with Gibbs phenomenon.
123 1 : struct UseHalo {
124 0 : using type = bool;
125 0 : static constexpr Options::String help = {
126 : "Use a width-one halo of FD elements around any troubled element."
127 : "\n"
128 : "This provides a buffer of FD subcells so that as a discontinuity "
129 : "moves from one element to another we do not get any Gibbs "
130 : "phenomenon."};
131 0 : using group = TroubledCellIndicator;
132 : };
133 :
134 : /// \brief The order of the FD derivative used.
135 : ///
136 : /// Must be one of 2, 4, 6, 8, or 10. If `Auto` then the derivative order is
137 : /// determined for you, typically the next-lowest even order compared with
138 : /// the reconstruction scheme order. E.g. for a 5th-order reconstruction we
139 : /// would use 4th order derivatives.
140 1 : struct FiniteDifferenceDerivativeOrder {
141 0 : using type = ::fd::DerivativeOrder;
142 0 : static constexpr Options::String help = {
143 : "The finite difference derivative order to use. If computed from the "
144 : "reconstruction, then the reconstruction method must support returning "
145 : "its reconstruction order."};
146 : };
147 :
148 0 : struct FdToDgTci {
149 0 : static constexpr Options::String help =
150 : "Options related to how quickly we switch from FD to DG.";
151 0 : using group = TroubledCellIndicator;
152 : };
153 : /// The number of time steps taken between calls to the TCI to check if we
154 : /// can go back to the DG grid. A value of `1` means every time step, while
155 : /// `2` means every other time step.
156 1 : struct NumberOfStepsBetweenTciCalls {
157 0 : static constexpr Options::String help{
158 : "The number of time steps taken between calls to the TCI to check if "
159 : "we can go back to the DG grid. A value of `1` means every time step, "
160 : "while `2` means every other time step."};
161 0 : using type = size_t;
162 0 : static constexpr type lower_bound() { return 1; }
163 0 : using group = FdToDgTci;
164 : };
165 : /// The number of time steps/TCI calls after a switch from DG to FD before we
166 : /// allow switching back to DG.
167 1 : struct MinTciCallsAfterRollback {
168 0 : static constexpr Options::String help{
169 : "The number of time steps/TCI calls after a switch from DG to FD "
170 : "before we allow switching back to DG."};
171 0 : using type = size_t;
172 0 : static constexpr type lower_bound() { return 1; }
173 0 : using group = FdToDgTci;
174 : };
175 : /// The number of time steps/TCI calls that the TCI needs to have decided
176 : /// switching to DG is fine before we actually do the switch.
177 1 : struct MinimumClearTcis {
178 0 : static constexpr Options::String help{
179 : "The number of time steps/TCI calls that the TCI needs to have decided "
180 : "switching to DG is fine before we actually do the switch."};
181 0 : using type = size_t;
182 0 : static constexpr type lower_bound() { return 1; }
183 0 : using group = FdToDgTci;
184 : };
185 :
186 : /// \brief FD to FD interpolation order.
187 : ///
188 1 : struct FdInterpolationOrder {
189 0 : using type = size_t;
190 0 : static constexpr Options::String help = {
191 : "The interpolation order to use when interpolating from FD block to FD "
192 : "block."};
193 0 : static constexpr type lower_bound() { return 1; }
194 0 : static constexpr type upper_bound() { return 3; }
195 : };
196 :
197 0 : using options =
198 : tmpl::list<PerssonExponent, PerssonNumHighestModes, RdmpDelta0,
199 : RdmpEpsilon, AlwaysUseSubcells, EnableExtensionDirections,
200 : SubcellToDgReconstructionMethod, UseHalo,
201 : FiniteDifferenceDerivativeOrder, NumberOfStepsBetweenTciCalls,
202 : MinTciCallsAfterRollback, MinimumClearTcis,
203 : FdInterpolationOrder>;
204 :
205 0 : static constexpr Options::String help{
206 : "System-agnostic options for the DG-subcell method."};
207 :
208 0 : SubcellOptions() = default;
209 0 : SubcellOptions(double persson_exponent, size_t persson_num_highest_modes,
210 : double rdmp_delta0, double rdmp_epsilon,
211 : bool always_use_subcells, bool enable_extension_directions,
212 : fd::ReconstructionMethod recons_method, bool use_halo,
213 : ::fd::DerivativeOrder finite_difference_derivative_order,
214 : size_t number_of_steps_between_tci_calls,
215 : size_t min_tci_calls_after_rollback,
216 : size_t min_clear_tci_before_dg,
217 : size_t fd_to_fd_interp_order = 1_st);
218 :
219 0 : void pup(PUP::er& p);
220 :
221 0 : double persson_exponent() const { return persson_exponent_; }
222 :
223 0 : size_t persson_num_highest_modes() const {
224 : return persson_num_highest_modes_;
225 : }
226 :
227 0 : double rdmp_delta0() const { return rdmp_delta0_; }
228 :
229 0 : double rdmp_epsilon() const { return rdmp_epsilon_; }
230 :
231 0 : bool always_use_subcells() const { return always_use_subcells_; }
232 :
233 0 : bool enable_extension_directions() const {
234 : return enable_extension_directions_;
235 : }
236 :
237 0 : fd::ReconstructionMethod reconstruction_method() const {
238 : return reconstruction_method_;
239 : }
240 :
241 0 : bool use_halo() const { return use_halo_; }
242 :
243 0 : size_t get_fd_to_fd_interp_order() const { return fd_to_fd_interp_order_; }
244 :
245 0 : ::fd::DerivativeOrder finite_difference_derivative_order() const {
246 : return finite_difference_derivative_order_;
247 : }
248 :
249 : /// The number of time steps between when we check the TCI to see if we can
250 : /// switch from FD back to DG.
251 1 : size_t number_of_steps_between_tci_calls() const {
252 : return number_of_steps_between_tci_calls_;
253 : }
254 :
255 : /// The number of time steps after a rollback before we check the TCI to see
256 : /// if we can switch from FD back to DG.
257 1 : size_t min_tci_calls_after_rollback() const {
258 : return min_tci_calls_after_rollback_;
259 : }
260 :
261 : /// The number of times the TCI must have flagged that we can switch back to
262 : /// DG before doing so.
263 : ///
264 : /// Note that if we only check the TCI every
265 : /// `number_of_steps_between_tci_calls()` then it takes
266 : /// `number_of_steps_between_tci_calls * min_clear_tci_before_dg` time steps
267 : /// before switching from FD back to DG.
268 : ///
269 : /// `0 means
270 1 : size_t min_clear_tci_before_dg() const { return min_clear_tci_before_dg_; }
271 :
272 : private:
273 0 : friend bool operator==(const SubcellOptions& lhs, const SubcellOptions& rhs);
274 :
275 0 : double persson_exponent_ = std::numeric_limits<double>::signaling_NaN();
276 0 : size_t persson_num_highest_modes_ =
277 : std::numeric_limits<size_t>::signaling_NaN();
278 0 : double rdmp_delta0_ = std::numeric_limits<double>::signaling_NaN();
279 0 : double rdmp_epsilon_ = std::numeric_limits<double>::signaling_NaN();
280 0 : bool always_use_subcells_ = false;
281 0 : bool enable_extension_directions_ = false;
282 0 : fd::ReconstructionMethod reconstruction_method_ =
283 : fd::ReconstructionMethod::AllDimsAtOnce;
284 0 : bool use_halo_{false};
285 0 : ::fd::DerivativeOrder finite_difference_derivative_order_{};
286 0 : size_t number_of_steps_between_tci_calls_{1};
287 0 : size_t min_tci_calls_after_rollback_{1};
288 0 : size_t min_clear_tci_before_dg_{0};
289 0 : size_t fd_to_fd_interp_order_{};
290 : };
291 :
292 0 : bool operator!=(const SubcellOptions& lhs, const SubcellOptions& rhs);
293 : } // namespace evolution::dg::subcell
|
