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