Line data Source code
1 0 : // Distributed under the MIT License.
2 : // See LICENSE.txt for details.
3 :
4 : #pragma once
5 :
6 : #include <algorithm>
7 : #include <array>
8 : #include <cstddef>
9 : #include <typeinfo>
10 : #include <utility>
11 :
12 : #include "DataStructures/DataBox/DataBox.hpp"
13 : #include "DataStructures/DataBox/ObservationBox.hpp"
14 : #include "Domain/Amr/Flag.hpp"
15 : #include "Domain/Amr/Helpers.hpp"
16 : #include "Domain/Structure/ElementId.hpp"
17 : #include "Domain/Tags.hpp"
18 : #include "IO/Logging/Tags.hpp"
19 : #include "IO/Logging/Verbosity.hpp"
20 : #include "NumericalAlgorithms/Spectral/IsValidDgMesh.hpp"
21 : #include "Options/String.hpp"
22 : #include "Parallel/GlobalCache.hpp"
23 : #include "Parallel/Printf/Printf.hpp"
24 : #include "ParallelAlgorithms/Amr/Criteria/Criterion.hpp"
25 : #include "ParallelAlgorithms/Amr/Criteria/Tags/Criteria.hpp"
26 : #include "ParallelAlgorithms/Amr/Policies/EnforcePolicies.hpp"
27 : #include "ParallelAlgorithms/Amr/Policies/Tags.hpp"
28 : #include "ParallelAlgorithms/Amr/Projectors/Mesh.hpp"
29 : #include "ParallelAlgorithms/Amr/Tags.hpp"
30 : #include "ParallelAlgorithms/EventsAndTriggers/Event.hpp"
31 : #include "Utilities/Algorithm.hpp"
32 : #include "Utilities/ErrorHandling/Assert.hpp"
33 : #include "Utilities/ErrorHandling/Error.hpp"
34 : #include "Utilities/Gsl.hpp"
35 : #include "Utilities/MakeArray.hpp"
36 : #include "Utilities/Serialization/CharmPupable.hpp"
37 : #include "Utilities/TMPL.hpp"
38 :
39 : /// \cond
40 : namespace PUP {
41 : class er;
42 : } // namespace PUP
43 : namespace Tags {
44 : struct TimeStepId;
45 : } // namespace Tags
46 : namespace evolution::dg::Tags {
47 : template <size_t Dim>
48 : struct MortarNextTemporalId;
49 : } // namespace evolution::dg::Tags
50 : /// \endcond
51 :
52 : namespace amr::Events {
53 : namespace detail {
54 : template <typename Criterion>
55 : struct get_tags {
56 : using type = typename Criterion::compute_tags_for_observation_box;
57 : };
58 :
59 : } // namespace detail
60 :
61 : /// \ingroup AmrGroup
62 : /// \brief Performs p-refinement on the domain
63 : ///
64 : /// \details
65 : /// - Loops over all p-refinement criteria specified in the input file.
66 : /// If no valid p-refinement decision is requested, no change is
67 : /// made to the Element.
68 : /// - Updates the Mesh and all return tags of Metavariables::amr::projectors
69 : ///
70 : /// \warning This does not communicate the new Mesh to its neighbors, nor does
71 : /// it update ::domain::Tags::NeighborMesh
72 1 : class RefineMesh : public Event {
73 : public:
74 : /// \cond
75 : explicit RefineMesh(CkMigrateMessage* m);
76 : using PUP::able::register_constructor;
77 : WRAPPED_PUPable_decl_template(RefineMesh); // NOLINT
78 : /// \endcond
79 :
80 0 : using options = tmpl::list<>;
81 0 : static constexpr Options::String help = {"Perform p-refinement"};
82 :
83 0 : RefineMesh();
84 :
85 0 : using compute_tags_for_observation_box = tmpl::list<>;
86 0 : using return_tags = tmpl::list<::Tags::DataBox>;
87 0 : using argument_tags = tmpl::list<>;
88 :
89 : template <typename DbTags, typename Metavariables, typename Component>
90 0 : void operator()(const gsl::not_null<db::DataBox<DbTags>*> box,
91 : Parallel::GlobalCache<Metavariables>& cache,
92 : const ElementId<Metavariables::volume_dim>& element_id,
93 : const Component* const /*meta*/,
94 : const ObservationValue& /*observation_value*/) const {
95 : constexpr size_t volume_dim = Metavariables::volume_dim;
96 : if (alg::any_of(
97 : db::get<evolution::dg::Tags::MortarNextTemporalId<volume_dim>>(
98 : *box),
99 : [&](const auto& mortar) {
100 : return mortar.second != db::get<::Tags::TimeStepId>(*box);
101 : })) {
102 : ERROR_NO_TRACE(
103 : "Cannot refine mesh when not temporally aligned with neighbors.");
104 : }
105 :
106 : // Evaluate AMR p-refinement criteria
107 : auto overall_decision = make_array<volume_dim>(amr::Flag::Undefined);
108 :
109 : using compute_tags = tmpl::remove_duplicates<tmpl::flatten<tmpl::transform<
110 : tmpl::at<typename Metavariables::factory_creation::factory_classes,
111 : Criterion>,
112 : detail::get_tags<tmpl::_1>>>>;
113 : auto observation_box = make_observation_box<compute_tags>(box);
114 :
115 : const auto& refinement_criteria =
116 : db::get<amr::Criteria::Tags::Criteria>(*box);
117 : for (const auto& criterion : refinement_criteria) {
118 : if (criterion->type() == amr::Criteria::Type::h) {
119 : continue;
120 : }
121 : auto decision = criterion->evaluate(observation_box, cache, element_id);
122 : ASSERT(alg::none_of(decision,
123 : [](amr::Flag flag) {
124 : return flag == amr::Flag::Split or
125 : flag == amr::Flag::Join;
126 : }),
127 : "The criterion '" << typeid(*criterion).name()
128 : << "' requested h-refinement, but claims to be "
129 : "for p-refinement.");
130 : const auto& topologies =
131 : db::get<::domain::Tags::Element<volume_dim>>(*box).topologies();
132 : enforce_p_refinement_topology_restrictions(make_not_null(&decision),
133 : topologies);
134 : for (size_t d = 0; d < volume_dim; ++d) {
135 : overall_decision[d] = std::max(overall_decision[d], decision[d]);
136 : }
137 : }
138 : // If no refinement criteria requested p-refinement, then set flag to
139 : // do nothing
140 : for (size_t d = 0; d < volume_dim; ++d) {
141 : if (overall_decision[d] == amr::Flag::Undefined) {
142 : overall_decision[d] = amr::Flag::DoNothing;
143 : }
144 : }
145 :
146 : // Now p-refine
147 : const auto old_mesh_and_element =
148 : std::make_pair(db::get<::domain::Tags::Mesh<volume_dim>>(*box),
149 : db::get<::domain::Tags::Element<volume_dim>>(*box));
150 : const auto& old_mesh = old_mesh_and_element.first;
151 : const auto& verbosity =
152 : db::get<logging::Tags::Verbosity<amr::OptionTags::AmrGroup>>(*box);
153 :
154 : // Enforce policies
155 : const auto& policies = db::get<amr::Tags::Policies>(*box);
156 : amr::enforce_policies(make_not_null(&overall_decision), policies,
157 : element_id, old_mesh);
158 :
159 : if (alg::any_of(overall_decision, [](amr::Flag flag) {
160 : return (flag == amr::Flag::IncreaseResolution or
161 : flag == amr::Flag::DecreaseResolution);
162 : })) {
163 : db::mutate<::domain::Tags::Mesh<volume_dim>>(
164 : [&old_mesh_and_element,
165 : &overall_decision](const gsl::not_null<Mesh<volume_dim>*> mesh) {
166 : *mesh = amr::projectors::p_refined_mesh(
167 : old_mesh_and_element.first, overall_decision,
168 : old_mesh_and_element.second.topologies());
169 : if (not Spectral::is_valid_dg_mesh(*mesh)) {
170 : ERROR("Invalid mesh: " << *mesh);
171 : }
172 : },
173 : box);
174 :
175 : if (verbosity >= Verbosity::Debug) {
176 : Parallel::printf(
177 : "Changing order of element %s: %s -> %s\n", element_id,
178 : old_mesh.extents(),
179 : db::get<::domain::Tags::Mesh<volume_dim>>(*box).extents());
180 : }
181 :
182 : // Run the projectors.
183 : tmpl::for_each<typename Metavariables::amr::projectors>(
184 : [&box, &old_mesh_and_element](auto projector_v) {
185 : using projector = typename decltype(projector_v)::type;
186 : try {
187 : db::mutate_apply<projector>(box, old_mesh_and_element);
188 : } catch (std::exception& e) {
189 : ERROR("Error in AMR projector '"
190 : << pretty_type::get_name<projector>() << "':\n"
191 : << e.what());
192 : }
193 : });
194 : }
195 : }
196 :
197 0 : using is_ready_argument_tags = tmpl::list<>;
198 :
199 : template <typename Metavariables, typename ArrayIndex, typename Component>
200 0 : bool is_ready(Parallel::GlobalCache<Metavariables>& /*cache*/,
201 : const ArrayIndex& /*array_index*/,
202 : const Component* const /*meta*/) const {
203 : return true;
204 : }
205 :
206 1 : bool needs_evolved_variables() const override { return true; }
207 :
208 : // NOLINTNEXTLINE(google-runtime-references)
209 0 : void pup(PUP::er& p) override;
210 : };
211 : } // namespace amr::Events
|