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 <cstdint>
8 : #include <optional>
9 : #include <pup.h>
10 : #include <type_traits>
11 :
12 : #include "DataStructures/MathWrapper.hpp"
13 : #include "DataStructures/TaggedVariant.hpp"
14 : #include "Time/History.hpp"
15 : #include "Time/StepperErrorEstimate.hpp"
16 : #include "Utilities/GenerateInstantiations.hpp"
17 : #include "Utilities/Gsl.hpp"
18 : #include "Utilities/Serialization/CharmPupable.hpp"
19 :
20 : /// \cond
21 : struct StepperErrorTolerances;
22 : class TimeDelta;
23 : class TimeStepId;
24 : /// \endcond
25 :
26 : /// \ingroup TimeSteppersGroup
27 : ///
28 : /// Holds classes that take time steps.
29 : namespace TimeSteppers {}
30 :
31 : /// \cond
32 : #define TIME_STEPPER_WRAPPED_TYPE(data) BOOST_PP_TUPLE_ELEM(0, data)
33 : #define TIME_STEPPER_DERIVED_CLASS(data) BOOST_PP_TUPLE_ELEM(1, data)
34 : #define TIME_STEPPER_DERIVED_CLASS_TEMPLATE(data) BOOST_PP_TUPLE_ELEM(2, data)
35 : /// \endcond
36 :
37 : namespace TimeSteppers {
38 : /// Minimum and maximum orders of a variable-order TimeStepper.
39 1 : struct VariableOrder {
40 0 : size_t minimum;
41 0 : size_t maximum;
42 :
43 0 : VariableOrder(const size_t minimum_in, const size_t maximum_in)
44 : : minimum(minimum_in), maximum(maximum_in) {}
45 : };
46 :
47 0 : bool operator==(const VariableOrder& a, const VariableOrder& b);
48 0 : bool operator!=(const VariableOrder& a, const VariableOrder& b);
49 :
50 0 : namespace Tags {
51 : /// Order of a fixed-order TimeStepper.
52 1 : struct FixedOrder {
53 0 : using type = size_t;
54 : };
55 :
56 : /// Minimum and maximum orders of a variable-order TimeStepper.
57 1 : struct VariableOrder {
58 0 : using type = TimeSteppers::VariableOrder;
59 : };
60 : } // namespace Tags
61 : } // namespace TimeSteppers
62 :
63 : /// \ingroup TimeSteppersGroup
64 : ///
65 : /// Abstract base class for TimeSteppers.
66 : ///
67 : /// Several of the member functions of this class are templated and
68 : /// perform type erasure before forwarding their arguments to the
69 : /// derived classes. This is implemented using the macros \ref
70 : /// TIME_STEPPER_DECLARE_OVERLOADS, which must be placed in a private
71 : /// section of the class body, and
72 : /// TIME_STEPPER_DEFINE_OVERLOADS(derived_class), which must be placed
73 : /// in the cpp file.
74 1 : class TimeStepper : public PUP::able {
75 : public:
76 0 : static constexpr bool imex = false;
77 0 : using provided_time_stepper_interfaces = tmpl::list<TimeStepper>;
78 :
79 0 : WRAPPED_PUPable_abstract(TimeStepper); // NOLINT
80 :
81 : /// \cond
82 : #define TIME_STEPPER_DECLARE_VIRTUALS_IMPL(_, data) \
83 : virtual void update_u_forward( \
84 : const gsl::not_null<TIME_STEPPER_WRAPPED_TYPE(data)*> u, \
85 : const TimeSteppers::ConstUntypedHistory<TIME_STEPPER_WRAPPED_TYPE( \
86 : data)>& history, \
87 : const TimeDelta& time_step) const = 0; \
88 : virtual std::optional<StepperErrorEstimate> update_u_forward( \
89 : const gsl::not_null<TIME_STEPPER_WRAPPED_TYPE(data)*> u, \
90 : const TimeSteppers::ConstUntypedHistory<TIME_STEPPER_WRAPPED_TYPE( \
91 : data)>& history, \
92 : const TimeDelta& time_step, \
93 : const StepperErrorTolerances& tolerances) const = 0; \
94 : virtual void clean_history_forward( \
95 : const TimeSteppers::MutableUntypedHistory<TIME_STEPPER_WRAPPED_TYPE( \
96 : data)>& history) const = 0; \
97 : virtual bool dense_update_u_forward( \
98 : const gsl::not_null<TIME_STEPPER_WRAPPED_TYPE(data)*> u, \
99 : const TimeSteppers::ConstUntypedHistory<TIME_STEPPER_WRAPPED_TYPE( \
100 : data)>& history, \
101 : double time) const = 0; \
102 : virtual bool can_change_step_size_forward( \
103 : const TimeStepId& time_id, \
104 : const TimeSteppers::ConstUntypedHistory<TIME_STEPPER_WRAPPED_TYPE( \
105 : data)>& history) const = 0;
106 :
107 : GENERATE_INSTANTIATIONS(TIME_STEPPER_DECLARE_VIRTUALS_IMPL,
108 : (MATH_WRAPPER_TYPES))
109 : #undef TIME_STEPPER_DECLARE_VIRTUALS_IMPL
110 : /// \endcond
111 :
112 : /// Set \p u to the value at the end of the current substep.
113 : ///
114 : /// Derived classes must implement this as a function with signature
115 : ///
116 : /// ```
117 : /// template <typename T>
118 : /// void update_u_impl(gsl::not_null<T*> u,
119 : /// const ConstUntypedHistory<T>& history,
120 : /// const TimeDelta& time_step) const;
121 : /// ```
122 : template <typename Vars>
123 1 : void update_u(const gsl::not_null<Vars*> u,
124 : const TimeSteppers::History<Vars>& history,
125 : const TimeDelta& time_step) const {
126 : return update_u_forward(&*make_math_wrapper(u), history.untyped(),
127 : time_step);
128 : }
129 :
130 : /// Set \p u to the value at the end of the current substep; report the error
131 : /// measure when available.
132 : ///
133 : /// For a substep method, the error measure will only be available on full
134 : /// steps. For a multistep method, the error measure will only be available
135 : /// when a sufficient number of steps are available in the `history` to
136 : /// compare two orders of step. Whenever the error measure is unavailable,
137 : /// the return value is empty
138 : ///
139 : /// If \p tolerances requests no estimates (e.g., a
140 : /// default-constructed object), no error measures are calculated,
141 : /// but any additional substeps necessary for error estimation are
142 : /// still taken. This is useful when a system integrates multiple
143 : /// variables in separate calls to the TimeStepper, but only uses
144 : /// estimates from some of them.
145 : ///
146 : /// Derived classes must implement this as a function with signature
147 : ///
148 : /// ```
149 : /// template <typename T>
150 : /// std::optional<StepperErrorEstimate> update_u_impl(
151 : /// gsl::not_null<T*> u,
152 : /// const ConstUntypedHistory<T>& history,
153 : /// const TimeDelta& time_step,
154 : /// const StepperErrorTolerances& tolerances) const;
155 : /// ```
156 : template <typename Vars>
157 1 : std::optional<StepperErrorEstimate> update_u(
158 : const gsl::not_null<Vars*> u, const TimeSteppers::History<Vars>& history,
159 : const TimeDelta& time_step,
160 : const StepperErrorTolerances& tolerances) const {
161 : return update_u_forward(&*make_math_wrapper(u), history.untyped(),
162 : time_step, tolerances);
163 : }
164 :
165 : /// Remove old entries from the history.
166 : ///
167 : /// This should be called after update_u and dense output.
168 : /// Afterward, the history will generally require a new entry to be
169 : /// added before it can be used by the TimeStepper.
170 : ///
171 : /// Derived classes must implement this as a function with signature
172 : ///
173 : /// ```
174 : /// template <typename T>
175 : /// void clean_history_impl(const MutableUntypedHistory<T>& history) const;
176 : /// ```
177 : template <typename Vars>
178 1 : void clean_history(
179 : const gsl::not_null<TimeSteppers::History<Vars>*> history) const {
180 : return clean_history_forward(history->untyped());
181 : }
182 :
183 : /// Compute the solution value at a time between steps. To evaluate
184 : /// at a time within a given step, call this method at the start of
185 : /// the step containing the time. The function returns true on
186 : /// success, otherwise the call should be retried after the next
187 : /// substep.
188 : ///
189 : /// The change from the partial step will be added to the initial
190 : /// value, so \p u should generally be initialized to
191 : /// `*history.step_start(time).value`. (TimeStepper
192 : /// implementations are required to keep this value in the history.)
193 : ///
194 : /// Derived classes must implement this as a function with signature
195 : ///
196 : /// ```
197 : /// template <typename T>
198 : /// bool dense_update_u_impl(
199 : /// gsl::not_null<T*> u, const ConstUntypedHistory<T>& history,
200 : /// double time) const;
201 : /// ```
202 : template <typename Vars>
203 1 : bool dense_update_u(const gsl::not_null<Vars*> u,
204 : const TimeSteppers::History<Vars>& history,
205 : const double time) const {
206 : return dense_update_u_forward(&*make_math_wrapper(u), history.untyped(),
207 : time);
208 : }
209 :
210 : /// The convergence order of the stepper
211 : virtual variants::TaggedVariant<TimeSteppers::Tags::FixedOrder,
212 : TimeSteppers::Tags::VariableOrder>
213 1 : order() const = 0;
214 :
215 : /// Number of substeps in this TimeStepper
216 1 : virtual uint64_t number_of_substeps() const = 0;
217 :
218 : /// Number of substeps in this TimeStepper when providing an error measure for
219 : /// adaptive time-stepping
220 : ///
221 : /// \details Certain substep methods (e.g. embedded RK4(3)) require additional
222 : /// steps when providing an error measure of the integration.
223 1 : virtual uint64_t number_of_substeps_for_error() const = 0;
224 :
225 : /// Number of past time steps needed for multi-step method
226 1 : virtual size_t number_of_past_steps() const = 0;
227 :
228 : /// Rough estimate of the maximum step size this method can take
229 : /// stably as a multiple of the step for Euler's method.
230 1 : virtual double stable_step() const = 0;
231 :
232 : /// Whether computational and temporal orderings of operations
233 : /// match.
234 : ///
235 : /// If this method returns true, then, for two time-stepper
236 : /// operations occurring at different simulation times, the
237 : /// temporally earlier operation will be performed first. These
238 : /// operations include RHS evaluation, dense output, and neighbor
239 : /// communication. In particular, dense output never requires
240 : /// performing a RHS evaluation later than the output time, so
241 : /// control systems measurements cannot cause deadlocks.
242 : ///
243 : /// \warning This guarantee only holds if the time steps themselves
244 : /// are monotonic, which can be violated during initialization.
245 1 : virtual bool monotonic() const = 0;
246 :
247 : /// Whether dense output makes use of the first-same-as-last
248 : /// optimization. Time steppers using this require the history to
249 : /// contain the RHS evaluation for the end of a step before dense
250 : /// output can be performed for that step. Time steppers not using
251 : /// this will delete any substep data after the final substep of
252 : /// each step.
253 1 : virtual bool dense_output_uses_fsal() const = 0;
254 :
255 : /// The TimeStepId after the current substep
256 1 : virtual TimeStepId next_time_id(const TimeStepId& current_id,
257 : const TimeDelta& time_step) const = 0;
258 :
259 : /// The TimeStepId after the current substep when providing an error measure
260 : /// for adaptive time-stepping.
261 : ///
262 : /// Certain substep methods (e.g. embedded RK4(3)) require additional
263 : /// steps when providing an error measure of the integration.
264 1 : virtual TimeStepId next_time_id_for_error(
265 : const TimeStepId& current_id, const TimeDelta& time_step) const = 0;
266 :
267 : /// Whether a change in the step size is allowed before taking
268 : /// a step. Step sizes can never be changed on a substep.
269 : ///
270 : /// Derived classes must implement this as a function with signature
271 : ///
272 : /// ```
273 : /// template <typename T>
274 : /// bool can_change_step_size_impl(
275 : /// const TimeStepId& time_id,
276 : /// const ConstUntypedHistory<T>& history) const;
277 : /// ```
278 : template <typename Vars>
279 1 : bool can_change_step_size(const TimeStepId& time_id,
280 : const TimeSteppers::History<Vars>& history) const {
281 : return can_change_step_size_forward(time_id, history.untyped());
282 : }
283 : };
284 :
285 : /// \cond
286 : #define TIME_STEPPER_DECLARE_OVERLOADS_IMPL(_, data) \
287 : void update_u_forward( \
288 : gsl::not_null<TIME_STEPPER_WRAPPED_TYPE(data)*> u, \
289 : const TimeSteppers::ConstUntypedHistory<TIME_STEPPER_WRAPPED_TYPE( \
290 : data)>& history, \
291 : const TimeDelta& time_step) const override; \
292 : std::optional<StepperErrorEstimate> update_u_forward( \
293 : gsl::not_null<TIME_STEPPER_WRAPPED_TYPE(data)*> u, \
294 : const TimeSteppers::ConstUntypedHistory<TIME_STEPPER_WRAPPED_TYPE( \
295 : data)>& history, \
296 : const TimeDelta& time_step, \
297 : const StepperErrorTolerances& tolerances) const override; \
298 : void clean_history_forward( \
299 : const TimeSteppers::MutableUntypedHistory<TIME_STEPPER_WRAPPED_TYPE( \
300 : data)>& history) const override; \
301 : bool dense_update_u_forward( \
302 : gsl::not_null<TIME_STEPPER_WRAPPED_TYPE(data)*> u, \
303 : const TimeSteppers::ConstUntypedHistory<TIME_STEPPER_WRAPPED_TYPE( \
304 : data)>& history, \
305 : double time) const override; \
306 : bool can_change_step_size_forward( \
307 : const TimeStepId& time_id, \
308 : const TimeSteppers::ConstUntypedHistory<TIME_STEPPER_WRAPPED_TYPE( \
309 : data)>& history) const override;
310 :
311 : #define TIME_STEPPER_DEFINE_OVERLOADS_IMPL(_, data) \
312 : TIME_STEPPER_DERIVED_CLASS_TEMPLATE(data) \
313 : void TIME_STEPPER_DERIVED_CLASS(data)::update_u_forward( \
314 : const gsl::not_null<TIME_STEPPER_WRAPPED_TYPE(data)*> u, \
315 : const TimeSteppers::ConstUntypedHistory<TIME_STEPPER_WRAPPED_TYPE( \
316 : data)>& history, \
317 : const TimeDelta& time_step) const { \
318 : update_u_impl(u, history, time_step); \
319 : } \
320 : TIME_STEPPER_DERIVED_CLASS_TEMPLATE(data) \
321 : std::optional<StepperErrorEstimate> \
322 : TIME_STEPPER_DERIVED_CLASS(data)::update_u_forward( \
323 : const gsl::not_null<TIME_STEPPER_WRAPPED_TYPE(data)*> u, \
324 : const TimeSteppers::ConstUntypedHistory<TIME_STEPPER_WRAPPED_TYPE( \
325 : data)>& history, \
326 : const TimeDelta& time_step, \
327 : const StepperErrorTolerances& tolerances) const { \
328 : return update_u_impl(u, history, time_step, tolerances); \
329 : } \
330 : TIME_STEPPER_DERIVED_CLASS_TEMPLATE(data) \
331 : void TIME_STEPPER_DERIVED_CLASS(data)::clean_history_forward( \
332 : const TimeSteppers::MutableUntypedHistory<TIME_STEPPER_WRAPPED_TYPE( \
333 : data)>& history) const { \
334 : return clean_history_impl(history); \
335 : } \
336 : TIME_STEPPER_DERIVED_CLASS_TEMPLATE(data) \
337 : bool TIME_STEPPER_DERIVED_CLASS(data)::dense_update_u_forward( \
338 : const gsl::not_null<TIME_STEPPER_WRAPPED_TYPE(data)*> u, \
339 : const TimeSteppers::ConstUntypedHistory<TIME_STEPPER_WRAPPED_TYPE( \
340 : data)>& history, \
341 : const double time) const { \
342 : return dense_update_u_impl(u, history, time); \
343 : } \
344 : TIME_STEPPER_DERIVED_CLASS_TEMPLATE(data) \
345 : bool TIME_STEPPER_DERIVED_CLASS(data)::can_change_step_size_forward( \
346 : const TimeStepId& time_id, \
347 : const TimeSteppers::ConstUntypedHistory<TIME_STEPPER_WRAPPED_TYPE( \
348 : data)>& history) const { \
349 : return can_change_step_size_impl(time_id, history); \
350 : }
351 : /// \endcond
352 :
353 : /// \ingroup TimeSteppersGroup
354 : /// Macro declaring overloaded detail methods in classes derived from
355 : /// TimeStepper. Must be placed in a private section of the class
356 : /// body.
357 1 : #define TIME_STEPPER_DECLARE_OVERLOADS \
358 : GENERATE_INSTANTIATIONS(TIME_STEPPER_DECLARE_OVERLOADS_IMPL, \
359 : (MATH_WRAPPER_TYPES))
360 :
361 : /// \ingroup TimeSteppersGroup
362 : /// Macro defining overloaded detail methods in classes derived from
363 : /// TimeStepper. Must be placed in the cpp file for the derived
364 : /// class.
365 : /// @{
366 1 : #define TIME_STEPPER_DEFINE_OVERLOADS(derived_class) \
367 : GENERATE_INSTANTIATIONS(TIME_STEPPER_DEFINE_OVERLOADS_IMPL, \
368 : (MATH_WRAPPER_TYPES), (derived_class), ())
369 1 : #define TIME_STEPPER_DEFINE_OVERLOADS_TEMPLATED(derived_class, template_args) \
370 : GENERATE_INSTANTIATIONS(TIME_STEPPER_DEFINE_OVERLOADS_IMPL, \
371 : (MATH_WRAPPER_TYPES), (derived_class), \
372 : (template <template_args>))
373 : /// @}
|
