Skip to content

Commit

Permalink
Add volume integrands for XCTS asymptotic quantities.
Browse files Browse the repository at this point in the history
  • Loading branch information
iago-mendes committed Sep 9, 2024
1 parent 8f7b3b6 commit 07af552
Show file tree
Hide file tree
Showing 11 changed files with 872 additions and 145 deletions.
9 changes: 6 additions & 3 deletions src/Elliptic/Systems/Xcts/Events/ObserveAdmIntegrals.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -131,6 +131,8 @@ void local_adm_integrals(
tenex::evaluate(mass_integrand(ti::I) * conformal_face_normal(ti::i));
const auto contracted_linear_momentum_integrand = tenex::evaluate<ti::I>(
linear_momentum_integrand(ti::I, ti::J) * conformal_face_normal(ti::j));
const auto contracted_center_of_mass_integrand = tenex::evaluate<ti::I>(
center_of_mass_integrand(ti::I, ti::J) * conformal_face_normal(ti::j));

// Take integrals
adm_mass->get() += definite_integral(
Expand All @@ -141,9 +143,10 @@ void local_adm_integrals(
definite_integral(contracted_linear_momentum_integrand.get(I) *
get(conformal_area_element),
face_mesh);
center_of_mass->get(I) += definite_integral(
center_of_mass_integrand.get(I) * get(conformal_area_element),
face_mesh);
center_of_mass->get(I) +=
definite_integral(contracted_center_of_mass_integrand.get(I) *
get(conformal_area_element),
face_mesh);
}
}
}
Expand Down
23 changes: 6 additions & 17 deletions src/PointwiseFunctions/Xcts/AdmLinearMomentum.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -3,19 +3,6 @@

#include "PointwiseFunctions/Xcts/AdmLinearMomentum.hpp"

#include <cmath>

#include "DataStructures/DataVector.hpp"
#include "DataStructures/Tensor/Slice.hpp"
#include "DataStructures/Tensor/Tensor.hpp"
#include "Domain/Structure/Direction.hpp"
#include "Domain/Structure/IndexToSliceAt.hpp"
#include "NumericalAlgorithms/LinearOperators/Divergence.hpp"
#include "NumericalAlgorithms/LinearOperators/PartialDerivatives.hpp"
#include "NumericalAlgorithms/Spectral/Mesh.hpp"
#include "Utilities/GenerateInstantiations.hpp"
#include "Utilities/Gsl.hpp"

namespace Xcts {

void adm_linear_momentum_surface_integrand(
Expand Down Expand Up @@ -47,34 +34,36 @@ void adm_linear_momentum_volume_integrand(
gsl::not_null<tnsr::I<DataVector, 3>*> result,
const tnsr::II<DataVector, 3>& surface_integrand,
const Scalar<DataVector>& conformal_factor,
const tnsr::i<DataVector, 3>& conformal_factor_deriv,
const tnsr::i<DataVector, 3>& deriv_conformal_factor,
const tnsr::ii<DataVector, 3>& conformal_metric,
const tnsr::II<DataVector, 3>& inv_conformal_metric,
const tnsr::Ijj<DataVector, 3>& conformal_christoffel_second_kind,
const tnsr::i<DataVector, 3>& conformal_christoffel_contracted) {
// Note: we can ignore the $1/(8\pi)$ term below because it is already
// included in `surface_integrand`.
tenex::evaluate<ti::I>(
result,
-(conformal_christoffel_second_kind(ti::I, ti::j, ti::k) *
surface_integrand(ti::J, ti::K) +
conformal_christoffel_contracted(ti::k) *
surface_integrand(ti::I, ti::K) -
2. * conformal_metric(ti::j, ti::k) * surface_integrand(ti::J, ti::K) *
inv_conformal_metric(ti::I, ti::L) * conformal_factor_deriv(ti::l) /
inv_conformal_metric(ti::I, ti::L) * deriv_conformal_factor(ti::l) /
conformal_factor()));
}

tnsr::I<DataVector, 3> adm_linear_momentum_volume_integrand(
const tnsr::II<DataVector, 3>& surface_integrand,
const Scalar<DataVector>& conformal_factor,
const tnsr::i<DataVector, 3>& conformal_factor_deriv,
const tnsr::i<DataVector, 3>& deriv_conformal_factor,
const tnsr::ii<DataVector, 3>& conformal_metric,
const tnsr::II<DataVector, 3>& inv_conformal_metric,
const tnsr::Ijj<DataVector, 3>& conformal_christoffel_second_kind,
const tnsr::i<DataVector, 3>& conformal_christoffel_contracted) {
tnsr::I<DataVector, 3> result;
adm_linear_momentum_volume_integrand(
make_not_null(&result), surface_integrand, conformal_factor,
conformal_factor_deriv, conformal_metric, inv_conformal_metric,
deriv_conformal_factor, conformal_metric, inv_conformal_metric,
conformal_christoffel_second_kind, conformal_christoffel_contracted);
return result;
}
Expand Down
50 changes: 23 additions & 27 deletions src/PointwiseFunctions/Xcts/AdmLinearMomentum.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -4,30 +4,26 @@
#pragma once

#include "DataStructures/DataVector.hpp"
#include "DataStructures/Tensor/Slice.hpp"
#include "DataStructures/Tensor/Tensor.hpp"
#include "Domain/Structure/Direction.hpp"
#include "Domain/Structure/IndexToSliceAt.hpp"
#include "NumericalAlgorithms/Spectral/Mesh.hpp"
#include "Utilities/Gsl.hpp"

#include <cmath>

#include "NumericalAlgorithms/LinearOperators/Divergence.hpp"
#include "NumericalAlgorithms/LinearOperators/PartialDerivatives.hpp"

namespace Xcts {

/// @{
/*!
* \brief Surface integrand for ADM linear momentum calculation defined as (see
* Eq. 20 in \cite Ossokine2015yla):
* \brief Surface integrand for the ADM linear momentum calculation.
*
* We define the ADM linear momentum integral as (see Eqs. 19-20 in
* \cite Ossokine2015yla):
*
* \begin{equation}
* \frac{1}{8\pi} \psi^{10} (K^{ij} - K \gamma^{ij})
* P_\text{ADM}^i = \frac{1}{8\pi}
* \oint_{S_\infty} \psi^10 \Big(
* K^{ij} - K \gamma^{ij}
* \Big) \, dS_j.
* \end{equation}
*
* \param result output buffer for the surface integrand
* \param result output pointer
* \param conformal_factor the conformal factor $\psi$
* \param inv_spatial_metric the inverse spatial metric $\gamma^{ij}$
* \param inv_extrinsic_curvature the inverse extrinsic curvature $K^{ij}$
Expand All @@ -54,23 +50,23 @@ tnsr::II<DataVector, 3> adm_linear_momentum_surface_integrand(
* Eq. 20 in \cite Ossokine2015yla):
*
* \begin{equation}
* - \frac{1}{8\pi} (
* \bar\Gamma^i_{jk} P^{jk}
* + \bar\Gamma^j_{jk} P^{jk}
* - 2 \bar\gamma_{jk} P^{jk} \bar\gamma^{il} \partial_l(\ln\psi)
* ),
* P_\text{ADM}^i = - \frac{1}{8\pi}
* \int_{V_\infty} \Big(
* \bar\Gamma^i_{jk} P^{jk}
* + \bar\Gamma^j_{jk} P^{jk}
* - 2 \bar\gamma_{jk} P^{jk} \bar\gamma^{il}
* \partial_l(\ln\psi)
* \Big) \, dV,
* \end{equation}
*
* where $\frac{1}{8\pi} P^{jk}$ is the result from
* where $1/(8\pi) P^{jk}$ is the result from
* `adm_linear_momentum_surface_integrand`.
*
* Note that we are including the negative sign in the integrand.
*
* \param result output buffer for the surface integrand
* \param surface_integrand the result of
* `adm_linear_momentum_surface_integrand`
* \param result output pointer
* \param surface_integrand the quantity $1/(8\pi) P^{ij}$ (result of
* `adm_linear_momentum_surface_integrand`)
* \param conformal_factor the conformal factor $\psi$
* \param conformal_factor_deriv the derivative of the conformal factor
* \param deriv_conformal_factor the gradient of the conformal factor
* $\partial_i\psi$
* \param conformal_metric the conformal metric $\bar\gamma_{ij}$
* \param inv_conformal_metric the inverse conformal metric $\bar\gamma^{ij}$
Expand All @@ -83,7 +79,7 @@ void adm_linear_momentum_volume_integrand(
gsl::not_null<tnsr::I<DataVector, 3>*> result,
const tnsr::II<DataVector, 3>& surface_integrand,
const Scalar<DataVector>& conformal_factor,
const tnsr::i<DataVector, 3>& conformal_factor_deriv,
const tnsr::i<DataVector, 3>& deriv_conformal_factor,
const tnsr::ii<DataVector, 3>& conformal_metric,
const tnsr::II<DataVector, 3>& inv_conformal_metric,
const tnsr::Ijj<DataVector, 3>& conformal_christoffel_second_kind,
Expand All @@ -93,7 +89,7 @@ void adm_linear_momentum_volume_integrand(
tnsr::I<DataVector, 3> adm_linear_momentum_volume_integrand(
const tnsr::II<DataVector, 3>& surface_integrand,
const Scalar<DataVector>& conformal_factor,
const tnsr::i<DataVector, 3>& conformal_factor_deriv,
const tnsr::i<DataVector, 3>& deriv_conformal_factor,
const tnsr::ii<DataVector, 3>& conformal_metric,
const tnsr::II<DataVector, 3>& inv_conformal_metric,
const tnsr::Ijj<DataVector, 3>& conformal_christoffel_second_kind,
Expand Down
52 changes: 52 additions & 0 deletions src/PointwiseFunctions/Xcts/AdmMass.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -33,4 +33,56 @@ tnsr::I<DataVector, 3> adm_mass_surface_integrand(
return result;
}

void adm_mass_volume_integrand(
gsl::not_null<Scalar<DataVector>*> result,
const Scalar<DataVector>& conformal_factor,
const Scalar<DataVector>& conformal_ricci_scalar,
const Scalar<DataVector>& trace_extrinsic_curvature,
const Scalar<DataVector>&
longitudinal_shift_minus_dt_conformal_metric_over_lapse_square,
const tnsr::II<DataVector, 3>& inv_conformal_metric,
const tnsr::iJK<DataVector, 3>& deriv_inv_conformal_metric,
const tnsr::Ijj<DataVector, 3>& conformal_christoffel_second_kind,
const tnsr::iJkk<DataVector, 3>& deriv_conformal_christoffel_second_kind) {
tenex::evaluate(
result,
1. / (16. * M_PI) *
(inv_conformal_metric(ti::J, ti::K) *
deriv_conformal_christoffel_second_kind(ti::i, ti::I, ti::j,
ti::k) +
conformal_christoffel_second_kind(ti::I, ti::j, ti::k) *
deriv_inv_conformal_metric(ti::i, ti::J, ti::K) -
inv_conformal_metric(ti::I, ti::J) *
deriv_conformal_christoffel_second_kind(ti::i, ti::K, ti::k,
ti::j) -
conformal_christoffel_second_kind(ti::K, ti::k, ti::j) *
deriv_inv_conformal_metric(ti::i, ti::I, ti::J) -
conformal_factor() * conformal_ricci_scalar() -
2. / 3. * pow<5>(conformal_factor()) *
square(trace_extrinsic_curvature()) +
pow<5>(conformal_factor()) / 4. *
longitudinal_shift_minus_dt_conformal_metric_over_lapse_square()));
}

Scalar<DataVector> adm_mass_volume_integrand(
const Scalar<DataVector>& conformal_factor,
const Scalar<DataVector>& conformal_ricci_scalar,
const Scalar<DataVector>& trace_extrinsic_curvature,
const Scalar<DataVector>&
longitudinal_shift_minus_dt_conformal_metric_over_lapse_square,
const tnsr::II<DataVector, 3>& inv_conformal_metric,
const tnsr::iJK<DataVector, 3>& deriv_inv_conformal_metric,
const tnsr::Ijj<DataVector, 3>& conformal_christoffel_second_kind,
const tnsr::iJkk<DataVector, 3>& deriv_conformal_christoffel_second_kind) {
Scalar<DataVector> result;
adm_mass_volume_integrand(
make_not_null(&result), conformal_factor, conformal_ricci_scalar,
trace_extrinsic_curvature,
longitudinal_shift_minus_dt_conformal_metric_over_lapse_square,
inv_conformal_metric, deriv_inv_conformal_metric,
conformal_christoffel_second_kind,
deriv_conformal_christoffel_second_kind);
return result;
}

} // namespace Xcts
92 changes: 81 additions & 11 deletions src/PointwiseFunctions/Xcts/AdmMass.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -16,20 +16,20 @@ namespace Xcts {
* We define the ADM mass integral as (see Eq. 3.139 in \cite BaumgarteShapiro):
*
* \begin{equation}
* M_{ADM}
* = \int_{S_\infty} \frac{1}{16\pi} (
* \bar\gamma^{jk} \bar\Gamma^i_{jk}
* - \bar\gamma^{ij} \bar\Gamma_{j}
* - 8 \bar\gamma^{ij} \partial_j \psi
* ) d\bar{S}_i.
* M_\text{ADM} = \frac{1}{16\pi}
* \oint_{S_\infty} \Big(
* \bar\gamma^{jk} \bar\Gamma^i_{jk}
* - \bar\gamma^{ij} \bar\Gamma_{j}
* - 8 \bar\gamma^{ij} \partial_j \psi
* \Big) dS_i.
* \end{equation}
*
* Note that we don't use the other versions presented in \cite BaumgarteShapiro
* of this integral because they make assumptions like $\bar\gamma = 1$,
* $\bar\Gamma^i_{ij} = 0$ and fast fall-off of the conformal metric.
* \note We don't use the other versions presented in \cite BaumgarteShapiro of
* this integral because they make assumptions like $\bar\gamma = 1$,
* $\bar\Gamma^i_{ij} = 0$ and faster fall-off of the conformal metric.
*
* \param result output buffer for the surface integrand
* \param deriv_conformal_factor the partial derivatives of the conformal factor
* \param result output pointer
* \param deriv_conformal_factor the gradient of the conformal factor
* $\partial_i \psi$
* \param inv_conformal_metric the inverse conformal metric $\bar\gamma^{ij}$
* \param conformal_christoffel_second_kind the conformal christoffel symbol
Expand All @@ -52,4 +52,74 @@ tnsr::I<DataVector, 3> adm_mass_surface_integrand(
const tnsr::i<DataVector, 3>& conformal_christoffel_contracted);
/// @}

/// @{
/*!
* \brief Volume integrand for the ADM mass calculation.
*
* We cast the ADM mass as an infinite volume integral by applying Gauss' law on
* the surface integral defined in `adm_mass_surface_integrand`:
*
* \begin{equation}
* M_\text{ADM} = \frac{1}{16\pi}
* \int_{V_\infty} \Big(
* \bar\gamma^{jk} \partial_i \bar\Gamma^i_{jk}
* + \bar\Gamma^i_{jk} \partial_i \bar\gamma^{jk}
* - \bar\gamma^{ij} \partial_i \bar\Gamma_j
* - \bar\Gamma_j \partial_i \bar\gamma^{ij}
* - 8 \bar D^2 \psi
* \Big) dV,
* \end{equation}
*
* where we can use the momentum constraint to replace $8 \bar D^2 \psi$ with
*
* \begin{equation}
* 8 \bar D^2 \psi = \psi \bar R + \frac{2}{3} \psi^5 K^2
* - \frac{1}{4} \psi^5 \frac{1}{\alpha^2}
* \Big[ (\bar L \beta)_{ij} - \bar u_{ij} \Big]
* \Big[ (\bar L \beta)^{ij} - \bar u^{ij} \Big].
* \end{equation}
*
* \note This is similar to Eq. 3.149 in \cite BaumgarteShapiro, except that
* here we don't assume $\bar\gamma = 1$.
*
* \param result output pointer
* \param conformal_factor the conformal factor
* \param conformal_ricci_scalar the conformal Ricci scalar $\bar R$
* \param trace_extrinsic_curvature the extrinsic curvature trace $K$
* \param longitudinal_shift_minus_dt_conformal_metric_over_lapse_square the
* quantity computed in
* `Xcts::Tags::LongitudinalShiftMinusDtConformalMetricOverLapseSquare`
* \param inv_conformal_metric the inverse conformal metric $\bar\gamma^{ij}$
* \param deriv_inv_conformal_metric the gradient of the inverse conformal
* metric $\partial_i \bar\gamma^{jk}$
* \param conformal_christoffel_second_kind the conformal christoffel symbol
* $\bar\Gamma^i_{jk}$
* \param deriv_conformal_christoffel_second_kind the gradient of the conformal
* christoffel symbol $\partial_i \bar\Gamma^j_{kl}$
*/
void adm_mass_volume_integrand(
gsl::not_null<Scalar<DataVector>*> result,
const Scalar<DataVector>& conformal_factor,
const Scalar<DataVector>& conformal_ricci_scalar,
const Scalar<DataVector>& trace_extrinsic_curvature,
const Scalar<DataVector>&
longitudinal_shift_minus_dt_conformal_metric_over_lapse_square,
const tnsr::II<DataVector, 3>& inv_conformal_metric,
const tnsr::iJK<DataVector, 3>& deriv_inv_conformal_metric,
const tnsr::Ijj<DataVector, 3>& conformal_christoffel_second_kind,
const tnsr::iJkk<DataVector, 3>& deriv_conformal_christoffel_second_kind);

/// Return-by-value overload
Scalar<DataVector> adm_mass_volume_integrand(
const Scalar<DataVector>& conformal_factor,
const Scalar<DataVector>& conformal_ricci_scalar,
const Scalar<DataVector>& trace_extrinsic_curvature,
const Scalar<DataVector>&
longitudinal_shift_minus_dt_conformal_metric_over_lapse_square,
const tnsr::II<DataVector, 3>& inv_conformal_metric,
const tnsr::iJK<DataVector, 3>& deriv_inv_conformal_metric,
const tnsr::Ijj<DataVector, 3>& conformal_christoffel_second_kind,
const tnsr::iJkk<DataVector, 3>& deriv_conformal_christoffel_second_kind);
/// @}

} // namespace Xcts
32 changes: 27 additions & 5 deletions src/PointwiseFunctions/Xcts/CenterOfMass.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -6,20 +6,42 @@
namespace Xcts {

void center_of_mass_surface_integrand(
gsl::not_null<tnsr::I<DataVector, 3>*> result,
gsl::not_null<tnsr::II<DataVector, 3>*> result,
const Scalar<DataVector>& conformal_factor,
const tnsr::I<DataVector, 3>& unit_normal) {
tenex::evaluate<ti::I>(result, 3. / (8. * M_PI) * pow<4>(conformal_factor()) *
unit_normal(ti::I));
tenex::evaluate<ti::I, ti::J>(result,
3. / (8. * M_PI) * pow<4>(conformal_factor()) *
unit_normal(ti::I) * unit_normal(ti::J));
}

tnsr::I<DataVector, 3> center_of_mass_surface_integrand(
tnsr::II<DataVector, 3> center_of_mass_surface_integrand(
const Scalar<DataVector>& conformal_factor,
const tnsr::I<DataVector, 3>& unit_normal) {
tnsr::I<DataVector, 3> result;
tnsr::II<DataVector, 3> result;
center_of_mass_surface_integrand(make_not_null(&result), conformal_factor,
unit_normal);
return result;
}

void center_of_mass_volume_integrand(
gsl::not_null<tnsr::I<DataVector, 3>*> result,
const Scalar<DataVector>& conformal_factor,
const tnsr::i<DataVector, 3, Frame::Inertial>& deriv_conformal_factor,
const tnsr::I<DataVector, 3>& unit_normal) {
tenex::evaluate<ti::I>(result, 3. / (8. * M_PI) * 4. *
pow<3>(conformal_factor()) *
deriv_conformal_factor(ti::j) *
unit_normal(ti::I) * unit_normal(ti::J));
}

tnsr::I<DataVector, 3> center_of_mass_volume_integrand(
const Scalar<DataVector>& conformal_factor,
const tnsr::i<DataVector, 3, Frame::Inertial>& deriv_conformal_factor,
const tnsr::I<DataVector, 3>& unit_normal) {
tnsr::I<DataVector, 3> result;
center_of_mass_volume_integrand(make_not_null(&result), conformal_factor,
deriv_conformal_factor, unit_normal);
return result;
}

} // namespace Xcts
Loading

0 comments on commit 07af552

Please sign in to comment.