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 :
8 : #include "DataStructures/Tensor/TypeAliases.hpp"
9 : #include "Domain/Tags.hpp"
10 : #include "Domain/TagsTimeDependent.hpp"
11 : #include "Evolution/Systems/GeneralizedHarmonic/ConstraintDamping/Tags.hpp"
12 : #include "Evolution/Systems/GeneralizedHarmonic/DuDtTempTags.hpp"
13 : #include "Evolution/Systems/GeneralizedHarmonic/GaugeSourceFunctions/Gauges.hpp"
14 : #include "Evolution/Systems/GeneralizedHarmonic/GaugeSourceFunctions/Tags/GaugeCondition.hpp"
15 : #include "Evolution/Systems/GeneralizedHarmonic/Tags.hpp"
16 : #include "PointwiseFunctions/GeneralRelativity/Tags.hpp"
17 : #include "PointwiseFunctions/GeneralRelativity/TagsDeclarations.hpp"
18 : #include "Utilities/TMPL.hpp"
19 :
20 : /// \cond
21 : class DataVector;
22 :
23 : namespace PUP {
24 : class er;
25 : } // namespace PUP
26 :
27 : namespace Tags {
28 : struct Time;
29 : } // namespace Tags
30 :
31 : namespace gsl {
32 : template <class T>
33 : class not_null;
34 : } // namespace gsl
35 :
36 : template <typename, typename, typename>
37 : class Tensor;
38 : /// \endcond
39 :
40 : namespace gh {
41 : /*!
42 : * \brief Compute the RHS of the Generalized Harmonic formulation of
43 : * Einstein's equations.
44 : *
45 : * The evolved variables are the spacetime metric \f$g_{ab}\f$, its spatial
46 : * derivative \f$\Phi_{iab}=\partial_i g_{ab}\f$, and conjugate momentum
47 : * \f$\Pi_{ab}=n^c\partial_c g_{ab}\f$, where \f$n^a\f$ is the spacetime
48 : * unit normal vector. The evolution equations are (Eqs. 35-57 of
49 : * \cite Lindblom2005qh)
50 : *
51 : * \f{align}{
52 : * \partial_t g_{ab}-
53 : * &\left(1+\gamma_1\right)\beta^k\partial_k g_{ab} - \gamma_1 v^i_g
54 : * \mathcal{C}_{iab} =
55 : * -\alpha \Pi_{ab}-\gamma_1\beta^i\Phi_{iab}, \\
56 : *
57 : * \partial_t\Pi_{ab}-
58 : * &\beta^k\partial_k\Pi_{ab} + \alpha \gamma^{ki}\partial_k\Phi_{iab}
59 : * - \gamma_1\gamma_2\beta^k\partial_kg_{ab} - \gamma_1 \gamma_2 v^i_g
60 : * \mathcal{C}_{iab} \notag\\
61 : * =&2\alpha g^{cd}\left(\gamma^{ij}\Phi_{ica}\Phi_{jdb}
62 : * - \Pi_{ca}\Pi_{db} - g^{ef}\Gamma_{ace}\Gamma_{bdf}\right) \notag \\
63 : * &-2\alpha \nabla_{(a}H_{b)}
64 : * - \frac{1}{2}\alpha n^c n^d\Pi_{cd}\Pi_{ab}
65 : * - \alpha n^c \Pi_{ci}\gamma^{ij}\Phi_{jab} \notag \\
66 : * &+\alpha \gamma_0\left(2\delta^c{}_{(a} n_{b)}
67 : * - (1 + \gamma_3)g_{ab}n^c\right)\mathcal{C}_c \notag \\
68 : * &+ 2 \gamma_4 \alpha \Pi_{ab} n^c \mathcal{C}_c \notag \\
69 : * &- \gamma_5\alpha n^c\mathcal{C}_c \left(\frac{\mathcal{C}_a\mathcal{C}_b
70 : * - \frac{1}{2} g_{ab} \mathcal{C}_d \mathcal{C}^d}
71 : * {\epsilon_{5} + 2 n^d \mathcal{C}_d n^e \mathcal{C}_e
72 : * + \mathcal{C}_d \mathcal{C}^d} \right) \notag \\
73 : * &-\gamma_1\gamma_2 \beta^i\Phi_{iab} \notag \\
74 : * &-16\pi \alpha \left(T_{ab} - \frac{1}{2}g_{ab}T^c{}_c\right),\\
75 : *
76 : * \partial_t\Phi_{iab}-
77 : * &\beta^k\partial_k\Phi_{iab} + \alpha \partial_i\Pi_{ab}
78 : * - \alpha \gamma_2\partial_ig_{ab} \notag \\
79 : * =&\frac{1}{2}\alpha n^c n^d\Phi_{icd}\Pi_{ab}
80 : * + \alpha \gamma^{jk}n^c\Phi_{ijc}\Phi_{kab} \notag \\
81 : * &-\alpha \gamma_2\Phi_{iab},
82 : * \f}
83 : *
84 : * where \f$H_a\f$ is the gauge source function,
85 : * \f$\mathcal{C}_a=H_a+\Gamma_a\f$ is the gauge constraint,
86 : * \f$\mathcal{C}_{iab}\f$ is the 3-index constraint, and
87 : * \f$v^i_g\f$ is the mesh velocity. The constraint
88 : * damping parameters \f$\gamma_0\f$ \f$\gamma_1\f$, \f$\gamma_2\f$,
89 : * \f$\gamma_3\f$, \f$\gamma_4\f$, and \f$\gamma_5\f$ have units of inverse time
90 : * and control the time scales on which the constraints are damped to zero.
91 : *
92 : * \note We have not coded up the constraint damping terms for \f$\gamma_3\f$,
93 : * \f$\gamma_4\f$, and \f$\gamma_5\f$. \f$\gamma_3\f$ was found to be essential
94 : * for evolutions of black strings by Pretorius and Lehner \cite Lehner2010pn.
95 : *
96 : * \note Here the only terms dependent on mesh velocity are constraint damping
97 : * terms added to this implementation of the generalized harmonic system.
98 : * Mesh-velocity corrections that are applicable to all systems are made in
99 : * `evolution::dg::Actions::detail::volume_terms()`.
100 : *
101 : * \warning When using harmonic gauge,
102 : * gr::Tags::SqrtDetSpatialMetric<DataVector> and
103 : * gr::Tags::SpacetimeChristoffelSecondKind<Dim, Frame::Inertial, DataVector>
104 : * are not computed. In Debug mode, they are filled with with signaling NaNs.
105 : */
106 : template <size_t Dim>
107 1 : struct TimeDerivative {
108 : public:
109 0 : using temporary_tags = tmpl::list<
110 : ::gh::ConstraintDamping::Tags::ConstraintGamma1,
111 : ::gh::ConstraintDamping::Tags::ConstraintGamma2,
112 : Tags::GaugeH<DataVector, Dim>,
113 : Tags::SpacetimeDerivGaugeH<DataVector, Dim>, Tags::Gamma1Gamma2,
114 : Tags::HalfPiTwoNormals, Tags::NormalDotOneIndexConstraint,
115 : Tags::Gamma1Plus1, Tags::PiOneNormal<Dim>,
116 : Tags::GaugeConstraint<DataVector, Dim>, Tags::HalfPhiTwoNormals<Dim>,
117 : Tags::ShiftDotThreeIndexConstraint<Dim>,
118 : Tags::MeshVelocityDotThreeIndexConstraint<Dim>, Tags::PhiOneNormal<Dim>,
119 : Tags::PiSecondIndexUp<Dim>, Tags::ThreeIndexConstraint<DataVector, Dim>,
120 : Tags::PhiFirstIndexUp<Dim>, Tags::PhiThirdIndexUp<Dim>,
121 : Tags::SpacetimeChristoffelFirstKindThirdIndexUp<Dim>,
122 : gr::Tags::Lapse<DataVector>, gr::Tags::Shift<DataVector, Dim>,
123 : gr::Tags::InverseSpatialMetric<DataVector, Dim>,
124 : gr::Tags::DetSpatialMetric<DataVector>,
125 : gr::Tags::SqrtDetSpatialMetric<DataVector>,
126 : gr::Tags::InverseSpacetimeMetric<DataVector, Dim>,
127 : gr::Tags::SpacetimeChristoffelFirstKind<DataVector, Dim>,
128 : gr::Tags::SpacetimeChristoffelSecondKind<DataVector, Dim>,
129 : gr::Tags::TraceSpacetimeChristoffelFirstKind<DataVector, Dim>,
130 : gr::Tags::SpacetimeNormalVector<DataVector, Dim>>;
131 0 : using argument_tags =
132 : tmpl::list<gr::Tags::SpacetimeMetric<DataVector, Dim>,
133 : Tags::Pi<DataVector, Dim>, Tags::Phi<DataVector, Dim>,
134 : ::gh::ConstraintDamping::Tags::ConstraintGamma0,
135 : ::gh::ConstraintDamping::Tags::ConstraintGamma1,
136 : ::gh::ConstraintDamping::Tags::ConstraintGamma2,
137 : gauges::Tags::GaugeCondition, domain::Tags::Mesh<Dim>,
138 : ::Tags::Time, domain::Tags::Coordinates<Dim, Frame::Inertial>,
139 : domain::Tags::InverseJacobian<Dim, Frame::ElementLogical,
140 : Frame::Inertial>,
141 : domain::Tags::MeshVelocity<Dim, Frame::Inertial>>;
142 :
143 0 : static void apply(
144 : gsl::not_null<tnsr::aa<DataVector, Dim>*> dt_spacetime_metric,
145 : gsl::not_null<tnsr::aa<DataVector, Dim>*> dt_pi,
146 : gsl::not_null<tnsr::iaa<DataVector, Dim>*> dt_phi,
147 : gsl::not_null<Scalar<DataVector>*> temp_gamma1,
148 : gsl::not_null<Scalar<DataVector>*> temp_gamma2,
149 : gsl::not_null<tnsr::a<DataVector, Dim>*> temp_gauge_function,
150 : gsl::not_null<tnsr::ab<DataVector, Dim>*>
151 : temp_spacetime_deriv_gauge_function,
152 : gsl::not_null<Scalar<DataVector>*> gamma1gamma2,
153 : gsl::not_null<Scalar<DataVector>*> half_half_pi_two_normals,
154 : gsl::not_null<Scalar<DataVector>*> normal_dot_gauge_constraint,
155 : gsl::not_null<Scalar<DataVector>*> gamma1_plus_1,
156 : gsl::not_null<tnsr::a<DataVector, Dim>*> pi_one_normal,
157 : gsl::not_null<tnsr::a<DataVector, Dim>*> gauge_constraint,
158 : gsl::not_null<tnsr::i<DataVector, Dim>*> half_phi_two_normals,
159 : gsl::not_null<tnsr::aa<DataVector, Dim>*>
160 : shift_dot_three_index_constraint,
161 : gsl::not_null<tnsr::aa<DataVector, Dim>*>
162 : mesh_velocity_dot_three_index_constraint,
163 : gsl::not_null<tnsr::ia<DataVector, Dim>*> phi_one_normal,
164 : gsl::not_null<tnsr::aB<DataVector, Dim>*> pi_2_up,
165 : gsl::not_null<tnsr::iaa<DataVector, Dim>*> three_index_constraint,
166 : gsl::not_null<tnsr::Iaa<DataVector, Dim>*> phi_1_up,
167 : gsl::not_null<tnsr::iaB<DataVector, Dim>*> phi_3_up,
168 : gsl::not_null<tnsr::abC<DataVector, Dim>*> christoffel_first_kind_3_up,
169 : gsl::not_null<Scalar<DataVector>*> lapse,
170 : gsl::not_null<tnsr::I<DataVector, Dim>*> shift,
171 : gsl::not_null<tnsr::II<DataVector, Dim>*> inverse_spatial_metric,
172 : gsl::not_null<Scalar<DataVector>*> det_spatial_metric,
173 : gsl::not_null<Scalar<DataVector>*> sqrt_det_spatial_metric,
174 : gsl::not_null<tnsr::AA<DataVector, Dim>*> inverse_spacetime_metric,
175 : gsl::not_null<tnsr::abb<DataVector, Dim>*> christoffel_first_kind,
176 : gsl::not_null<tnsr::Abb<DataVector, Dim>*> christoffel_second_kind,
177 : gsl::not_null<tnsr::a<DataVector, Dim>*> trace_christoffel,
178 : gsl::not_null<tnsr::A<DataVector, Dim>*> normal_spacetime_vector,
179 : const tnsr::iaa<DataVector, Dim>& d_spacetime_metric,
180 : const tnsr::iaa<DataVector, Dim>& d_pi,
181 : const tnsr::ijaa<DataVector, Dim>& d_phi,
182 : const tnsr::aa<DataVector, Dim>& spacetime_metric,
183 : const tnsr::aa<DataVector, Dim>& pi,
184 : const tnsr::iaa<DataVector, Dim>& phi, const Scalar<DataVector>& gamma0,
185 : const Scalar<DataVector>& gamma1, const Scalar<DataVector>& gamma2,
186 : const gauges::GaugeCondition& gauge_condition, const Mesh<Dim>& mesh,
187 : double time,
188 : const tnsr::I<DataVector, Dim, Frame::Inertial>& inertial_coords,
189 : const InverseJacobian<DataVector, Dim, Frame::ElementLogical,
190 : Frame::Inertial>& inverse_jacobian,
191 : const std::optional<tnsr::I<DataVector, Dim, Frame::Inertial>>&
192 : mesh_velocity);
193 : };
194 : } // namespace gh
|