BoundaryConditions.hpp
1 // Distributed under the MIT License.
2 // See LICENSE.txt for details.
3 
4 #pragma once
5 
7 
8 #include <array>
9 #include <memory>
10 #include <optional>
11 #include <random>
12 #include <regex>
13 #include <stdexcept>
14 #include <string>
15 
17 #include "DataStructures/DataBox/PrefixHelpers.hpp"
19 #include "DataStructures/DataVector.hpp"
20 #include "DataStructures/Index.hpp"
21 #include "DataStructures/Tensor/EagerMath/Magnitude.hpp"
24 #include "Domain/BoundaryConditions/BoundaryCondition.hpp"
25 #include "Domain/BoundaryConditions/Periodic.hpp"
26 #include "Domain/Tags.hpp"
27 #include "Evolution/BoundaryConditions/Type.hpp"
28 #include "Evolution/DiscontinuousGalerkin/Actions/ComputeTimeDerivativeHelpers.hpp"
29 #include "Evolution/DiscontinuousGalerkin/Actions/NormalCovectorAndMagnitude.hpp"
30 #include "Framework/Pypp.hpp"
32 #include "Framework/TestCreation.hpp"
35 #include "Helpers/Evolution/DiscontinuousGalerkin/NormalVectors.hpp"
36 #include "Helpers/Evolution/DiscontinuousGalerkin/Range.hpp"
38 #include "Utilities/Gsl.hpp"
39 #include "Utilities/NoSuchType.hpp"
40 #include "Utilities/TMPL.hpp"
41 #include "Utilities/TaggedTuple.hpp"
42 
43 namespace TestHelpers::evolution::dg {
44 /// Tags for testing DG code
45 namespace Tags {
46 /// Tag for a `TaggedTuple` that holds the name of the python function that
47 /// gives the result for computing what `Tag` should be.
48 ///
49 /// If a Tag is only needed for or different for a specific boundary correction,
50 /// then the `BoundaryCorrection` template parameter must be set.
51 ///
52 /// The `Tag` template parameter is an input to the `dg_package_data` function.
53 /// That is, the tag is one of the evolved variables, fluxes, or any other
54 /// `dg_package_data_temporary_tags` that the boundary correction needs.
55 template <typename Tag, typename BoundaryCorrection = NoSuchType>
57  using tag = Tag;
58  using boundary_correction = BoundaryCorrection;
59  using type = std::string;
60 };
61 
62 /// The name of the python function that returns the error message.
63 ///
64 /// If `BoundaryCorrection` is `NoSuchType` then the same function will be used
65 /// for all boundary corrections.
66 ///
67 /// The python function must return `None` if there shouldn't be an error
68 /// message.
69 template <typename BoundaryCorrection = NoSuchType>
71  using boundary_correction = BoundaryCorrection;
72  using type = std::string;
73 };
74 } // namespace Tags
75 
76 namespace detail {
77 template <typename ConversionClassList, typename... Args>
78 static std::optional<std::string> call_for_error_message(
79  const std::string& module_name, const std::string& function_name,
80  const Args&... t) {
81  static_assert(sizeof...(Args) > 0,
82  "Call to python which returns a Tensor of DataVectors must "
83  "pass at least one argument");
84  using ReturnType = std::optional<std::string>;
85 
86  PyObject* python_module = PyImport_ImportModule(module_name.c_str());
87  if (python_module == nullptr) {
88  PyErr_Print();
89  throw std::runtime_error{std::string("Could not find python module.\n") +
90  module_name};
91  }
92  PyObject* func = PyObject_GetAttrString(python_module, function_name.c_str());
93  if (func == nullptr or not PyCallable_Check(func)) {
94  if (PyErr_Occurred()) {
95  PyErr_Print();
96  }
97  throw std::runtime_error{"Could not find python function in module.\n"};
98  }
99 
100  const std::array<size_t, sizeof...(Args)> arg_sizes{
101  {pypp::detail::ContainerPackAndUnpack<
102  Args, ConversionClassList>::get_size(t)...}};
103  const size_t npts = *std::max_element(arg_sizes.begin(), arg_sizes.end());
104  for (size_t i = 0; i < arg_sizes.size(); ++i) {
105  if (gsl::at(arg_sizes, i) != 1 and gsl::at(arg_sizes, i) != npts) {
106  ERROR("Each argument must return size 1 or "
107  << npts
108  << " (the number of points in the DataVector), but argument number "
109  << i << " has size " << gsl::at(arg_sizes, i));
110  }
111  }
112  ReturnType result{};
113 
114  for (size_t s = 0; s < npts; ++s) {
115  PyObject* args = pypp::make_py_tuple(
116  pypp::detail::ContainerPackAndUnpack<Args, ConversionClassList>::unpack(
117  t, s)...);
118  auto ret =
119  pypp::detail::call_work<typename pypp::detail::ContainerPackAndUnpack<
120  ReturnType, ConversionClassList>::unpacked_container>(python_module,
121  func, args);
122  if (ret.has_value()) {
123  result = std::move(ret);
124  break;
125  }
126  }
127  Py_DECREF(func); // NOLINT
128  Py_DECREF(python_module); // NOLINT
129  return result;
130 }
131 
132 template <typename BoundaryCorrection, typename... PythonFunctionNameTags>
133 const std::string& get_python_error_message_function(
135  python_boundary_condition_functions) {
136  return tuples::get<tmpl::conditional_t<
137  tmpl::list_contains_v<tmpl::list<PythonFunctionNameTags...>,
138  Tags::PythonFunctionForErrorMessage<NoSuchType>>,
139  Tags::PythonFunctionForErrorMessage<NoSuchType>,
140  Tags::PythonFunctionForErrorMessage<BoundaryCorrection>>>(
141  python_boundary_condition_functions);
142 }
143 
144 template <typename Tag, typename BoundaryCorrection,
145  typename... PythonFunctionNameTags>
146 const std::string& get_python_tag_function(
148  python_boundary_condition_functions) {
149  return tuples::get<tmpl::conditional_t<
150  tmpl::list_contains_v<tmpl::list<PythonFunctionNameTags...>,
151  Tags::PythonFunctionName<Tag, NoSuchType>>,
152  Tags::PythonFunctionName<Tag, NoSuchType>,
153  Tags::PythonFunctionName<Tag, BoundaryCorrection>>>(
154  python_boundary_condition_functions);
155 }
156 
157 template <typename BoundaryConditionHelper, typename AllTagsOnFaceList,
158  typename... TagsFromFace, typename... VolumeArgs>
159 void apply_boundary_condition_impl(
160  BoundaryConditionHelper& boundary_condition_helper,
161  const Variables<AllTagsOnFaceList>& fields_on_interior_face,
162  tmpl::list<TagsFromFace...> /*meta*/,
163  const VolumeArgs&... volume_args) noexcept {
164  boundary_condition_helper(get<TagsFromFace>(fields_on_interior_face)...,
165  volume_args...);
166 }
167 
168 template <typename System, typename ConversionClassList,
169  typename BoundaryCorrection, typename... PythonFunctionNameTags,
170  typename BoundaryCondition, size_t FaceDim, typename DbTagsList,
171  typename... RangeTags, typename... EvolvedVariablesTags,
172  typename... BoundaryCorrectionPackagedDataInputTags,
173  typename... BoundaryConditionVolumeTags>
174 void test_boundary_condition_with_python_impl(
175  const gsl::not_null<std::mt19937*> generator,
176  const std::string& python_module,
178  python_boundary_condition_functions,
179  const BoundaryCondition& boundary_condition,
180  const Index<FaceDim>& face_points,
181  const db::DataBox<DbTagsList>& box_of_volume_data,
182  const tuples::TaggedTuple<Tags::Range<RangeTags>...>& ranges,
183  const bool use_moving_mesh, tmpl::list<EvolvedVariablesTags...> /*meta*/,
184  tmpl::list<BoundaryCorrectionPackagedDataInputTags...> /*meta*/,
185  tmpl::list<BoundaryConditionVolumeTags...> /*meta*/, const double epsilon) {
186  CAPTURE(FaceDim);
187  CAPTURE(python_module);
188  CAPTURE(face_points);
189  CAPTURE(use_moving_mesh);
190  CAPTURE(epsilon);
191  CAPTURE(pretty_type::short_name<BoundaryCorrection>());
192  const size_t number_of_points_on_face = face_points.product();
193 
194  using variables_tag = typename System::variables_tag;
195  using variables_tags = typename variables_tag::tags_list;
196  using flux_variables = typename System::flux_variables;
197  using dt_variables_tags = db::wrap_tags_in<::Tags::dt, variables_tags>;
198 
199  constexpr bool uses_ghost =
200  BoundaryCondition::bc_type ==
201  ::evolution::BoundaryConditions::Type::Ghost or
202  BoundaryCondition::bc_type ==
203  ::evolution::BoundaryConditions::Type::GhostAndTimeDerivative;
204  constexpr bool uses_time_derivative_condition =
205  BoundaryCondition::bc_type ==
206  ::evolution::BoundaryConditions::Type::TimeDerivative or
207  BoundaryCondition::bc_type ==
208  ::evolution::BoundaryConditions::Type::GhostAndTimeDerivative;
209 
210  // List that holds the inverse spatial metric if it's needed
211  using inverse_spatial_metric_list =
212  ::evolution::dg::Actions::detail::inverse_spatial_metric_tag<System>;
213  constexpr bool has_inv_spatial_metric =
214  ::evolution::dg::Actions::detail::has_inverse_spatial_metric_tag_v<
215  System>;
216 
217  // Set up tags for boundary conditions
218  using bcondition_interior_temp_tags =
219  typename BoundaryCondition::dg_interior_temporary_tags;
220  using bcondition_interior_prim_tags =
221  typename ::evolution::dg::Actions::detail::get_primitive_vars<
222  System::has_primitive_and_conservative_vars>::
223  template boundary_condition_interior_tags<BoundaryCondition>;
224  using bcondition_interior_evolved_vars_tags =
225  typename BoundaryCondition::dg_interior_evolved_variables_tags;
226  using bcondition_interior_dt_evolved_vars_tags =
227  ::evolution::dg::Actions::detail::get_dt_vars_from_boundary_condition<
228  BoundaryCondition>;
229  using bcondition_interior_deriv_evolved_vars_tags =
230  ::evolution::dg::Actions::detail::get_deriv_vars_from_boundary_condition<
231  BoundaryCondition>;
232  using bcondition_interior_tags = tmpl::remove_duplicates<tmpl::append<
233  tmpl::conditional_t<has_inv_spatial_metric,
234  tmpl::list<::evolution::dg::Actions::detail::
235  NormalVector<FaceDim + 1>>,
236  tmpl::list<>>,
237  bcondition_interior_evolved_vars_tags, bcondition_interior_prim_tags,
238  bcondition_interior_temp_tags, bcondition_interior_dt_evolved_vars_tags,
239  bcondition_interior_deriv_evolved_vars_tags>>;
240 
241  std::uniform_real_distribution<> dist(-1., 1.);
242 
243  // Fill all fields with random values in [-1,1), then, for each tag with a
244  // specified range, overwrite with new random values in [min,max)
245  Variables<tmpl::remove_duplicates<
246  tmpl::append<bcondition_interior_tags, inverse_spatial_metric_list>>>
247  interior_face_fields{number_of_points_on_face};
248  fill_with_random_values(make_not_null(&interior_face_fields), generator,
249  make_not_null(&dist));
250  tmpl::for_each<tmpl::list<RangeTags...>>([&generator, &interior_face_fields,
251  &ranges](auto tag_v) {
252  using tag = tmpl::type_from<decltype(tag_v)>;
253  const std::array<double, 2>& range = tuples::get<Tags::Range<tag>>(ranges);
254  std::uniform_real_distribution<> local_dist(range[0], range[1]);
255  fill_with_random_values(make_not_null(&get<tag>(interior_face_fields)),
256  generator, make_not_null(&local_dist));
257  });
258 
259  auto interior_normal_covector = // Unnormalized at this point
261  tnsr::i<DataVector, FaceDim + 1, Frame::Inertial>>(
262  generator, make_not_null(&dist), number_of_points_on_face);
263  if constexpr (has_inv_spatial_metric) {
264  auto& inv_spatial_metric =
265  get<tmpl::front<inverse_spatial_metric_list>>(interior_face_fields);
266  detail::adjust_spatial_metric_or_inverse(
267  make_not_null(&inv_spatial_metric));
268  auto& normal_vector =
269  get<::evolution::dg::Actions::detail::NormalVector<FaceDim + 1>>(
270  interior_face_fields);
271  detail::normalize_vector_and_covector(
272  make_not_null(&interior_normal_covector), make_not_null(&normal_vector),
273  inv_spatial_metric);
274  } else {
275  const Scalar<DataVector> magnitude = ::magnitude(interior_normal_covector);
276  for (DataVector& component : interior_normal_covector) {
277  component /= get(magnitude);
278  }
279  }
280 
281  // Set a random mesh velocity. We allow for velocities above 1 to test "faster
282  // than light" mesh movement.
284  if (use_moving_mesh) {
285  std::uniform_real_distribution<> local_dist(-2., 2.);
286  face_mesh_velocity =
287  make_with_random_values<tnsr::I<DataVector, FaceDim + 1>>(
288  generator, make_not_null(&local_dist), number_of_points_on_face);
289  }
290 
291  if constexpr (BoundaryCondition::bc_type ==
292  ::evolution::BoundaryConditions::Type::Outflow) {
293  // Outflow boundary conditions only check that all characteristic speeds
294  // are directed out of the element. If there are any inward directed
295  // fields then the boundary condition should error.
296  const auto apply_bc =
297  [&boundary_condition, &face_mesh_velocity, &interior_normal_covector,
298  &python_boundary_condition_functions,
299  &python_module](const auto&... face_and_volume_args) noexcept {
300  const std::optional<std::string> error_msg =
301  boundary_condition.dg_outflow(face_mesh_velocity,
302  interior_normal_covector,
303  face_and_volume_args...);
304  const std::string& python_error_msg_function =
305  get_python_error_message_function<BoundaryCorrection>(
306  python_boundary_condition_functions);
307  const auto python_error_message =
308  call_for_error_message<ConversionClassList>(
309  python_module, python_error_msg_function, face_mesh_velocity,
310  interior_normal_covector, face_and_volume_args...);
311  CAPTURE(python_error_msg_function);
312  CAPTURE(python_error_message.value_or(""));
313  CAPTURE(error_msg.value_or(""));
314  REQUIRE(python_error_message.has_value() == error_msg.has_value());
315  if (python_error_message.has_value() and error_msg.has_value()) {
316  std::smatch matcher{};
317  CHECK(std::regex_search(*error_msg, matcher,
318  std::regex{*python_error_message}));
319  }
320  };
321  apply_boundary_condition_impl(
322  apply_bc, interior_face_fields, bcondition_interior_tags{},
323  db::get<BoundaryConditionVolumeTags>(box_of_volume_data)...);
324  }
325 
326  if constexpr (uses_time_derivative_condition) {
327  Variables<dt_variables_tags> time_derivative_correction{
328  number_of_points_on_face};
329  auto apply_bc = [&boundary_condition, epsilon, &face_mesh_velocity,
330  &interior_normal_covector,
331  &python_boundary_condition_functions, &python_module,
332  &time_derivative_correction](
333  const auto&... interior_face_and_volume_args) {
334  const std::optional<std::string> error_msg =
335  boundary_condition.dg_time_derivative(
337  time_derivative_correction))...,
338  face_mesh_velocity, interior_normal_covector,
339  interior_face_and_volume_args...);
340 
341  const std::string& python_error_msg_function =
342  get_python_error_message_function<BoundaryCorrection>(
343  python_boundary_condition_functions);
344  const auto python_error_message =
345  call_for_error_message<ConversionClassList>(
346  python_module, python_error_msg_function, face_mesh_velocity,
347  interior_normal_covector, interior_face_and_volume_args...);
348  CAPTURE(python_error_msg_function);
349  CAPTURE(python_error_message.value_or(""));
350  CAPTURE(error_msg.value_or(""));
351  REQUIRE(python_error_message.has_value() == error_msg.has_value());
352  if (python_error_message.has_value() and error_msg.has_value()) {
353  std::smatch matcher{};
354  CHECK(std::regex_search(*error_msg, matcher,
355  std::regex{*python_error_message}));
356  }
357 
358  // Check that the values were set correctly
359  tmpl::for_each<dt_variables_tags>([&](auto dt_var_tag_v) {
360  using DtVarTag = tmpl::type_from<decltype(dt_var_tag_v)>;
361  // Use NoSuchType for the BoundaryCorrection since dt-type corrections
362  // should be boundary correction agnostic.
363  const std::string& python_tag_function =
364  get_python_tag_function<DtVarTag, NoSuchType>(
365  python_boundary_condition_functions);
366  CAPTURE(python_tag_function);
367  CAPTURE(pretty_type::short_name<DtVarTag>());
368  typename DtVarTag::type python_result{};
369  try {
370  python_result =
371  pypp::call<typename DtVarTag::type, ConversionClassList>(
372  python_module, python_tag_function, face_mesh_velocity,
373  interior_normal_covector, interior_face_and_volume_args...);
374  } catch (const std::exception& e) {
375  INFO("Failed python call with '" << e.what() << "'");
376  // Use REQUIRE(false) to print all the CAPTURE variables
377  REQUIRE(false);
378  }
380  get<DtVarTag>(time_derivative_correction), python_result,
381  Approx::custom().epsilon(epsilon).scale(1.0));
382  });
383  };
384  apply_boundary_condition_impl(
385  apply_bc, interior_face_fields, bcondition_interior_tags{},
386  db::get<BoundaryConditionVolumeTags>(box_of_volume_data)...);
387  }
388 
389  if constexpr (uses_ghost) {
390  using fluxes_tags =
391  db::wrap_tags_in<::Tags::Flux, flux_variables,
392  tmpl::size_t<FaceDim + 1>, Frame::Inertial>;
393  using correction_temp_tags =
394  typename BoundaryCorrection::dg_package_data_temporary_tags;
395  using correction_prim_tags =
396  typename ::evolution::dg::Actions::detail::get_primitive_vars<
397  System::has_primitive_and_conservative_vars>::
398  template f<BoundaryCorrection>;
399  using tags_on_exterior_face = tmpl::append<
400  variables_tags, fluxes_tags, correction_temp_tags, correction_prim_tags,
401  inverse_spatial_metric_list,
402  tmpl::list<
403  ::evolution::dg::Actions::detail::OneOverNormalVectorMagnitude,
404  ::evolution::dg::Actions::detail::NormalVector<FaceDim + 1>,
406  Variables<tags_on_exterior_face> exterior_face_fields{
407  number_of_points_on_face};
408  auto apply_bc = [&boundary_condition, epsilon, &exterior_face_fields,
409  &face_mesh_velocity, &interior_normal_covector,
410  &python_boundary_condition_functions, &python_module](
411  const auto&... interior_face_and_volume_args) {
412  std::optional<std::string> error_msg{};
413  if constexpr (has_inv_spatial_metric) {
414  error_msg = boundary_condition.dg_ghost(
415  make_not_null(&get<BoundaryCorrectionPackagedDataInputTags>(
416  exterior_face_fields))...,
417  make_not_null(&get<typename System::inverse_spatial_metric_tag>(
418  exterior_face_fields)),
419  face_mesh_velocity, interior_normal_covector,
420  interior_face_and_volume_args...);
421  } else {
422  error_msg = boundary_condition.dg_ghost(
423  make_not_null(&get<BoundaryCorrectionPackagedDataInputTags>(
424  exterior_face_fields))...,
425  face_mesh_velocity, interior_normal_covector,
426  interior_face_and_volume_args...);
427  }
428 
429  const std::string& python_error_msg_function =
430  get_python_error_message_function<BoundaryCorrection>(
431  python_boundary_condition_functions);
432  const auto python_error_message =
433  call_for_error_message<ConversionClassList>(
434  python_module, python_error_msg_function, face_mesh_velocity,
435  interior_normal_covector, interior_face_and_volume_args...);
436  CAPTURE(python_error_msg_function);
437  CAPTURE(python_error_message.value_or(""));
438  CAPTURE(error_msg.value_or(""));
439  REQUIRE(python_error_message.has_value() == error_msg.has_value());
440  if (python_error_message.has_value() and error_msg.has_value()) {
441  std::smatch matcher{};
442  CHECK(std::regex_search(*error_msg, matcher,
443  std::regex{*python_error_message}));
444  }
445  if (python_error_message.has_value()) {
446  return;
447  }
448 
449  // Check that the values were set correctly
450  tmpl::for_each<tmpl::list<BoundaryCorrectionPackagedDataInputTags...>>(
451  [&](auto boundary_correction_tag_v) {
452  using BoundaryCorrectionTag =
453  tmpl::type_from<decltype(boundary_correction_tag_v)>;
454  const std::string& python_tag_function =
455  get_python_tag_function<BoundaryCorrectionTag,
456  BoundaryCorrection>(
457  python_boundary_condition_functions);
458  CAPTURE(python_tag_function);
459  CAPTURE(pretty_type::short_name<BoundaryCorrectionTag>());
460  typename BoundaryCorrectionTag::type python_result{};
461  try {
462  python_result = pypp::call<typename BoundaryCorrectionTag::type,
463  ConversionClassList>(
464  python_module, python_tag_function, face_mesh_velocity,
465  interior_normal_covector, interior_face_and_volume_args...);
466  } catch (const std::exception& e) {
467  INFO("Failed python call with '" << e.what() << "'");
468  // Use REQUIRE(false) to print all the CAPTURE variables
469  REQUIRE(false);
470  }
472  get<BoundaryCorrectionTag>(exterior_face_fields), python_result,
473  Approx::custom().epsilon(epsilon).scale(1.0));
474  });
475  };
476  // Since any of the variables in `exterior_face_fields` could be mutated
477  // from their projected state, we don't pass in the interior tags explicitly
478  // since that will likely give a false sense of "no aliasing"
479  apply_boundary_condition_impl(
480  apply_bc, interior_face_fields, bcondition_interior_tags{},
481  db::get<BoundaryConditionVolumeTags>(box_of_volume_data)...);
482  }
483 }
484 
485 template <typename T, typename = std::void_t<>>
486 struct get_boundary_conditions_impl {
487  using type = tmpl::list<>;
488 };
489 
490 template <typename T>
491 struct get_boundary_conditions_impl<
492  T, std::void_t<typename T::boundary_conditions_base>> {
493  using type = typename T::boundary_conditions_base::creatable_classes;
494 };
495 
496 template <typename T>
497 using get_boundary_conditions = typename get_boundary_conditions_impl<T>::type;
498 } // namespace detail
499 
500 /*!
501  * \brief Test a boundary condition against python code and that it satisfies
502  * the required interface.
503  *
504  * The boundary conditions return a `std::optional<std::string>` that is the
505  * error message. For ghost cell boundary conditions they must also return the
506  * arguments needed by the boundary correction's `dg_package_data` function by
507  * `gsl::not_null`. Time derivative boundary conditions return the correction
508  * added to the time derivatives by `gsl::not_null`, while outflow boundary
509  * conditions should only check that the boundary is actually an outflow
510  * boundary. Therefore, the comparison implementation in python must have a
511  * function for each of these. Which function is called is specified using a
512  * `python_boundary_condition_functions` and `Tags::PythonFunctionName`.
513  *
514  * The specific boundary condition and system need to be explicitly given.
515  * Once the boundary condition and boundary correction base classes are
516  * listed/available in the `System` class we can remove needing to specify
517  * those. The `ConversionClassList` is forwarded to `pypp` to allow custom
518  * conversions of classes to python, such as analytic solutions or equations of
519  * state.
520  *
521  * - The random number generator is passed in so that the seed can be easily
522  * controlled externally. Use, e.g. `MAKE_GENERATOR(gen)` to create a
523  * generator.
524  * - The python function names are given in a `TaggedTuple`. A
525  * `TestHelpers::evolution::dg::Tags::PythonFunctionForErrorMessage` tag must
526  * be given for specifying the error message that the boundary condition
527  * returns. In many cases this function will simply return `None`. The tag can
528  * optionally specify the boundary correction for which it is to be used, so
529  * that a different error message could be printed if a different boundary
530  * correction is used.
531  * The `TestHelpers::evolution::dg::Tags::PythonFunctionName` tag is used to
532  * give the name of the python function for each return argument. The tags are
533  * the input to the `package_data` function for Ghost boundary conditions, and
534  * `::Tags::dt<evolved_var_tag>` for time derivative boundary conditions.
535  * - `factory_string` is a string used to create the boundary condition from the
536  * factory
537  * - `face_points` is the grid points on the interface. Generally 5 grid points
538  * per dimension in 2d and 3d is recommended to catch indexing errors. In 1d
539  * there is only ever one point on the interface
540  * - `box_of_volume_data` is a `db::DataBox` that contains all of the
541  * `dg_gridless_tags` of the boundary condition. This is not a `TaggedTuple`
542  * so that all the different types of tags, like base tags, can be supported
543  * and easily tested.
544  * - `ranges` is a `TaggedTuple` of
545  * `TestHelpers::evolution::dg::Tags::Range<tag>` specifying a custom range in
546  * which to generate the random values. This can be used for ensuring that
547  * positive quantities are randomly generated on the interval
548  * `[lower_bound,upper_bound)`, choosing `lower_bound` to be `0` or some small
549  * number. The default interval if a tag is not listed is `[-1,1)`.
550  * - `epsilon` is the relative tolerance to use in the random tests
551  */
552 template <typename BoundaryCondition, typename BoundaryConditionBase,
553  typename System, typename BoundaryCorrectionsList,
554  typename ConversionClassList = tmpl::list<>, size_t FaceDim,
555  typename DbTagsList, typename... RangeTags,
556  typename... PythonFunctionNameTags>
557 void test_boundary_condition_with_python(
558  const gsl::not_null<std::mt19937*> generator,
559  const std::string& python_module,
561  python_boundary_condition_functions,
562  const std::string& factory_string, const Index<FaceDim>& face_points,
563  const db::DataBox<DbTagsList>& box_of_volume_data,
564  const tuples::TaggedTuple<Tags::Range<RangeTags>...>& ranges,
565  const double epsilon = 1.0e-12) {
566  PUPable_reg(BoundaryCondition);
567  static_assert(std::is_final_v<std::decay_t<BoundaryCondition>>,
568  "All boundary condition classes must be marked `final`.");
569  using variables_tags = typename System::variables_tag::tags_list;
570  using flux_variables = typename System::flux_variables;
571  using fluxes_tags =
575  boundary_condition =
576  TestHelpers::test_creation<std::unique_ptr<BoundaryConditionBase>>(
577  factory_string);
578 
579  REQUIRE_FALSE(
580  domain::BoundaryConditions::is_periodic(boundary_condition->get_clone()));
581 
582  tmpl::for_each<BoundaryCorrectionsList>(
583  [&boundary_condition, &box_of_volume_data, epsilon, &face_points,
584  &generator, &python_boundary_condition_functions, &python_module,
585  &ranges](auto boundary_correction_v) {
586  using BoundaryCorrection =
587  tmpl::type_from<decltype(boundary_correction_v)>;
588  using package_data_input_tags = tmpl::append<
589  variables_tags, fluxes_tags,
590  typename BoundaryCorrection::dg_package_data_temporary_tags,
591  typename ::evolution::dg::Actions::detail::get_primitive_vars<
592  System::has_primitive_and_conservative_vars>::
593  template f<BoundaryCorrection>>;
594  using boundary_condition_dg_gridless_tags =
595  typename BoundaryCondition::dg_gridless_tags;
596  for (const auto use_moving_mesh : {false, true}) {
597  detail::test_boundary_condition_with_python_impl<
598  System, ConversionClassList, BoundaryCorrection>(
599  generator, python_module, python_boundary_condition_functions,
600  dynamic_cast<const BoundaryCondition&>(*boundary_condition),
601  face_points, box_of_volume_data, ranges, use_moving_mesh,
602  variables_tags{}, package_data_input_tags{},
603  boundary_condition_dg_gridless_tags{}, epsilon);
604  // Now serialize and deserialize and test again
605  INFO("Test boundary condition after serialization.");
606  const auto deserialized_bc =
607  serialize_and_deserialize(boundary_condition);
608  detail::test_boundary_condition_with_python_impl<
609  System, ConversionClassList, BoundaryCorrection>(
610  generator, python_module, python_boundary_condition_functions,
611  dynamic_cast<const BoundaryCondition&>(*deserialized_bc),
612  face_points, box_of_volume_data, ranges, use_moving_mesh,
613  variables_tags{}, package_data_input_tags{},
614  boundary_condition_dg_gridless_tags{}, epsilon);
615  }
616  });
617 }
618 
619 /// Test that a boundary condition is correctly identified as periodic.
620 template <typename BoundaryCondition, typename BoundaryConditionBase>
621 void test_periodic_condition(const std::string& factory_string) {
622  PUPable_reg(BoundaryCondition);
623  const auto boundary_condition =
624  TestHelpers::test_creation<std::unique_ptr<BoundaryConditionBase>>(
625  factory_string);
626  REQUIRE(typeid(*boundary_condition.get()) == typeid(BoundaryCondition));
627  const auto bc_clone = boundary_condition->get_clone();
629  const auto bc_deserialized = serialize_and_deserialize(bc_clone);
630  CHECK(domain::BoundaryConditions::is_periodic(bc_deserialized));
631 }
632 } // namespace TestHelpers::evolution::dg
gsl::at
constexpr T & at(std::array< T, N > &arr, Size index)
Retrieve a entry from a container, with checks in Debug mode that the index being retrieved is valid.
Definition: Gsl.hpp:125
regex
std::string
Frame::Inertial
Definition: IndexType.hpp:44
std::exception
Tags.hpp
TestingFramework.hpp
random
fill_with_random_values
void fill_with_random_values(const gsl::not_null< T * > data, const gsl::not_null< UniformRandomBitGenerator * > generator, const gsl::not_null< RandomNumberDistribution * > distribution) noexcept
Fill an existing data structure with random values.
Definition: MakeWithRandomValues.hpp:165
MakeWithRandomValues.hpp
db::get
const auto & get(const DataBox< TagList > &box) noexcept
Retrieve the item with tag Tag from the DataBox.
Definition: DataBox.hpp:791
Index
Definition: Index.hpp:31
CHECK_ITERABLE_CUSTOM_APPROX
#define CHECK_ITERABLE_CUSTOM_APPROX(a, b, appx)
Same as CHECK_ITERABLE_APPROX with user-defined Approx. The third argument should be of type Approx.
Definition: TestingFramework.hpp:134
Tags::Flux
Prefix indicating a flux.
Definition: Prefixes.hpp:40
TestHelpers::evolution::dg::Tags::PythonFunctionName
Tag for a TaggedTuple that holds the name of the python function that gives the result for computing ...
Definition: BoundaryConditions.hpp:56
get_size
decltype(auto) get_size(const T &t, SizeFunction size=GetContainerSize{}) noexcept
Retrieve the size of t if t.size() is a valid expression, otherwise if T is fundamental or a std::com...
Definition: ContainerHelpers.hpp:145
std::uniform_real_distribution
pypp::call
ReturnType call(const std::string &module_name, const std::string &function_name, const Args &... t)
Calls a Python function from a module/file with given parameters.
Definition: Pypp.hpp:494
TestHelpers.hpp
stdexcept
ERROR
#define ERROR(m)
prints an error message to the standard error stream and aborts the program.
Definition: Error.hpp:37
Pypp.hpp
DataBox.hpp
PyppFundamentals.hpp
array
Index.hpp
std::runtime_error
tuples::TaggedTuple
An associative container that is indexed by structs.
Definition: TaggedTuple.hpp:271
DataVector
Stores a collection of function values.
Definition: DataVector.hpp:46
std::regex
memory
serialize_and_deserialize
T serialize_and_deserialize(const T &t)
Serializes and deserializes an object t of type T
Definition: TestHelpers.hpp:45
std::decay_t
Variables.hpp
pypp::make_py_tuple
PyObject * make_py_tuple(const Args &... t)
Create a python tuple from Args.
Definition: PyppFundamentals.hpp:72
Tags::dt
Prefix indicating a time derivative.
Definition: Prefixes.hpp:29
Scalar
Tensor< T, Symmetry<>, index_list<> > Scalar
Definition: TypeAliases.hpp:21
make_with_random_values
ReturnType make_with_random_values(const gsl::not_null< UniformRandomBitGenerator * > generator, const gsl::not_null< RandomNumberDistribution * > distribution, const T &used_for_size) noexcept
Make a data structure and fill it with random values.
Definition: MakeWithRandomValues.hpp:185
Gsl.hpp
domain::BoundaryConditions::is_periodic
bool is_periodic(const std::unique_ptr< BoundaryCondition > &boundary_condition) noexcept
Check if a boundary condition inherits from MarkAsPeriodic, which constitutes as it being marked as a...
Tensor.hpp
db::wrap_tags_in
tmpl::transform< TagList, tmpl::bind< Wrapper, tmpl::_1, tmpl::pin< Args >... > > wrap_tags_in
Create a new tmpl::list of tags by wrapping each tag in TagList in Wrapper<_, Args....
Definition: PrefixHelpers.hpp:30
magnitude
Scalar< DataType > magnitude(const Tensor< DataType, Symmetry< 1 >, index_list< Index >> &vector) noexcept
Compute the Euclidean magnitude of a rank-1 tensor.
Definition: Magnitude.hpp:27
optional
std::smatch
evolution::dg::Tags::NormalCovector
The normal covector to the interface.
Definition: NormalVectorTags.hpp:24
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
std::unique_ptr< domain::BoundaryConditions::BoundaryCondition >
Prefixes.hpp
Index::product
constexpr size_t product() const noexcept
The product of the indices. If Dim = 0, the product is defined as 1.
Definition: Index.hpp:72
TMPL.hpp
TestHelpers::evolution::dg::Tags::PythonFunctionForErrorMessage
The name of the python function that returns the error message.
Definition: BoundaryConditions.hpp:70
Mesh.hpp
gsl::not_null
Require a pointer to not be a nullptr
Definition: ReadSpecPiecewisePolynomial.hpp:13
string