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 <pup.h>
8 :
9 : #include "DataStructures/ComplexDataVector.hpp"
10 : #include "DataStructures/DataVector.hpp"
11 : #include "DataStructures/Tensor/Tensor.hpp"
12 : #include "DataStructures/Variables.hpp"
13 : #include "DataStructures/VariablesTag.hpp"
14 : #include "Elliptic/Systems/SelfForce/GeneralRelativity/Tags.hpp"
15 : #include "NumericalAlgorithms/DiscontinuousGalerkin/Tags.hpp"
16 : #include "PointwiseFunctions/GeneralRelativity/Tags.hpp"
17 : #include "Utilities/Gsl.hpp"
18 : #include "Utilities/MakeWithValue.hpp"
19 : #include "Utilities/TMPL.hpp"
20 :
21 : namespace GrSelfForce {
22 :
23 : /// @{
24 : /// We're working with 4D tensors to represent the 10 independent components
25 : /// we're solving for, but we only take 2D spatial derivatives, so we define
26 : /// these mixed-dimension tensors for gradients and fluxes.
27 1 : using GradTensorType =
28 : TensorMetafunctions::prepend_spatial_index<tnsr::aa<ComplexDataVector, 3>,
29 : 2, UpLo::Lo, Frame::Inertial>;
30 1 : using FluxTensorType =
31 : TensorMetafunctions::prepend_spatial_index<tnsr::aa<ComplexDataVector, 3>,
32 : 2, UpLo::Up, Frame::Inertial>;
33 : /// @}
34 :
35 : /*!
36 : * \brief The first-order flux $F^i=\{\partial_{r_\star}, \alpha
37 : * \partial_\theta\}\Psi_m$.
38 : */
39 1 : void fluxes(gsl::not_null<FluxTensorType*> flux,
40 : const Scalar<ComplexDataVector>& alpha,
41 : const GradTensorType& field_gradient);
42 :
43 : /*!
44 : * \brief The first-order flux on an element face
45 : * $F^i=\{n_{r_\star}, \alpha n_\theta\}\Psi_m$.
46 : */
47 1 : void fluxes_on_face(gsl::not_null<FluxTensorType*> flux,
48 : const Scalar<ComplexDataVector>& alpha,
49 : const tnsr::I<DataVector, 2>& face_normal_vector,
50 : const tnsr::aa<ComplexDataVector, 3>& field);
51 :
52 : /*!
53 : * \brief The source term $\beta_{ab}^{cd} (\Psi_m)_{cd} + \gamma_{iab}^{cd}
54 : * F^i_{cd}$.
55 : */
56 1 : void add_sources(gsl::not_null<tnsr::aa<ComplexDataVector, 3>*> source,
57 : const tnsr::aaBB<ComplexDataVector, 3>& beta,
58 : const tnsr::aaBB<ComplexDataVector, 3>& gamma_rstar,
59 : const tnsr::aaBB<ComplexDataVector, 3>& gamma_theta,
60 : const tnsr::aa<ComplexDataVector, 3>& field,
61 : const FluxTensorType& flux);
62 :
63 : /// Fluxes $F^i$ for the gravitational self-force system.
64 : /// \see GrSelfForce::FirstOrderSystem
65 1 : struct Fluxes {
66 0 : using argument_tags = tmpl::list<Tags::Alpha>;
67 0 : using volume_tags = tmpl::list<>;
68 0 : using const_global_cache_tags = tmpl::list<>;
69 0 : static constexpr bool is_trivial = false;
70 0 : static constexpr bool is_discontinuous = false;
71 0 : static void apply(gsl::not_null<FluxTensorType*> flux,
72 : const Scalar<ComplexDataVector>& alpha,
73 : const tnsr::aa<ComplexDataVector, 3>& /*field*/,
74 : const GradTensorType& field_gradient);
75 0 : static void apply(gsl::not_null<FluxTensorType*> flux,
76 : const Scalar<ComplexDataVector>& alpha,
77 : const tnsr::i<DataVector, 2>& /*face_normal*/,
78 : const tnsr::I<DataVector, 2>& face_normal_vector,
79 : const tnsr::aa<ComplexDataVector, 3>& field);
80 : };
81 :
82 : /// Source terms for the gravitational self-force system.
83 : /// \see GrSelfForce::FirstOrderSystem
84 1 : struct Sources {
85 0 : using argument_tags =
86 : tmpl::list<Tags::Beta, Tags::GammaRstar, Tags::GammaTheta>;
87 0 : using const_global_cache_tags = tmpl::list<>;
88 0 : static void apply(
89 : gsl::not_null<tnsr::aa<ComplexDataVector, 3>*> scalar_equation,
90 : const tnsr::aaBB<ComplexDataVector, 3>& beta,
91 : const tnsr::aaBB<ComplexDataVector, 3>& gamma_rstar,
92 : const tnsr::aaBB<ComplexDataVector, 3>& gamma_theta,
93 : const tnsr::aa<ComplexDataVector, 3>& field,
94 : const GradTensorType& /*field_gradient*/, const FluxTensorType& flux);
95 : };
96 :
97 : /*!
98 : * \brief Adds or subtracts the singular field to/from the received data on
99 : * element boundaries.
100 : *
101 : * In the regularized region we solve for the regularized field
102 : * \begin{equation}
103 : * \Psi_m^R = \Psi_m - \Psi_m^P
104 : * \text{,}
105 : * \end{equation}
106 : * so we subtract the singular field on the regularized side (where
107 : * `field_is_regularized` is true) and add it on the other side of the boundary
108 : * (where `field_is_regularized` is false). We do the same for the received
109 : * normal dot flux $n_i F^i$, but with an extra minus sign because this quantity
110 : * is defined with the face normal from the perspective of the sending element
111 : * (see `elliptic::protocols::FirstOrderSystem`).
112 : */
113 1 : struct ModifyBoundaryData {
114 : private:
115 0 : static constexpr size_t Dim = 2;
116 0 : using singular_vars_on_mortars_tag =
117 : ::Tags::Variables<tmpl::list<Tags::SingularField,
118 : ::Tags::NormalDotFlux<Tags::SingularField>>>;
119 :
120 : public:
121 0 : using argument_tags =
122 : tmpl::list<Tags::FieldIsRegularized,
123 : ::Tags::Mortars<Tags::FieldIsRegularized, Dim>,
124 : ::Tags::Mortars<singular_vars_on_mortars_tag, Dim>>;
125 0 : static void apply(
126 : gsl::not_null<tnsr::aa<ComplexDataVector, 3>*> field,
127 : gsl::not_null<tnsr::aa<ComplexDataVector, 3>*> n_dot_flux,
128 : const DirectionalId<Dim>& mortar_id, bool field_is_regularized,
129 : const DirectionalIdMap<Dim, bool>& neighbors_field_is_regularized,
130 : const DirectionalIdMap<Dim, typename singular_vars_on_mortars_tag::type>&
131 : singular_vars_on_mortars);
132 : };
133 :
134 : } // namespace GrSelfForce
|