17 #include "DataStructures/DataBox/PrefixHelpers.hpp"
19 #include "DataStructures/DataVector.hpp"
21 #include "DataStructures/Tensor/EagerMath/Magnitude.hpp"
24 #include "Domain/BoundaryConditions/BoundaryCondition.hpp"
25 #include "Domain/BoundaryConditions/Periodic.hpp"
27 #include "Evolution/BoundaryConditions/Type.hpp"
28 #include "Evolution/DiscontinuousGalerkin/Actions/ComputeTimeDerivativeHelpers.hpp"
29 #include "Evolution/DiscontinuousGalerkin/Actions/NormalCovectorAndMagnitude.hpp"
32 #include "Framework/TestCreation.hpp"
35 #include "Helpers/Evolution/DiscontinuousGalerkin/NormalVectors.hpp"
38 #include "Utilities/NoSuchType.hpp"
40 #include "Utilities/TaggedTuple.hpp"
42 namespace TestHelpers::evolution::dg {
54 template <
typename Tag,
typename BoundaryCorrection = NoSuchType>
57 using boundary_correction = BoundaryCorrection;
68 template <
typename BoundaryCorrection = NoSuchType>
70 using boundary_correction = BoundaryCorrection;
76 template <
typename Tag>
84 template <
typename ConversionClassList,
typename... Args>
88 static_assert(
sizeof...(Args) > 0,
89 "Call to python which returns a Tensor of DataVectors must "
90 "pass at least one argument");
93 PyObject* python_module = PyImport_ImportModule(module_name.c_str());
94 if (python_module ==
nullptr) {
99 PyObject* func = PyObject_GetAttrString(python_module, function_name.c_str());
100 if (func ==
nullptr or not PyCallable_Check(func)) {
101 if (PyErr_Occurred()) {
107 const std::array<size_t,
sizeof...(Args)> arg_sizes{
108 {pypp::detail::ContainerPackAndUnpack<
110 const size_t npts = *std::max_element(arg_sizes.begin(), arg_sizes.end());
111 for (
size_t i = 0; i < arg_sizes.size(); ++i) {
112 if (
gsl::at(arg_sizes, i) != 1 and
gsl::at(arg_sizes, i) != npts) {
113 ERROR(
"Each argument must return size 1 or "
115 <<
" (the number of points in the DataVector), but argument number "
116 << i <<
" has size " <<
gsl::at(arg_sizes, i));
121 for (
size_t s = 0; s < npts; ++s) {
123 pypp::detail::ContainerPackAndUnpack<Args, ConversionClassList>::unpack(
126 pypp::detail::call_work<
typename pypp::detail::ContainerPackAndUnpack<
127 ReturnType, ConversionClassList>::unpacked_container>(python_module,
129 if (ret.has_value()) {
130 result = std::move(ret);
135 Py_DECREF(python_module);
139 template <
typename BoundaryCorrection,
typename... PythonFunctionNameTags>
140 const std::string& get_python_error_message_function(
142 python_boundary_condition_functions) {
144 tmpl::list_contains_v<tmpl::list<PythonFunctionNameTags...>,
145 Tags::PythonFunctionForErrorMessage<NoSuchType>>,
146 Tags::PythonFunctionForErrorMessage<NoSuchType>,
147 Tags::PythonFunctionForErrorMessage<BoundaryCorrection>>>(
148 python_boundary_condition_functions);
151 template <
typename Tag,
typename BoundaryCorrection,
152 typename... PythonFunctionNameTags>
155 python_boundary_condition_functions) {
157 tmpl::list_contains_v<tmpl::list<PythonFunctionNameTags...>,
158 Tags::PythonFunctionName<Tag, NoSuchType>>,
159 Tags::PythonFunctionName<Tag, NoSuchType>,
160 Tags::PythonFunctionName<Tag, BoundaryCorrection>>>(
161 python_boundary_condition_functions);
164 template <
typename BoundaryConditionHelper,
typename AllTagsOnFaceList,
165 typename... TagsFromFace,
typename... VolumeArgs>
166 void apply_boundary_condition_impl(
167 BoundaryConditionHelper& boundary_condition_helper,
168 const Variables<AllTagsOnFaceList>& fields_on_interior_face,
169 tmpl::list<TagsFromFace...> ,
170 const VolumeArgs&... volume_args) noexcept {
171 boundary_condition_helper(get<TagsFromFace>(fields_on_interior_face)...,
175 template <
typename System,
typename ConversionClassList,
176 typename BoundaryCorrection,
typename... PythonFunctionNameTags,
177 typename BoundaryCondition,
size_t FaceDim,
typename DbTagsList,
178 typename... RangeTags,
typename... EvolvedVariablesTags,
179 typename... BoundaryCorrectionPackagedDataInputTags,
180 typename... BoundaryConditionVolumeTags>
181 void test_boundary_condition_with_python_impl(
185 python_boundary_condition_functions,
186 const BoundaryCondition& boundary_condition,
188 const db::DataBox<DbTagsList>& box_of_volume_data,
190 const bool use_moving_mesh, tmpl::list<EvolvedVariablesTags...> ,
191 tmpl::list<BoundaryCorrectionPackagedDataInputTags...> ,
192 tmpl::list<BoundaryConditionVolumeTags...> ,
const double epsilon) {
194 CAPTURE(python_module);
195 CAPTURE(face_points);
196 CAPTURE(use_moving_mesh);
198 CAPTURE(pretty_type::short_name<BoundaryCorrection>());
199 const size_t number_of_points_on_face = face_points.
product();
201 using variables_tag =
typename System::variables_tag;
202 using variables_tags =
typename variables_tag::tags_list;
203 using flux_variables =
typename System::flux_variables;
206 constexpr
bool uses_ghost =
207 BoundaryCondition::bc_type ==
208 ::evolution::BoundaryConditions::Type::Ghost or
209 BoundaryCondition::bc_type ==
210 ::evolution::BoundaryConditions::Type::GhostAndTimeDerivative;
211 constexpr
bool uses_time_derivative_condition =
212 BoundaryCondition::bc_type ==
213 ::evolution::BoundaryConditions::Type::TimeDerivative or
214 BoundaryCondition::bc_type ==
215 ::evolution::BoundaryConditions::Type::GhostAndTimeDerivative;
218 using inverse_spatial_metric_list =
219 ::evolution::dg::Actions::detail::inverse_spatial_metric_tag<System>;
220 constexpr
bool has_inv_spatial_metric =
221 ::evolution::dg::Actions::detail::has_inverse_spatial_metric_tag_v<
225 using bcondition_interior_temp_tags =
226 typename BoundaryCondition::dg_interior_temporary_tags;
227 using bcondition_interior_prim_tags =
228 typename ::evolution::dg::Actions::detail::get_primitive_vars<
229 System::has_primitive_and_conservative_vars>::
230 template boundary_condition_interior_tags<BoundaryCondition>;
231 using bcondition_interior_evolved_vars_tags =
232 typename BoundaryCondition::dg_interior_evolved_variables_tags;
233 using bcondition_interior_dt_evolved_vars_tags =
234 ::evolution::dg::Actions::detail::get_dt_vars_from_boundary_condition<
236 using bcondition_interior_deriv_evolved_vars_tags =
237 ::evolution::dg::Actions::detail::get_deriv_vars_from_boundary_condition<
239 using bcondition_interior_tags = tmpl::remove_duplicates<tmpl::append<
240 tmpl::conditional_t<has_inv_spatial_metric,
241 tmpl::list<::evolution::dg::Actions::detail::
242 NormalVector<FaceDim + 1>>,
244 bcondition_interior_evolved_vars_tags, bcondition_interior_prim_tags,
245 bcondition_interior_temp_tags, bcondition_interior_dt_evolved_vars_tags,
246 bcondition_interior_deriv_evolved_vars_tags>>;
250 Variables<tmpl::remove_duplicates<
251 tmpl::append<bcondition_interior_tags, inverse_spatial_metric_list>>>
252 interior_face_fields{number_of_points_on_face};
256 tmpl::for_each<tmpl::list<RangeTags...>>([&generator, &interior_face_fields,
257 &ranges](
auto tag_v) {
258 using tag = tmpl::type_from<decltype(tag_v)>;
265 auto interior_normal_covector =
267 tnsr::i<DataVector, FaceDim + 1, Frame::Inertial>>(
269 if constexpr (has_inv_spatial_metric) {
270 auto& inv_spatial_metric =
271 get<tmpl::front<inverse_spatial_metric_list>>(interior_face_fields);
272 detail::adjust_inverse_spatial_metric(
make_not_null(&inv_spatial_metric));
273 auto& normal_vector =
274 get<::evolution::dg::Actions::detail::NormalVector<FaceDim + 1>>(
275 interior_face_fields);
276 detail::normalize_vector_and_covector(
281 for (
DataVector& component : interior_normal_covector) {
289 if (use_moving_mesh) {
292 make_with_random_values<tnsr::I<DataVector, FaceDim + 1>>(
293 generator,
make_not_null(&local_dist), number_of_points_on_face);
296 if constexpr (BoundaryCondition::bc_type ==
297 ::evolution::BoundaryConditions::Type::Outflow) {
301 const auto apply_bc =
302 [&boundary_condition, &face_mesh_velocity, &interior_normal_covector,
303 &python_boundary_condition_functions,
304 &python_module](
const auto&... face_and_volume_args) noexcept {
306 boundary_condition.dg_outflow(face_mesh_velocity,
307 interior_normal_covector,
308 face_and_volume_args...);
310 get_python_error_message_function<BoundaryCorrection>(
311 python_boundary_condition_functions);
312 const auto python_error_message =
313 call_for_error_message<ConversionClassList>(
314 python_module, python_error_msg_function, face_mesh_velocity,
315 interior_normal_covector, face_and_volume_args...);
316 CAPTURE(python_error_msg_function);
317 CAPTURE(python_error_message.value_or(
""));
318 CAPTURE(error_msg.value_or(
""));
319 REQUIRE(python_error_message.has_value() == error_msg.has_value());
320 if (python_error_message.has_value() and error_msg.has_value()) {
322 CHECK(std::regex_search(*error_msg, matcher,
326 apply_boundary_condition_impl(
327 apply_bc, interior_face_fields, bcondition_interior_tags{},
328 db::get<BoundaryConditionVolumeTags>(box_of_volume_data)...);
331 if constexpr (uses_time_derivative_condition) {
332 Variables<dt_variables_tags> time_derivative_correction{
333 number_of_points_on_face};
334 auto apply_bc = [&boundary_condition, epsilon, &face_mesh_velocity,
335 &interior_normal_covector,
336 &python_boundary_condition_functions, &python_module,
337 &time_derivative_correction](
338 const auto&... interior_face_and_volume_args) {
340 boundary_condition.dg_time_derivative(
342 time_derivative_correction))...,
343 face_mesh_velocity, interior_normal_covector,
344 interior_face_and_volume_args...);
347 get_python_error_message_function<BoundaryCorrection>(
348 python_boundary_condition_functions);
349 const auto python_error_message =
350 call_for_error_message<ConversionClassList>(
351 python_module, python_error_msg_function, face_mesh_velocity,
352 interior_normal_covector, interior_face_and_volume_args...);
353 CAPTURE(python_error_msg_function);
354 CAPTURE(python_error_message.value_or(
""));
355 CAPTURE(error_msg.value_or(
""));
356 REQUIRE(python_error_message.has_value() == error_msg.has_value());
357 if (python_error_message.has_value() and error_msg.has_value()) {
359 CHECK(std::regex_search(*error_msg, matcher,
364 tmpl::for_each<dt_variables_tags>([&](
auto dt_var_tag_v) {
365 using DtVarTag = tmpl::type_from<decltype(dt_var_tag_v)>;
369 get_python_tag_function<DtVarTag, NoSuchType>(
370 python_boundary_condition_functions);
371 CAPTURE(python_tag_function);
372 CAPTURE(pretty_type::short_name<DtVarTag>());
373 typename DtVarTag::type python_result{};
376 pypp::call<typename DtVarTag::type, ConversionClassList>(
377 python_module, python_tag_function, face_mesh_velocity,
378 interior_normal_covector, interior_face_and_volume_args...);
380 INFO(
"Failed python call with '" << e.what() <<
"'");
385 get<DtVarTag>(time_derivative_correction), python_result,
386 Approx::custom().epsilon(epsilon).scale(1.0));
389 apply_boundary_condition_impl(
390 apply_bc, interior_face_fields, bcondition_interior_tags{},
391 db::get<BoundaryConditionVolumeTags>(box_of_volume_data)...);
394 if constexpr (uses_ghost) {
398 using correction_temp_tags =
399 typename BoundaryCorrection::dg_package_data_temporary_tags;
400 using correction_prim_tags =
401 typename ::evolution::dg::Actions::detail::get_primitive_vars<
402 System::has_primitive_and_conservative_vars>::
403 template f<BoundaryCorrection>;
404 using tags_on_exterior_face = tmpl::append<
405 variables_tags, fluxes_tags, correction_temp_tags, correction_prim_tags,
406 inverse_spatial_metric_list,
408 ::evolution::dg::Actions::detail::OneOverNormalVectorMagnitude,
409 ::evolution::dg::Actions::detail::NormalVector<FaceDim + 1>,
411 Variables<tags_on_exterior_face> exterior_face_fields{
412 number_of_points_on_face};
413 auto apply_bc = [&boundary_condition, epsilon, &exterior_face_fields,
414 &face_mesh_velocity, &interior_normal_covector,
415 &python_boundary_condition_functions, &python_module](
416 const auto&... interior_face_and_volume_args) {
419 exterior_face_fields))...,
420 face_mesh_velocity, interior_normal_covector,
421 interior_face_and_volume_args...);
424 get_python_error_message_function<BoundaryCorrection>(
425 python_boundary_condition_functions);
426 const auto python_error_message =
427 call_for_error_message<ConversionClassList>(
428 python_module, python_error_msg_function, face_mesh_velocity,
429 interior_normal_covector, interior_face_and_volume_args...);
430 CAPTURE(python_error_msg_function);
431 CAPTURE(python_error_message.value_or(
""));
432 CAPTURE(error_msg.value_or(
""));
433 REQUIRE(python_error_message.has_value() == error_msg.has_value());
434 if (python_error_message.has_value() and error_msg.has_value()) {
436 CHECK(std::regex_search(*error_msg, matcher,
439 if (python_error_message.has_value()) {
444 tmpl::for_each<tmpl::list<BoundaryCorrectionPackagedDataInputTags...>>(
445 [&](
auto boundary_correction_tag_v) {
446 using BoundaryCorrectionTag =
447 tmpl::type_from<decltype(boundary_correction_tag_v)>;
449 get_python_tag_function<BoundaryCorrectionTag,
451 python_boundary_condition_functions);
452 CAPTURE(python_tag_function);
453 CAPTURE(pretty_type::short_name<BoundaryCorrectionTag>());
454 typename BoundaryCorrectionTag::type python_result{};
456 python_result =
pypp::call<
typename BoundaryCorrectionTag::type,
457 ConversionClassList>(
458 python_module, python_tag_function, face_mesh_velocity,
459 interior_normal_covector, interior_face_and_volume_args...);
461 INFO(
"Failed python call with '" << e.what() <<
"'");
466 get<BoundaryCorrectionTag>(exterior_face_fields), python_result,
467 Approx::custom().epsilon(epsilon).scale(1.0));
473 apply_boundary_condition_impl(
474 apply_bc, interior_face_fields, bcondition_interior_tags{},
475 db::get<BoundaryConditionVolumeTags>(box_of_volume_data)...);
479 template <
typename T,
typename = std::
void_t<>>
480 struct get_boundary_conditions_impl {
481 using type = tmpl::list<>;
484 template <
typename T>
485 struct get_boundary_conditions_impl<
486 T, std::void_t<typename T::boundary_conditions_base>> {
487 using type =
typename T::boundary_conditions_base::creatable_classes;
490 template <
typename T>
491 using get_boundary_conditions =
typename get_boundary_conditions_impl<T>::type;
546 template <
typename BoundaryCondition,
typename BoundaryConditionBase,
547 typename System,
typename BoundaryCorrectionsList,
548 typename ConversionClassList = tmpl::list<>,
size_t FaceDim,
549 typename DbTagsList,
typename... RangeTags,
550 typename... PythonFunctionNameTags>
551 void test_boundary_condition_with_python(
555 python_boundary_condition_functions,
557 const db::DataBox<DbTagsList>& box_of_volume_data,
559 const double epsilon = 1.0e-12) {
560 PUPable_reg(BoundaryCondition);
562 "All boundary condition classes must be marked `final`.");
563 using variables_tags =
typename System::variables_tag::tags_list;
564 using flux_variables =
typename System::flux_variables;
570 TestHelpers::test_factory_creation<BoundaryConditionBase>(
576 tmpl::for_each<BoundaryCorrectionsList>(
577 [&boundary_condition, &box_of_volume_data, epsilon, &face_points,
578 &generator, &python_boundary_condition_functions, &python_module,
579 &ranges](
auto boundary_correction_v) {
580 using BoundaryCorrection =
581 tmpl::type_from<decltype(boundary_correction_v)>;
582 using package_data_input_tags = tmpl::append<
583 variables_tags, fluxes_tags,
584 typename BoundaryCorrection::dg_package_data_temporary_tags,
585 typename ::evolution::dg::Actions::detail::get_primitive_vars<
586 System::has_primitive_and_conservative_vars>::
587 template f<BoundaryCorrection>>;
588 using boundary_condition_dg_gridless_tags =
589 typename BoundaryCondition::dg_gridless_tags;
590 for (
const auto use_moving_mesh : {
false,
true}) {
591 detail::test_boundary_condition_with_python_impl<
592 System, ConversionClassList, BoundaryCorrection>(
593 generator, python_module, python_boundary_condition_functions,
594 dynamic_cast<const BoundaryCondition&
>(*boundary_condition),
595 face_points, box_of_volume_data, ranges, use_moving_mesh,
596 variables_tags{}, package_data_input_tags{},
597 boundary_condition_dg_gridless_tags{}, epsilon);
599 INFO(
"Test boundary condition after serialization.");
600 const auto deserialized_bc =
602 detail::test_boundary_condition_with_python_impl<
603 System, ConversionClassList, BoundaryCorrection>(
604 generator, python_module, python_boundary_condition_functions,
605 dynamic_cast<const BoundaryCondition&
>(*deserialized_bc),
606 face_points, box_of_volume_data, ranges, use_moving_mesh,
607 variables_tags{}, package_data_input_tags{},
608 boundary_condition_dg_gridless_tags{}, epsilon);
614 template <
typename BoundaryCondition,
typename BoundaryConditionBase>
615 void test_periodic_condition(
const std::string& factory_string) {
616 PUPable_reg(BoundaryCondition);
617 const auto boundary_condition =
618 TestHelpers::test_factory_creation<BoundaryConditionBase>(factory_string);
619 REQUIRE(
typeid(*boundary_condition.get()) ==
typeid(BoundaryCondition));
620 const auto bc_clone = boundary_condition->get_clone();