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 <optional>
11 : #include <pup.h>
12 : #include <utility>
13 :
14 : #include "DataStructures/DataBox/Prefixes.hpp"
15 : #include "DataStructures/Tensor/TypeAliases.hpp"
16 : #include "DataStructures/VariablesTag.hpp"
17 : #include "Domain/Structure/DirectionalIdMap.hpp"
18 : #include "Domain/Tags.hpp"
19 : #include "Evolution/DgSubcell/Tags/GhostDataForReconstruction.hpp"
20 : #include "Evolution/DgSubcell/Tags/Mesh.hpp"
21 : #include "Evolution/DiscontinuousGalerkin/Actions/NormalCovectorAndMagnitude.hpp"
22 : #include "Evolution/Systems/ForceFree/FiniteDifference/Reconstructor.hpp"
23 : #include "Evolution/Systems/ForceFree/FiniteDifference/Tags.hpp"
24 : #include "Evolution/Systems/ForceFree/Tags.hpp"
25 : #include "NumericalAlgorithms/FiniteDifference/FallbackReconstructorType.hpp"
26 : #include "Options/Auto.hpp"
27 : #include "Options/Context.hpp"
28 : #include "Options/String.hpp"
29 : #include "PointwiseFunctions/GeneralRelativity/TagsDeclarations.hpp"
30 : #include "Utilities/TMPL.hpp"
31 :
32 : /// \cond
33 : class DataVector;
34 : template <size_t dim>
35 : class Direction;
36 : template <size_t dim>
37 : class Element;
38 : template <size_t dim>
39 : class ElementId;
40 : template <size_t dim>
41 : class Mesh;
42 : template <typename recons_tags>
43 : class Variables;
44 : namespace gsl {
45 : template <typename>
46 : class not_null;
47 : } // namespace gsl
48 : namespace PUP {
49 : class er;
50 : } // namespace PUP
51 : /// \endcond
52 :
53 1 : namespace ForceFree::fd {
54 :
55 : /*!
56 : * \brief Adaptive order FD reconstruction. See
57 : * ::fd::reconstruction::positivity_preserving_adaptive_order() for details.
58 : * Note that in the ForceFree evolution system no variable needs to be kept
59 : * positive.
60 : *
61 : */
62 1 : class AdaptiveOrder : public Reconstructor {
63 : private:
64 0 : using TildeE = ForceFree::Tags::TildeE;
65 0 : using TildeB = ForceFree::Tags::TildeB;
66 0 : using TildePsi = ForceFree::Tags::TildePsi;
67 0 : using TildePhi = ForceFree::Tags::TildePhi;
68 0 : using TildeQ = ForceFree::Tags::TildeQ;
69 0 : using TildeJ = ForceFree::Tags::TildeJ;
70 :
71 0 : using volume_vars_tags =
72 : tmpl::list<TildeE, TildeB, TildePsi, TildePhi, TildeQ>;
73 :
74 0 : using recons_tags = ForceFree::fd::tags_list_for_reconstruction;
75 :
76 0 : using FallbackReconstructorType =
77 : ::fd::reconstruction::FallbackReconstructorType;
78 :
79 : public:
80 0 : static constexpr size_t dim = 3;
81 :
82 0 : struct Alpha5 {
83 0 : using type = double;
84 0 : static constexpr Options::String help = {
85 : "The alpha parameter in the Persson convergence measurement. 4 is the "
86 : "right value, but anything in the range of 3-5 is 'reasonable'. "
87 : "Smaller values allow for more oscillations."};
88 : };
89 0 : struct Alpha7 {
90 0 : using type = Options::Auto<double, Options::AutoLabel::None>;
91 0 : static constexpr Options::String help = {
92 : "The alpha parameter in the Persson convergence measurement. 4 is the "
93 : "right value, but anything in the range of 3-5 is 'reasonable'. "
94 : "Smaller values allow for more oscillations. If not specified then "
95 : "7th-order reconstruction is not used."};
96 : };
97 0 : struct Alpha9 {
98 0 : using type = Options::Auto<double, Options::AutoLabel::None>;
99 0 : static constexpr Options::String help = {
100 : "The alpha parameter in the Persson convergence measurement. 4 is the "
101 : "right value, but anything in the range of 3-5 is 'reasonable'. "
102 : "Smaller values allow for more oscillations. If not specified then "
103 : "9th-order reconstruction is not used."};
104 : };
105 0 : struct LowOrderReconstructor {
106 0 : using type = FallbackReconstructorType;
107 0 : static constexpr Options::String help = {
108 : "The 2nd/3rd-order reconstruction scheme to use if unlimited 5th-order "
109 : "isn't okay."};
110 : };
111 :
112 0 : using options = tmpl::list<Alpha5, Alpha7, Alpha9, LowOrderReconstructor>;
113 :
114 0 : static constexpr Options::String help{"Adaptive-order reconstruction."};
115 0 : AdaptiveOrder() = default;
116 0 : AdaptiveOrder(AdaptiveOrder&&) = default;
117 0 : AdaptiveOrder& operator=(AdaptiveOrder&&) = default;
118 0 : AdaptiveOrder(const AdaptiveOrder&) = default;
119 0 : AdaptiveOrder& operator=(const AdaptiveOrder&) = default;
120 0 : ~AdaptiveOrder() override = default;
121 :
122 0 : AdaptiveOrder(double alpha_5, std::optional<double> alpha_7,
123 : std::optional<double> alpha_9,
124 : FallbackReconstructorType low_order_reconstructor,
125 : const Options::Context& context = {});
126 :
127 0 : explicit AdaptiveOrder(CkMigrateMessage* msg);
128 :
129 0 : WRAPPED_PUPable_decl_base_template(Reconstructor, AdaptiveOrder);
130 :
131 0 : auto get_clone() const -> std::unique_ptr<Reconstructor> override;
132 :
133 0 : static constexpr bool use_adaptive_order = true;
134 0 : bool supports_adaptive_order() const override { return use_adaptive_order; }
135 :
136 0 : void pup(PUP::er& p) override;
137 :
138 0 : size_t ghost_zone_size() const override {
139 : return eight_to_the_alpha_9_.has_value()
140 : ? 5
141 : : (six_to_the_alpha_7_.has_value() ? 4 : 3);
142 : }
143 :
144 0 : using reconstruction_argument_tags =
145 : tmpl::list<::Tags::Variables<volume_vars_tags>, TildeJ,
146 : domain::Tags::Element<dim>,
147 : evolution::dg::subcell::Tags::GhostDataForReconstruction<dim>,
148 : evolution::dg::subcell::Tags::Mesh<dim>>;
149 :
150 0 : void reconstruct(
151 : gsl::not_null<std::array<Variables<recons_tags>, dim>*>
152 : vars_on_lower_face,
153 : gsl::not_null<std::array<Variables<recons_tags>, dim>*>
154 : vars_on_upper_face,
155 : const Variables<volume_vars_tags>& volume_vars,
156 : const tnsr::I<DataVector, 3, Frame::Inertial>& tilde_j,
157 : const Element<dim>& element,
158 : const DirectionalIdMap<dim, evolution::dg::subcell::GhostData>&
159 : ghost_data,
160 : const Mesh<dim>& subcell_mesh) const;
161 :
162 0 : void reconstruct_fd_neighbor(
163 : gsl::not_null<Variables<recons_tags>*> vars_on_face,
164 : const Variables<volume_vars_tags>& volume_vars,
165 : const tnsr::I<DataVector, 3, Frame::Inertial>& tilde_j,
166 : const Element<dim>& element,
167 : const DirectionalIdMap<dim, evolution::dg::subcell::GhostData>&
168 : ghost_data,
169 : const Mesh<dim>& subcell_mesh,
170 : const Direction<dim> direction_to_reconstruct) const;
171 :
172 : private:
173 : // NOLINTNEXTLINE(readability-redundant-declaration)
174 0 : friend bool operator==(const AdaptiveOrder& lhs, const AdaptiveOrder& rhs);
175 0 : friend bool operator!=(const AdaptiveOrder& lhs, const AdaptiveOrder& rhs);
176 0 : void set_function_pointers();
177 :
178 0 : double four_to_the_alpha_5_ = std::numeric_limits<double>::signaling_NaN();
179 0 : std::optional<double> six_to_the_alpha_7_{};
180 0 : std::optional<double> eight_to_the_alpha_9_{};
181 0 : FallbackReconstructorType low_order_reconstructor_ =
182 : FallbackReconstructorType::None;
183 :
184 0 : using PointerReconsOrder = void (*)(
185 : gsl::not_null<std::array<gsl::span<double>, dim>*>,
186 : gsl::not_null<std::array<gsl::span<double>, dim>*>,
187 : gsl::not_null<std::optional<std::array<gsl::span<std::uint8_t>, dim>>*>,
188 : const gsl::span<const double>&,
189 : const DirectionMap<dim, gsl::span<const double>>&, const Index<dim>&,
190 : size_t, double, double, double);
191 0 : using PointerRecons =
192 : void (*)(gsl::not_null<std::array<gsl::span<double>, dim>*>,
193 : gsl::not_null<std::array<gsl::span<double>, dim>*>,
194 : const gsl::span<const double>&,
195 : const DirectionMap<dim, gsl::span<const double>>&,
196 : const Index<dim>&, size_t, double, double, double);
197 0 : PointerRecons reconstruct_ = nullptr;
198 0 : PointerReconsOrder pp_reconstruct_ = nullptr;
199 :
200 0 : using PointerNeighbor = void (*)(gsl::not_null<DataVector*>,
201 : const DataVector&, const DataVector&,
202 : const Index<dim>&, const Index<dim>&,
203 : const Direction<dim>&, const double&,
204 : const double&, const double&);
205 0 : PointerNeighbor reconstruct_lower_neighbor_ = nullptr;
206 0 : PointerNeighbor reconstruct_upper_neighbor_ = nullptr;
207 0 : PointerNeighbor pp_reconstruct_lower_neighbor_ = nullptr;
208 0 : PointerNeighbor pp_reconstruct_upper_neighbor_ = nullptr;
209 : };
210 :
211 : } // namespace ForceFree::fd
|
