Line data Source code
1 0 : // Distributed under the MIT License.
2 : // See LICENSE.txt for details.
3 :
4 : #pragma once
5 :
6 : #include <cstddef>
7 : #include <memory>
8 : #include <optional>
9 : #include <pup.h>
10 : #include <string>
11 : #include <utility>
12 :
13 : #include "DataStructures/DataBox/PrefixHelpers.hpp"
14 : #include "DataStructures/DataVector.hpp"
15 : #include "DataStructures/Tensor/Tensor.hpp"
16 : #include "Evolution/BoundaryConditions/Type.hpp"
17 : #include "Evolution/DiscontinuousGalerkin/Actions/ComputeTimeDerivativeHelpers.hpp"
18 : #include "Evolution/Systems/CurvedScalarWave/BoundaryConditions/Factory.hpp"
19 : #include "Evolution/Systems/CurvedScalarWave/System.hpp"
20 : #include "Evolution/Systems/GeneralizedHarmonic/BoundaryConditions/Factory.hpp"
21 : #include "Evolution/Systems/GeneralizedHarmonic/System.hpp"
22 : #include "Evolution/Systems/ScalarTensor/BoundaryConditions/BoundaryCondition.hpp"
23 : #include "Evolution/Systems/ScalarTensor/Tags.hpp"
24 : #include "NumericalAlgorithms/DiscontinuousGalerkin/Formulation.hpp"
25 : #include "Options/String.hpp"
26 : #include "PointwiseFunctions/GeneralRelativity/Lapse.hpp"
27 : #include "PointwiseFunctions/GeneralRelativity/Shift.hpp"
28 : #include "PointwiseFunctions/GeneralRelativity/SpatialMetric.hpp"
29 : #include "Utilities/Gsl.hpp"
30 : #include "Utilities/Serialization/CharmPupable.hpp"
31 : #include "Utilities/TMPL.hpp"
32 : #include "Utilities/TaggedTuple.hpp"
33 :
34 : namespace ScalarTensor::BoundaryConditions {
35 : namespace detail {
36 :
37 : template <evolution::BoundaryConditions::Type GhBcType,
38 : evolution::BoundaryConditions::Type ScalarBcType>
39 : struct UnionOfBcTypes {
40 : static_assert(
41 : (GhBcType == evolution::BoundaryConditions::Type::TimeDerivative) or
42 : (ScalarBcType ==
43 : evolution::BoundaryConditions::Type::TimeDerivative) or
44 : (GhBcType ==
45 : evolution::BoundaryConditions::Type::GhostAndTimeDerivative) or
46 : (ScalarBcType ==
47 : evolution::BoundaryConditions::Type::GhostAndTimeDerivative),
48 : "TimeDerivative boundary conditions are not yet supported.");
49 : static constexpr evolution::BoundaryConditions::Type bc_type =
50 : evolution::BoundaryConditions::Type::Ghost;
51 : };
52 :
53 : template <>
54 : struct UnionOfBcTypes<evolution::BoundaryConditions::Type::Ghost,
55 : evolution::BoundaryConditions::Type::Ghost> {
56 : static constexpr evolution::BoundaryConditions::Type bc_type =
57 : evolution::BoundaryConditions::Type::Ghost;
58 : };
59 :
60 : template <evolution::BoundaryConditions::Type GhBcType>
61 : struct UnionOfBcTypes<
62 : GhBcType, evolution::BoundaryConditions::Type::DemandOutgoingCharSpeeds> {
63 : static_assert(
64 : GhBcType == evolution::BoundaryConditions::Type::DemandOutgoingCharSpeeds,
65 : "If either boundary condition in `ProductOfConditions` has "
66 : "`Type::DemandOutgoingCharSpeeds`, both must have "
67 : "`Type::DemandOutgoingCharSpeeds`");
68 : };
69 :
70 : template <evolution::BoundaryConditions::Type ScalarBcType>
71 : struct UnionOfBcTypes<
72 : evolution::BoundaryConditions::Type::DemandOutgoingCharSpeeds,
73 : ScalarBcType> {
74 : static_assert(
75 : ScalarBcType ==
76 : evolution::BoundaryConditions::Type::DemandOutgoingCharSpeeds,
77 : "If either boundary condition in `ProductOfConditions` has "
78 : "`Type::DemandOutgoingCharSpeeds`, both must have "
79 : "`Type::DemandOutgoingCharSpeeds`");
80 : };
81 :
82 : template <>
83 : struct UnionOfBcTypes<
84 : evolution::BoundaryConditions::Type::DemandOutgoingCharSpeeds,
85 : evolution::BoundaryConditions::Type::DemandOutgoingCharSpeeds> {
86 : static constexpr evolution::BoundaryConditions::Type bc_type =
87 : evolution::BoundaryConditions::Type::DemandOutgoingCharSpeeds;
88 : };
89 :
90 : } // namespace detail
91 :
92 : /*!
93 : * \brief Apply a boundary condition to the combined Generalized Harmonic (GH)
94 : * and CurvedScalarWave (CSW) system using the boundary conditions defined
95 : * separately for the GH and CSW systems.
96 : */
97 : template <typename DerivedGhCondition, typename DerivedScalarCondition>
98 1 : class ProductOfConditions final : public BoundaryCondition {
99 : public:
100 0 : static constexpr size_t dim = 3;
101 0 : static constexpr evolution::BoundaryConditions::Type bc_type =
102 : detail::UnionOfBcTypes<DerivedGhCondition::bc_type,
103 : DerivedScalarCondition::bc_type>::bc_type;
104 :
105 0 : using dg_interior_evolved_variables_tags =
106 : tmpl::remove_duplicates<tmpl::append<
107 : typename DerivedGhCondition::dg_interior_evolved_variables_tags,
108 : typename DerivedScalarCondition::dg_interior_evolved_variables_tags>>;
109 :
110 0 : using dg_interior_temporary_tags = tmpl::remove_duplicates<tmpl::append<
111 : typename DerivedGhCondition::dg_interior_temporary_tags,
112 : typename DerivedScalarCondition::dg_interior_temporary_tags>>;
113 :
114 0 : using dg_gridless_tags = tmpl::remove_duplicates<
115 : tmpl::append<typename DerivedGhCondition::dg_gridless_tags,
116 : typename DerivedScalarCondition::dg_gridless_tags>>;
117 :
118 0 : using dg_interior_dt_vars_tags = tmpl::append<
119 : evolution::dg::Actions::detail::get_dt_vars_from_boundary_condition<
120 : DerivedGhCondition>,
121 : evolution::dg::Actions::detail::get_dt_vars_from_boundary_condition<
122 : DerivedScalarCondition>>;
123 :
124 0 : using dg_interior_deriv_vars_tags = tmpl::append<
125 : evolution::dg::Actions::detail::get_deriv_vars_from_boundary_condition<
126 : DerivedGhCondition>,
127 : evolution::dg::Actions::detail::get_deriv_vars_from_boundary_condition<
128 : DerivedScalarCondition>>;
129 :
130 0 : static std::string name() {
131 : return "Product" + pretty_type::name<DerivedGhCondition>() + "And" +
132 : pretty_type::name<DerivedScalarCondition>();
133 : }
134 :
135 0 : struct GhCondition {
136 0 : using type = DerivedGhCondition;
137 0 : static std::string name() {
138 : return "GeneralizedHarmonic" + pretty_type::name<DerivedGhCondition>();
139 : }
140 0 : static constexpr Options::String help{
141 : "The Generalized Harmonic part of the product boundary condition"};
142 : };
143 :
144 0 : struct ScalarCondition {
145 0 : using type = DerivedScalarCondition;
146 0 : static std::string name() {
147 : return "Scalar" + pretty_type::name<DerivedScalarCondition>();
148 : }
149 0 : static constexpr Options::String help{
150 : "The Scalar part of the product boundary condition"};
151 : };
152 :
153 0 : using options = tmpl::list<GhCondition, ScalarCondition>;
154 :
155 0 : static constexpr Options::String help = {
156 : "Direct product of a GH and CurvedScalarWave boundary conditions. "
157 : "See the documentation for the two individual boundary conditions for "
158 : "further details."};
159 :
160 0 : ProductOfConditions() = default;
161 0 : ProductOfConditions(DerivedGhCondition gh_condition,
162 : DerivedScalarCondition scalar_condition)
163 : : derived_gh_condition_{gh_condition},
164 : derived_scalar_condition_{scalar_condition} {}
165 0 : ProductOfConditions(const ProductOfConditions&) = default;
166 0 : ProductOfConditions& operator=(const ProductOfConditions&) = default;
167 0 : ProductOfConditions(ProductOfConditions&&) = default;
168 0 : ProductOfConditions& operator=(ProductOfConditions&&) = default;
169 0 : ~ProductOfConditions() override = default;
170 :
171 : /// \cond
172 : explicit ProductOfConditions(CkMigrateMessage* msg)
173 : : BoundaryCondition(msg) {}
174 : using PUP::able::register_constructor;
175 : WRAPPED_PUPable_decl_base_template(
176 : domain::BoundaryConditions::BoundaryCondition, ProductOfConditions);
177 : /// \endcond
178 :
179 0 : void pup(PUP::er& p) override;
180 :
181 0 : auto get_clone() const -> std::unique_ptr<
182 : domain::BoundaryConditions::BoundaryCondition> override;
183 :
184 0 : std::optional<std::string> dg_demand_outgoing_char_speeds(
185 : // GH arguments
186 : const std::optional<tnsr::I<DataVector, dim, Frame::Inertial>>&
187 : face_mesh_velocity,
188 : const tnsr::i<DataVector, dim, Frame::Inertial>&
189 : outward_directed_normal_covector,
190 : const tnsr::I<DataVector, dim, Frame::Inertial>&
191 : outward_directed_normal_vector,
192 : const Scalar<DataVector>& gamma_1, const Scalar<DataVector>& lapse,
193 : const tnsr::I<DataVector, dim, Frame::Inertial>& shift,
194 : // Scalar arguments
195 : const Scalar<DataVector>& gamma1_scalar) const {
196 : // DemandOutgoingCharSpeeds condition is only valid if both boundary
197 : // conditions are DemandOutgoingCharSpeeds, so we directly apply both. A
198 : // static_assert elsewhere is triggered if only one boundary condition is
199 : // DemandOutgoingCharSpeeds.
200 : auto gh_string = derived_gh_condition_.dg_demand_outgoing_char_speeds(
201 : face_mesh_velocity, outward_directed_normal_covector,
202 : outward_directed_normal_vector, gamma_1, lapse, shift);
203 :
204 : auto scalar_string =
205 : derived_scalar_condition_.dg_demand_outgoing_char_speeds(
206 : face_mesh_velocity, outward_directed_normal_covector,
207 : outward_directed_normal_vector, gamma1_scalar, lapse, shift);
208 :
209 : if (not gh_string.has_value()) {
210 : return scalar_string;
211 : }
212 : if (not scalar_string.has_value()) {
213 : return gh_string;
214 : }
215 : return gh_string.value() + ";" + scalar_string.value();
216 : }
217 :
218 : // We overload dg_ghost for the different analytic boundary conditions
219 :
220 : // Boundary conditions for DirichletAnalytic and AnalyticConstant
221 0 : std::optional<std::string> dg_ghost(
222 : // GH evolved variables
223 : const gsl::not_null<tnsr::aa<DataVector, dim, Frame::Inertial>*>
224 : spacetime_metric,
225 : const gsl::not_null<tnsr::aa<DataVector, dim, Frame::Inertial>*> pi,
226 : const gsl::not_null<tnsr::iaa<DataVector, dim, Frame::Inertial>*> phi,
227 : // Scalar evolved variables
228 : const gsl::not_null<Scalar<DataVector>*> psi_scalar,
229 : const gsl::not_null<Scalar<DataVector>*> pi_scalar,
230 : const gsl::not_null<tnsr::i<DataVector, dim, Frame::Inertial>*>
231 : phi_scalar,
232 : // GH temporary variables
233 : const gsl::not_null<Scalar<DataVector>*> gamma1,
234 : const gsl::not_null<Scalar<DataVector>*> gamma2,
235 : const gsl::not_null<Scalar<DataVector>*> lapse,
236 : const gsl::not_null<tnsr::I<DataVector, dim, Frame::Inertial>*> shift,
237 : // Scalar temporary variables
238 : const gsl::not_null<Scalar<DataVector>*> gamma1_scalar,
239 : const gsl::not_null<Scalar<DataVector>*> gamma2_scalar,
240 : // Inverse metric
241 : const gsl::not_null<tnsr::II<DataVector, dim, Frame::Inertial>*>
242 : inv_spatial_metric,
243 : // Mesh variables
244 : const std::optional<tnsr::I<DataVector, dim, Frame::Inertial>>&
245 : face_mesh_velocity,
246 : const tnsr::i<DataVector, dim, Frame::Inertial>& normal_covector,
247 : const tnsr::I<DataVector, dim, Frame::Inertial>& normal_vector,
248 : // GH interior variables
249 : const tnsr::I<DataVector, dim, Frame::Inertial>& coords,
250 : const Scalar<DataVector>& interior_gamma1,
251 : const Scalar<DataVector>& interior_gamma2,
252 : // Scalar interior variables
253 : const tnsr::II<DataVector, dim, Frame::Inertial>&
254 : inverse_spatial_metric_interior,
255 : const Scalar<DataVector>& gamma1_interior_scalar,
256 : const Scalar<DataVector>& gamma2_interior_scalar,
257 : const Scalar<DataVector>& lapse_interior,
258 : const tnsr::I<DataVector, dim>& shift_interior, const double time) const {
259 : // For gh::BoundaryConditions::DirichletAnalytic
260 : auto gh_string = derived_gh_condition_.dg_ghost(
261 : spacetime_metric, pi, phi, gamma1, gamma2, lapse, shift,
262 : inv_spatial_metric, face_mesh_velocity, normal_covector, normal_vector,
263 : coords, interior_gamma1, interior_gamma2, time);
264 :
265 : // For CurvedScalarWave::BoundaryConditions::AnalyticConstant
266 : auto scalar_string = derived_scalar_condition_.dg_ghost(
267 : psi_scalar, pi_scalar, phi_scalar, lapse, shift, gamma1_scalar,
268 : gamma2_scalar, inv_spatial_metric, face_mesh_velocity, normal_covector,
269 : normal_vector, inverse_spatial_metric_interior, gamma1_interior_scalar,
270 : gamma2_interior_scalar, lapse_interior, shift_interior);
271 :
272 : if (not gh_string.has_value()) {
273 : return scalar_string;
274 : }
275 : if (not scalar_string.has_value()) {
276 : return gh_string;
277 : }
278 : return gh_string.value() + ";" + scalar_string.value();
279 : }
280 :
281 : // Boundary conditions for DirichletMinkowski and AnalyticConstant
282 0 : std::optional<std::string> dg_ghost(
283 : // GH evolved variables
284 : const gsl::not_null<tnsr::aa<DataVector, dim, Frame::Inertial>*>
285 : spacetime_metric,
286 : const gsl::not_null<tnsr::aa<DataVector, dim, Frame::Inertial>*> pi,
287 : const gsl::not_null<tnsr::iaa<DataVector, dim, Frame::Inertial>*> phi,
288 : // Scalar evolved variables
289 : const gsl::not_null<Scalar<DataVector>*> psi_scalar,
290 : const gsl::not_null<Scalar<DataVector>*> pi_scalar,
291 : const gsl::not_null<tnsr::i<DataVector, dim, Frame::Inertial>*>
292 : phi_scalar,
293 : // GH temporary variables
294 : const gsl::not_null<Scalar<DataVector>*> gamma1,
295 : const gsl::not_null<Scalar<DataVector>*> gamma2,
296 : const gsl::not_null<Scalar<DataVector>*> lapse,
297 : const gsl::not_null<tnsr::I<DataVector, dim, Frame::Inertial>*> shift,
298 : // Scalar temporary variables
299 : const gsl::not_null<Scalar<DataVector>*> gamma1_scalar,
300 : const gsl::not_null<Scalar<DataVector>*> gamma2_scalar,
301 : // Inverse metric
302 : const gsl::not_null<tnsr::II<DataVector, dim, Frame::Inertial>*>
303 : inv_spatial_metric,
304 : // Mesh variables
305 : const std::optional<tnsr::I<DataVector, dim, Frame::Inertial>>&
306 : face_mesh_velocity,
307 : const tnsr::i<DataVector, dim, Frame::Inertial>& normal_covector,
308 : const tnsr::I<DataVector, dim, Frame::Inertial>& normal_vector,
309 : // GH interior variables
310 : const Scalar<DataVector>& interior_gamma1,
311 : const Scalar<DataVector>& interior_gamma2,
312 : // Scalar interior variables
313 : const tnsr::II<DataVector, dim, Frame::Inertial>&
314 : inverse_spatial_metric_interior,
315 : const Scalar<DataVector>& gamma1_interior_scalar,
316 : const Scalar<DataVector>& gamma2_interior_scalar,
317 : const Scalar<DataVector>& lapse_interior,
318 : const tnsr::I<DataVector, dim>& shift_interior) const {
319 : // For gh::BoundaryConditions::DirichletMinkowski
320 : auto gh_string = derived_gh_condition_.dg_ghost(
321 : spacetime_metric, pi, phi, gamma1, gamma2, lapse, shift,
322 : inv_spatial_metric, face_mesh_velocity, normal_covector, normal_vector,
323 : interior_gamma1, interior_gamma2);
324 :
325 : // For CurvedScalarWave::BoundaryConditions::AnalyticConstant
326 : auto scalar_string = derived_scalar_condition_.dg_ghost(
327 : psi_scalar, pi_scalar, phi_scalar, lapse, shift, gamma1_scalar,
328 : gamma2_scalar, inv_spatial_metric, face_mesh_velocity, normal_covector,
329 : normal_vector, inverse_spatial_metric_interior, gamma1_interior_scalar,
330 : gamma2_interior_scalar, lapse_interior, shift_interior);
331 :
332 : if (not gh_string.has_value()) {
333 : return scalar_string;
334 : }
335 : if (not scalar_string.has_value()) {
336 : return gh_string;
337 : }
338 : return gh_string.value() + ";" + scalar_string.value();
339 : }
340 :
341 : private:
342 0 : DerivedGhCondition derived_gh_condition_;
343 0 : DerivedScalarCondition derived_scalar_condition_;
344 : };
345 :
346 : template <typename DerivedGhCondition, typename DerivedScalarCondition>
347 : void ProductOfConditions<DerivedGhCondition, DerivedScalarCondition>::pup(
348 : PUP::er& p) {
349 : BoundaryCondition::pup(p);
350 : p | derived_gh_condition_;
351 : p | derived_scalar_condition_;
352 : }
353 :
354 : template <typename DerivedGhCondition, typename DerivedScalarCondition>
355 : auto ProductOfConditions<DerivedGhCondition,
356 : DerivedScalarCondition>::get_clone() const
357 : -> std::unique_ptr<domain::BoundaryConditions::BoundaryCondition> {
358 : return std::make_unique<ProductOfConditions>(*this);
359 : }
360 :
361 : /// \cond
362 : template <typename DerivedGhCondition, typename DerivedScalarCondition>
363 : PUP::able::PUP_ID ProductOfConditions<DerivedGhCondition,
364 : DerivedScalarCondition>::my_PUP_ID =
365 : 0; // NOLINT
366 : /// \endcond
367 :
368 : template <typename DerivedGhCondition, typename DerivedScalarCondition>
369 0 : ProductOfConditions(DerivedGhCondition gh_condition,
370 : DerivedScalarCondition scalar_condition)
371 : -> ProductOfConditions<DerivedGhCondition, DerivedScalarCondition>;
372 :
373 : } // namespace ScalarTensor::BoundaryConditions
|
