ElementActions.hpp
1 // Distributed under the MIT License.
2 // See LICENSE.txt for details.
3 
4 #pragma once
5 
6 #include <tuple>
7 
9 #include "DataStructures/DataBox/PrefixHelpers.hpp"
10 #include "IO/Observer/Actions.hpp"
11 #include "IO/Observer/Helpers.hpp"
12 #include "IO/Observer/ObservationId.hpp"
13 #include "IO/Observer/ObserverComponent.hpp"
14 #include "IO/Observer/ReductionActions.hpp"
15 #include "IO/Observer/TypeOfObservation.hpp"
16 #include "Informer/Verbosity.hpp"
17 #include "NumericalAlgorithms/Convergence/HasConverged.hpp"
18 #include "NumericalAlgorithms/LinearSolver/Gmres.hpp"
21 #include "Parallel/Invoke.hpp"
22 #include "Parallel/Printf.hpp"
23 #include "Parallel/Reduction.hpp"
24 #include "ParallelAlgorithms/Initialization/MergeIntoDataBox.hpp"
26 #include "Utilities/GetOutput.hpp"
27 #include "Utilities/Gsl.hpp"
28 #include "Utilities/Requires.hpp"
29 #include "Utilities/TMPL.hpp"
30 
31 /// \cond
32 namespace tuples {
33 template <typename...>
34 class TaggedTuple;
35 } // namespace tuples
36 /// \endcond
37 
38 /// Functionality shared between parallel linear solvers that have no global
39 /// synchronization points
41 
42 using reduction_data = Parallel::ReductionData<
43  // Iteration
45  // Residual
47 
48 template <typename OptionsGroup>
50 
51 template <typename FieldsTag, typename OptionsGroup, typename DbTagsList,
52  typename Metavariables, typename ArrayIndex>
53 void contribute_to_residual_observation(
54  const db::DataBox<DbTagsList>& box,
56  const ArrayIndex& array_index) noexcept {
57  using residual_tag =
59  using residual_magnitude_square_tag =
61 
62  const size_t iteration_id =
63  get<LinearSolver::Tags::IterationId<OptionsGroup>>(box);
64  const double residual_magnitude_square =
65  get<residual_magnitude_square_tag>(box);
66  auto& local_observer =
67  *Parallel::get_parallel_component<observers::Observer<Metavariables>>(
68  cache)
69  .ckLocalBranch();
70  Parallel::simple_action<observers::Actions::ContributeReductionData>(
71  local_observer,
72  observers::ObservationId(iteration_id,
74  std::string{"/" + option_name<OptionsGroup>() + "Residuals"},
75  std::vector<std::string>{"Iteration", "Residual"},
76  reduction_data{iteration_id, residual_magnitude_square});
77  if (UNLIKELY(static_cast<int>(
79  static_cast<int>(::Verbosity::Verbose))) {
80  if (iteration_id == 0) {
82  "Linear solver '" + option_name<OptionsGroup>() +
83  "' initialized on element %s. Remaining local residual: %e\n",
84  get_output(array_index), sqrt(residual_magnitude_square));
85  } else {
86  Parallel::printf("Linear solver '" + option_name<OptionsGroup>() +
87  "' iteration %zu done on element %s. Remaining "
88  "local residual: %e\n",
89  iteration_id, get_output(array_index),
90  sqrt(residual_magnitude_square));
91  }
92  }
93 }
94 
95 template <typename FieldsTag, typename OptionsGroup, typename SourceTag>
97  private:
98  using fields_tag = FieldsTag;
99  using operator_applied_to_fields_tag =
101  using source_tag = SourceTag;
102  using residual_tag =
104  using residual_magnitude_square_tag =
106 
107  public:
108  using const_global_cache_tags =
109  tmpl::list<LinearSolver::Tags::Iterations<OptionsGroup>>;
110 
111  template <typename DbTagsList, typename... InboxTags, typename Metavariables,
112  typename ArrayIndex, typename ActionList,
113  typename ParallelComponent>
114  static auto apply(db::DataBox<DbTagsList>& box,
115  const tuples::TaggedTuple<InboxTags...>& /*inboxes*/,
117  const ArrayIndex& /*array_index*/,
118  const ActionList /*meta*/,
119  const ParallelComponent* const /*meta*/) noexcept {
120  return std::make_tuple(
124  residual_magnitude_square_tag,
125  operator_applied_to_fields_tag>,
129  OptionsGroup>>>(
130  std::move(box),
131  // The `PrepareSolve` action populates these tags with initial
132  // values, except for `operator_applied_to_fields_tag` which is
133  // expected to be updated in every iteration of the algorithm
137  }
138 };
139 
140 template <typename OptionsGroup>
142  template <typename ParallelComponent, typename DbTagsList,
143  typename ArrayIndex>
145  register_info(const db::DataBox<DbTagsList>& box,
146  const ArrayIndex& /*array_index*/) noexcept {
147  return {
150  static_cast<double>(
151  db::get<LinearSolver::Tags::IterationId<OptionsGroup>>(box)),
153  }
154 };
155 
156 template <typename FieldsTag, typename OptionsGroup, typename SourceTag>
157 using RegisterElement =
159 
160 template <typename FieldsTag, typename OptionsGroup, typename SourceTag>
161 struct PrepareSolve {
162  private:
163  using fields_tag = FieldsTag;
164  using residual_tag =
166  using residual_magnitude_square_tag =
168 
169  public:
170  using const_global_cache_tags =
171  tmpl::list<LinearSolver::Tags::Verbosity<OptionsGroup>>;
172 
173  template <typename DbTagsList, typename... InboxTags, typename Metavariables,
174  typename ArrayIndex, typename ActionList,
175  typename ParallelComponent>
178  const tuples::TaggedTuple<InboxTags...>& /*inboxes*/,
180  const ArrayIndex& array_index, const ActionList /*meta*/,
181  const ParallelComponent* const /*meta*/) noexcept {
182  db::mutate<LinearSolver::Tags::IterationId<OptionsGroup>,
183  residual_magnitude_square_tag>(
184  make_not_null(&box),
185  [](const gsl::not_null<size_t*> iteration_id,
186  const gsl::not_null<double*> residual_magnitude_square,
187  const auto& residual) noexcept {
188  *iteration_id = 0;
189  *residual_magnitude_square = inner_product(residual, residual);
190  },
191  get<residual_tag>(box));
192  // Observe the initial residual even if no steps are going to be performed
193  contribute_to_residual_observation<FieldsTag, OptionsGroup>(box, cache,
194  array_index);
195  return {std::move(box)};
196  }
197 };
198 
199 template <typename FieldsTag, typename OptionsGroup, typename SourceTag>
200 struct CompleteStep {
201  private:
202  using fields_tag = FieldsTag;
203  using residual_tag =
205  using residual_magnitude_square_tag =
207 
208  public:
209  using const_global_cache_tags =
210  tmpl::list<LinearSolver::Tags::Verbosity<OptionsGroup>>;
211 
212  template <typename DbTagsList, typename... InboxTags, typename Metavariables,
213  typename ArrayIndex, typename ActionList,
214  typename ParallelComponent>
219  const ArrayIndex& array_index, const ActionList /*meta*/,
220  const ParallelComponent* const /*meta*/) noexcept {
221  // Update and observe element-local residual magnitude
222  db::mutate<residual_magnitude_square_tag,
224  make_not_null(&box),
225  [](const gsl::not_null<double*> residual_magnitude_square,
226  const gsl::not_null<size_t*> iteration_id,
227  const auto& residual) noexcept {
228  *residual_magnitude_square = inner_product(residual, residual);
229  ++(*iteration_id);
230  },
231  get<residual_tag>(box));
232  contribute_to_residual_observation<FieldsTag, OptionsGroup>(box, cache,
233  array_index);
234  return {std::move(box)};
235  }
236 };
237 
238 } // namespace LinearSolver::async_solvers
Parallel::ConstGlobalCache
Definition: ElementReceiveInterpPoints.hpp:16
observers::ObservationId
A type-erased identifier that combines the identifier's type and hash used to uniquely identify an ob...
Definition: ObservationId.hpp:42
std::string
UNLIKELY
#define UNLIKELY(x)
Definition: Gsl.hpp:73
get
constexpr Tag::type & get(Variables< TagList > &v) noexcept
Return Tag::type pointing into the contiguous array.
Definition: Variables.hpp:689
Parallel::ReductionDatum
The data to be reduced, and invokables to be called whenever two reduction messages are combined and ...
Definition: Reduction.hpp:64
InnerProduct.hpp
std::pair
std::vector< std::string >
Parallel::printf
void printf(const std::string &format, Args &&... args)
Print an atomic message to stdout with C printf usage.
Definition: Printf.hpp:103
LinearSolver::Tags::IterationId
Holds an IterationId that identifies a step in the linear solver algorithm.
Definition: Tags.hpp:84
db::AddComputeTags
tmpl::flatten< tmpl::list< Tags... > > AddComputeTags
List of Compute Item Tags to add to the DataBox.
Definition: DataBox.hpp:1171
tuple
LinearSolver::async_solvers::RegisterObservers
Definition: ElementActions.hpp:141
LinearSolver::async_solvers::InitializeElement
Definition: ElementActions.hpp:96
db::add_tag_prefix
Prefix< DataBox_detail::dispatch_add_tag_prefix_impl< Prefix, Tag, Args... >, Args... > add_tag_prefix
Definition: PrefixHelpers.hpp:145
Tags.hpp
db::mutate
void mutate(const gsl::not_null< DataBox< TagList > * > box, Invokable &&invokable, Args &&... args) noexcept
Allows changing the state of one or more non-computed elements in the DataBox.
Definition: DataBox.hpp:1020
std::sqrt
T sqrt(T... args)
get_output
std::string get_output(const T &t) noexcept
Get the streamed output of t as a std::string
Definition: GetOutput.hpp:14
funcl::Sqrt
Functional for computing sqrt on an object.
Definition: Functional.hpp:273
db::apply
constexpr auto apply(F &&f, const DataBox< BoxTags > &box, Args &&... args) noexcept
Apply the invokable f with argument Tags TagsList from DataBox box
Definition: DataBox.hpp:1443
Printf.hpp
LinearSolver::inner_product
double inner_product(const Lhs &lhs, const Rhs &rhs) noexcept
The local part of the Euclidean inner product on the vector space w.r.t. which the addition and scala...
Definition: InnerProduct.hpp:71
DataBox.hpp
LinearSolver::async_solvers::ElementObservationType
Definition: ElementActions.hpp:49
db::item_type
typename DataBox_detail::item_type_impl< TagList, Tag >::type item_type
Get the type that can be written to the Tag. If it is a base tag then a TagList must be passed as a s...
Definition: DataBoxTag.hpp:246
tuples::TaggedTuple
An associative container that is indexed by structs.
Definition: TaggedTuple.hpp:272
db::AddSimpleTags
tmpl::flatten< tmpl::list< Tags... > > AddSimpleTags
List of Tags to add to the DataBox.
Definition: DataBox.hpp:1150
LinearSolver::Tags::Verbosity
Definition: Tags.hpp:386
Gsl.hpp
LinearSolver::async_solvers
Functionality shared between parallel linear solvers that have no global synchronization points.
Definition: ElementActions.hpp:40
Requires.hpp
observers::TypeOfObservation::Reduction
@ Reduction
The sender will only perform reduction observations.
Initialization::merge_into_databox
auto merge_into_databox(db::DataBox< DbTagsList > &&box, Args &&... args) noexcept
Add tags that are not yet in the DataBox.
Definition: MergeIntoDataBox.hpp:133
LinearSolver::Tags::HasConvergedByIterationsCompute
Employs the LinearSolver::Tags::Iterations to determine the linear solver has converged once it has c...
Definition: Tags.hpp:284
LinearSolver::async_solvers::PrepareSolve
Definition: ElementActions.hpp:161
LinearSolver::Tags::ResidualCompute
Compute the residual from the SourceTag and the db::add_tag_prefix<LinearSolver::Tags::OperatorAppl...
Definition: Tags.hpp:110
make_not_null
gsl::not_null< T * > make_not_null(T *ptr) noexcept
Construct a not_null from a pointer. Often this will be done as an implicit conversion,...
Definition: Gsl.hpp:880
db::DataBox
Definition: InterpolationTargetWedgeSectionTorus.hpp:24
std::numeric_limits
TMPL.hpp
ConstGlobalCache.hpp
gsl::not_null
Require a pointer to not be a nullptr
Definition: Gsl.hpp:183
observers::Actions::RegisterWithObservers
Registers itself with the local observer parallel component so the observer knows to expect data from...
Definition: Actions.hpp:232
LinearSolver::async_solvers::CompleteStep
Definition: ElementActions.hpp:200