Line data Source code
1 0 : // Distributed under the MIT License.
2 : // See LICENSE.txt for details.
3 :
4 : #pragma once
5 :
6 : #include <array>
7 : #include <cstddef>
8 : #include <memory>
9 : #include <optional>
10 : #include <pup.h>
11 : #include <string>
12 : #include <type_traits>
13 :
14 : #include "DataStructures/DataBox/Prefixes.hpp"
15 : #include "DataStructures/DataVector.hpp"
16 : #include "DataStructures/Tensor/Tensor.hpp"
17 : #include "DataStructures/Variables.hpp"
18 : #include "Evolution/BoundaryConditions/Type.hpp"
19 : #include "Evolution/Systems/GeneralizedHarmonic/BoundaryConditions/BoundaryCondition.hpp"
20 : #include "Evolution/Systems/GeneralizedHarmonic/ConstraintDamping/Tags.hpp"
21 : #include "Evolution/Systems/GeneralizedHarmonic/Tags.hpp"
22 : #include "Options/Options.hpp"
23 : #include "Options/String.hpp"
24 : #include "PointwiseFunctions/AnalyticData/Tags.hpp"
25 : #include "PointwiseFunctions/AnalyticSolutions/AnalyticSolution.hpp"
26 : #include "PointwiseFunctions/GeneralRelativity/Tags.hpp"
27 : #include "Utilities/Gsl.hpp"
28 : #include "Utilities/Serialization/CharmPupable.hpp"
29 : #include "Utilities/TMPL.hpp"
30 :
31 : /// \cond
32 : namespace domain::Tags {
33 : template <size_t Dim, typename Frame>
34 : struct Coordinates;
35 : } // namespace domain::Tags
36 : /// \endcond
37 :
38 1 : namespace gh::BoundaryConditions::detail {
39 : enum class ConstraintPreservingBjorhusType {
40 : ConstraintPreserving,
41 : ConstraintPreservingPhysical
42 : };
43 :
44 : ConstraintPreservingBjorhusType
45 : convert_constraint_preserving_bjorhus_type_from_yaml(
46 : const Options::Option& options);
47 : } // namespace gh::BoundaryConditions::detail
48 :
49 : namespace gh::BoundaryConditions {
50 : /*!
51 : * \brief Sets constraint preserving boundary conditions using the Bjorhus
52 : * method.
53 : *
54 : * \details Boundary conditions for the generalized harmonic evolution system
55 : * can be divided in to three parts, constraint-preserving, physical and gauge
56 : * boundary conditions.
57 : *
58 : * The generalized harmonic (GH) evolution system is a first-order reduction of
59 : * Einstein equations brought about by the imposition of GH gauge. This
60 : * introduces constraints on the free (evolved) variables in addition to the
61 : * standard Hamiltonian and momentum constraints. The constraint-preserving
62 : * portion of the boundary conditions is designed to prevent the influx of
63 : * constraint violations from external faces of the evolution domain, by damping
64 : * them away on a controlled and short time-scale. These conditions are imposed
65 : * as corrections to the characteristic projections of the right-hand-sides of
66 : * the GH evolution equations (i.e. using Bjorhus' method \cite Bjorhus1995),
67 : * as written down in Eq. (63) - (65) of \cite Lindblom2005qh . In addition to
68 : * these equations, the fourth projection is simply frozen in the unlikely case
69 : * its coordinate speed becomes negative, i.e. \f$d_t u^{\hat{1}+}{}_{ab}=0\f$
70 : * (in the notation of \cite Lindblom2005qh). The gauge degrees
71 : * of freedom are controlled by imposing a Sommerfeld-type condition (\f$L=0\f$
72 : * member of the hierarchy derived in \cite BaylissTurkel) that allow gauge
73 : * perturbations to pass through the boundary without strong reflections. These
74 : * assume a spherical outer boundary, and can be written down as in Eq. (25) of
75 : * \cite Rinne2007ui . Finally, the physical boundary conditions control the
76 : * influx of inward propagating gravitational-wave solutions from the external
77 : * boundaries. These are derived by considering the evolution system of the Weyl
78 : * curvature tensor, and controlling the inward propagating characteristics of
79 : * the system that are proportional to the Newman-Penrose curvature spinor
80 : * components \f$\Psi_4\f$ and \f$\Psi_0\f$. Here we use Eq. (68) of
81 : * \cite Lindblom2005qh to disallow any incoming waves. It is to be noted that
82 : * all the above conditions are also imposed on characteristic modes with speeds
83 : * exactly zero.
84 : *
85 : * This class provides two choices of combinations of the above corrections:
86 : * - `ConstraintPreserving` : this imposes the constraint-preserving and
87 : * gauge-controlling corrections;
88 : * - `ConstraintPreservingPhysical` : this additionally restricts the influx of
89 : * any physical gravitational waves from the outer boundary, in addition to
90 : * preventing the influx of constraint violations and gauge perturbations.
91 : *
92 : * We refer to `Bjorhus::constraint_preserving_bjorhus_corrections_dt_v_psi()`,
93 : * `Bjorhus::constraint_preserving_bjorhus_corrections_dt_v_zero()`,
94 : * `Bjorhus::constraint_preserving_bjorhus_corrections_dt_v_minus()`, and
95 : * `Bjorhus::constraint_preserving_physical_bjorhus_corrections_dt_v_minus()`
96 : * for the further details on implementation.
97 : *
98 : * \note These boundary conditions assume a spherical outer boundary. Also, we
99 : * do not yet have an option to inject incoming gravitational waves at the outer
100 : * boundary.
101 : */
102 : template <size_t Dim>
103 1 : class ConstraintPreservingBjorhus final : public BoundaryCondition<Dim> {
104 : public:
105 0 : struct TypeOptionTag {
106 0 : using type = detail::ConstraintPreservingBjorhusType;
107 0 : static std::string name() { return "Type"; }
108 0 : static constexpr Options::String help{
109 : "Whether to impose ConstraintPreserving, with or without physical "
110 : "terms for VMinus."};
111 : };
112 :
113 0 : using options = tmpl::list<TypeOptionTag>;
114 0 : static constexpr Options::String help{
115 : "ConstraintPreservingBjorhus boundary conditions setting the value of the"
116 : "time derivatives of the spacetime metric, Phi and Pi to expressions that"
117 : "prevent the influx of constraint violations and reflections."};
118 0 : static std::string name() { return "ConstraintPreservingBjorhus"; }
119 :
120 0 : explicit ConstraintPreservingBjorhus(
121 : detail::ConstraintPreservingBjorhusType type);
122 :
123 0 : ConstraintPreservingBjorhus() = default;
124 : /// \cond
125 : ConstraintPreservingBjorhus(ConstraintPreservingBjorhus&&) = default;
126 : ConstraintPreservingBjorhus& operator=(ConstraintPreservingBjorhus&&) =
127 : default;
128 : ConstraintPreservingBjorhus(const ConstraintPreservingBjorhus&) = default;
129 : ConstraintPreservingBjorhus& operator=(const ConstraintPreservingBjorhus&) =
130 : default;
131 : /// \endcond
132 0 : ~ConstraintPreservingBjorhus() override = default;
133 :
134 0 : explicit ConstraintPreservingBjorhus(CkMigrateMessage* msg);
135 :
136 0 : WRAPPED_PUPable_decl_base_template(
137 : domain::BoundaryConditions::BoundaryCondition,
138 : ConstraintPreservingBjorhus);
139 :
140 0 : auto get_clone() const -> std::unique_ptr<
141 : domain::BoundaryConditions::BoundaryCondition> override;
142 :
143 0 : static constexpr evolution::BoundaryConditions::Type bc_type =
144 : evolution::BoundaryConditions::Type::TimeDerivative;
145 :
146 0 : void pup(PUP::er& p) override;
147 :
148 0 : using dg_interior_evolved_variables_tags =
149 : tmpl::list<gr::Tags::SpacetimeMetric<DataVector, Dim>,
150 : Tags::Pi<DataVector, Dim>, Tags::Phi<DataVector, Dim>>;
151 0 : using dg_interior_temporary_tags =
152 : tmpl::list<domain::Tags::Coordinates<Dim, Frame::Inertial>,
153 : ConstraintDamping::Tags::ConstraintGamma1,
154 : ConstraintDamping::Tags::ConstraintGamma2,
155 : gr::Tags::Lapse<DataVector>, gr::Tags::Shift<DataVector, Dim>,
156 : gr::Tags::InverseSpacetimeMetric<DataVector, Dim>,
157 : gr::Tags::SpacetimeNormalVector<DataVector, Dim>,
158 : Tags::ThreeIndexConstraint<DataVector, Dim>,
159 : Tags::GaugeH<DataVector, Dim>,
160 : Tags::SpacetimeDerivGaugeH<DataVector, Dim>>;
161 0 : using dg_interior_dt_vars_tags =
162 : tmpl::list<::Tags::dt<gr::Tags::SpacetimeMetric<DataVector, Dim>>,
163 : ::Tags::dt<Tags::Pi<DataVector, Dim>>,
164 : ::Tags::dt<Tags::Phi<DataVector, Dim>>>;
165 0 : using dg_interior_deriv_vars_tags =
166 : tmpl::list<::Tags::deriv<gr::Tags::SpacetimeMetric<DataVector, Dim>,
167 : tmpl::size_t<Dim>, Frame::Inertial>,
168 : ::Tags::deriv<Tags::Pi<DataVector, Dim>, tmpl::size_t<Dim>,
169 : Frame::Inertial>,
170 : ::Tags::deriv<Tags::Phi<DataVector, Dim>, tmpl::size_t<Dim>,
171 : Frame::Inertial>>;
172 0 : using dg_gridless_tags = tmpl::list<>;
173 :
174 0 : std::optional<std::string> dg_time_derivative(
175 : gsl::not_null<tnsr::aa<DataVector, Dim, Frame::Inertial>*>
176 : dt_spacetime_metric_correction,
177 : gsl::not_null<tnsr::aa<DataVector, Dim, Frame::Inertial>*>
178 : dt_pi_correction,
179 : gsl::not_null<tnsr::iaa<DataVector, Dim, Frame::Inertial>*>
180 : dt_phi_correction,
181 :
182 : const std::optional<tnsr::I<DataVector, Dim, Frame::Inertial>>&
183 : face_mesh_velocity,
184 : const tnsr::i<DataVector, Dim, Frame::Inertial>& normal_covector,
185 : const tnsr::I<DataVector, Dim, Frame::Inertial>& /*normal_vector*/,
186 : // c.f. dg_interior_evolved_variables_tags
187 : const tnsr::aa<DataVector, Dim, Frame::Inertial>& spacetime_metric,
188 : const tnsr::aa<DataVector, Dim, Frame::Inertial>& pi,
189 : const tnsr::iaa<DataVector, Dim, Frame::Inertial>& phi,
190 : // c.f. dg_interior_temporary_tags
191 : const tnsr::I<DataVector, Dim, Frame::Inertial>& coords,
192 : const Scalar<DataVector>& gamma1, const Scalar<DataVector>& gamma2,
193 : const Scalar<DataVector>& lapse,
194 : const tnsr::I<DataVector, Dim, Frame::Inertial>& shift,
195 : const tnsr::AA<DataVector, Dim, Frame::Inertial>&
196 : inverse_spacetime_metric,
197 : const tnsr::A<DataVector, Dim, Frame::Inertial>&
198 : spacetime_unit_normal_vector,
199 : const tnsr::iaa<DataVector, Dim, Frame::Inertial>& three_index_constraint,
200 : const tnsr::a<DataVector, Dim, Frame::Inertial>& gauge_source,
201 : const tnsr::ab<DataVector, Dim, Frame::Inertial>&
202 : spacetime_deriv_gauge_source,
203 : // c.f. dg_interior_dt_vars_tags
204 : const tnsr::aa<DataVector, Dim, Frame::Inertial>&
205 : logical_dt_spacetime_metric,
206 : const tnsr::aa<DataVector, Dim, Frame::Inertial>& logical_dt_pi,
207 : const tnsr::iaa<DataVector, Dim, Frame::Inertial>& logical_dt_phi,
208 : // c.f. dg_interior_deriv_vars_tags
209 : const tnsr::iaa<DataVector, Dim, Frame::Inertial>& d_spacetime_metric,
210 : const tnsr::iaa<DataVector, Dim, Frame::Inertial>& d_pi,
211 : const tnsr::ijaa<DataVector, Dim, Frame::Inertial>& d_phi) const;
212 :
213 : private:
214 0 : void compute_intermediate_vars(
215 : gsl::not_null<tnsr::I<DataVector, Dim, Frame::Inertial>*>
216 : unit_interface_normal_vector,
217 : gsl::not_null<tnsr::iaa<DataVector, Dim, Frame::Inertial>*>
218 : four_index_constraint,
219 : gsl::not_null<tnsr::II<DataVector, Dim, Frame::Inertial>*>
220 : inverse_spatial_metric,
221 : gsl::not_null<tnsr::ii<DataVector, Dim, Frame::Inertial>*>
222 : extrinsic_curvature,
223 : gsl::not_null<tnsr::a<DataVector, Dim, Frame::Inertial>*>
224 : incoming_null_one_form,
225 : gsl::not_null<tnsr::a<DataVector, Dim, Frame::Inertial>*>
226 : outgoing_null_one_form,
227 : gsl::not_null<tnsr::A<DataVector, Dim, Frame::Inertial>*>
228 : incoming_null_vector,
229 : gsl::not_null<tnsr::A<DataVector, Dim, Frame::Inertial>*>
230 : outgoing_null_vector,
231 : gsl::not_null<tnsr::aa<DataVector, Dim, Frame::Inertial>*> projection_ab,
232 : gsl::not_null<tnsr::Ab<DataVector, Dim, Frame::Inertial>*> projection_Ab,
233 : gsl::not_null<tnsr::AA<DataVector, Dim, Frame::Inertial>*> projection_AB,
234 : gsl::not_null<tnsr::aa<DataVector, Dim, Frame::Inertial>*>
235 : char_projected_rhs_dt_v_psi,
236 : gsl::not_null<tnsr::iaa<DataVector, Dim, Frame::Inertial>*>
237 : char_projected_rhs_dt_v_zero,
238 : gsl::not_null<tnsr::aa<DataVector, Dim, Frame::Inertial>*>
239 : char_projected_rhs_dt_v_plus,
240 : gsl::not_null<tnsr::aa<DataVector, Dim, Frame::Inertial>*>
241 : char_projected_rhs_dt_v_minus,
242 : gsl::not_null<tnsr::a<DataVector, Dim, Frame::Inertial>*>
243 : constraint_char_zero_plus,
244 : gsl::not_null<tnsr::a<DataVector, Dim, Frame::Inertial>*>
245 : constraint_char_zero_minus,
246 : gsl::not_null<std::array<DataVector, 4>*> char_speeds,
247 :
248 : const std::optional<tnsr::I<DataVector, Dim, Frame::Inertial>>&
249 : face_mesh_velocity,
250 : const tnsr::i<DataVector, Dim, Frame::Inertial>& normal_covector,
251 : const tnsr::aa<DataVector, Dim, Frame::Inertial>& pi,
252 : const tnsr::iaa<DataVector, Dim, Frame::Inertial>& phi,
253 : const tnsr::aa<DataVector, Dim, Frame::Inertial>& spacetime_metric,
254 : const tnsr::I<DataVector, Dim, Frame::Inertial>& coords,
255 : const Scalar<DataVector>& gamma1, const Scalar<DataVector>& gamma2,
256 : const Scalar<DataVector>& lapse,
257 : const tnsr::I<DataVector, Dim, Frame::Inertial>& shift,
258 : const tnsr::AA<DataVector, Dim, Frame::Inertial>&
259 : inverse_spacetime_metric,
260 : const tnsr::A<DataVector, Dim, Frame::Inertial>&
261 : spacetime_unit_normal_vector,
262 : const tnsr::a<DataVector, Dim, Frame::Inertial>&
263 : spacetime_unit_normal_one_form,
264 : const tnsr::iaa<DataVector, Dim, Frame::Inertial>& three_index_constraint,
265 : const tnsr::a<DataVector, Dim, Frame::Inertial>& gauge_source,
266 : const tnsr::ab<DataVector, Dim, Frame::Inertial>&
267 : spacetime_deriv_gauge_source,
268 : const tnsr::aa<DataVector, Dim, Frame::Inertial>& dt_pi,
269 : const tnsr::iaa<DataVector, Dim, Frame::Inertial>& dt_phi,
270 : const tnsr::aa<DataVector, Dim, Frame::Inertial>& dt_spacetime_metric,
271 : const tnsr::iaa<DataVector, Dim, Frame::Inertial>& d_pi,
272 : const tnsr::ijaa<DataVector, Dim, Frame::Inertial>& d_phi,
273 : const tnsr::iaa<DataVector, Dim, Frame::Inertial>& d_spacetime_metric)
274 : const;
275 :
276 0 : detail::ConstraintPreservingBjorhusType type_{
277 : detail::ConstraintPreservingBjorhusType::ConstraintPreservingPhysical};
278 : };
279 : } // namespace gh::BoundaryConditions
280 :
281 : template <>
282 : struct Options::create_from_yaml<
283 : gh::BoundaryConditions::detail::ConstraintPreservingBjorhusType> {
284 : template <typename Metavariables>
285 : static
286 : typename gh::BoundaryConditions::detail::ConstraintPreservingBjorhusType
287 : create(const Options::Option& options) {
288 : return gh::BoundaryConditions::detail::
289 : convert_constraint_preserving_bjorhus_type_from_yaml(options);
290 : }
291 : };
|