SpECTRE  v2024.09.29
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Modules Pages
Xcts Namespace Reference

Items related to solving the Extended Conformal Thin Sandwich (XCTS) decomposition of the Einstein constraint equations. More...

Namespaces

namespace  AnalyticData
 Analytic data for the XCTS equations, i.e. field configurations that do not solve the equations but are used as background or initial guess.
 
namespace  Solutions
 Analytic solutions of the XCTS equations.
 
namespace  Tags
 Tags related to the XCTS equations.
 

Classes

struct  FirstOrderSystem
 The Extended Conformal Thin Sandwich (XCTS) decomposition of the Einstein constraint equations, formulated as a set of coupled first-order partial differential equations. More...
 
struct  Fluxes
 The fluxes Fi for the first-order formulation of the XCTS equations. More...
 
struct  LinearizedSources
 The linearization of the sources S for the first-order formulation of the XCTS equations. More...
 
struct  Sources
 The sources S for the first-order formulation of the XCTS equations. More...
 
struct  SpacetimeQuantitiesComputer
 CachedTempBuffer computer class for 3+1 quantities from XCTS variables. See Xcts::SpacetimeQuantities. More...
 

Typedefs

using SpacetimeQuantities = CachedTempBuffer< Tags::ConformalFactor< DataVector >, Tags::LapseTimesConformalFactor< DataVector >, ::Tags::deriv< Tags::ConformalFactor< DataVector >, tmpl::size_t< 3 >, Frame::Inertial >, ::Tags::deriv<::Tags::deriv< Tags::ConformalFactor< DataVector >, tmpl::size_t< 3 >, Frame::Inertial >, tmpl::size_t< 3 >, Frame::Inertial >, detail::ConformalLaplacianOfConformalFactor< DataVector >, ::Tags::deriv< Tags::LapseTimesConformalFactor< DataVector >, tmpl::size_t< 3 >, Frame::Inertial >, ::Tags::deriv< Tags::ShiftExcess< DataVector, 3, Frame::Inertial >, tmpl::size_t< 3 >, Frame::Inertial >, Xcts::Tags::ShiftStrain< DataVector, 3, Frame::Inertial >, Xcts::Tags::LongitudinalShiftExcess< DataVector, 3, Frame::Inertial >, ::Tags::div< Tags::LongitudinalShiftExcess< DataVector, 3, Frame::Inertial > >, detail::LongitudinalShiftMinusDtConformalMetric< DataVector >, gr::Tags::SpatialMetric< DataVector, 3 >, gr::Tags::InverseSpatialMetric< DataVector, 3 >, gr::Tags::Lapse< DataVector >, gr::Tags::Shift< DataVector, 3 >, gr::Tags::ExtrinsicCurvature< DataVector, 3 >, gr::Tags::SpatialChristoffelSecondKind< DataVector, 3 >, gr::Tags::SpatialRicci< DataVector, 3 >, gr::Tags::HamiltonianConstraint< DataVector >, gr::Tags::MomentumConstraint< DataVector, 3 > >
 General-relativistic 3+1 quantities computed from XCTS variables.
 

Enumerations

enum class  Equations { Hamiltonian , HamiltonianAndLapse , HamiltonianLapseAndShift }
 Indicates a subset of the XCTS equations. More...
 
enum class  Geometry { FlatCartesian , Curved }
 Types of conformal geometries for the XCTS equations. More...
 

Functions

template<int ConformalMatterScale>
void add_hamiltonian_sources (gsl::not_null< Scalar< DataVector > * > hamiltonian_constraint, const Scalar< DataVector > &conformal_energy_density, const Scalar< DataVector > &extrinsic_curvature_trace, const Scalar< DataVector > &conformal_factor_minus_one)
 Add the nonlinear source to the Hamiltonian constraint on a flat conformal background in Cartesian coordinates and with u¯ij=0=βi. More...
 
template<int ConformalMatterScale>
void add_linearized_hamiltonian_sources (gsl::not_null< Scalar< DataVector > * > linearized_hamiltonian_constraint, const Scalar< DataVector > &conformal_energy_density, const Scalar< DataVector > &extrinsic_curvature_trace, const Scalar< DataVector > &conformal_factor_minus_one, const Scalar< DataVector > &conformal_factor_correction)
 The linearization of add_hamiltonian_sources
 
void add_distortion_hamiltonian_sources (gsl::not_null< Scalar< DataVector > * > hamiltonian_constraint, const Scalar< DataVector > &longitudinal_shift_minus_dt_conformal_metric_over_lapse_square, const Scalar< DataVector > &conformal_factor_minus_one)
 Add the "distortion" source term to the Hamiltonian constraint. More...
 
void add_linearized_distortion_hamiltonian_sources (gsl::not_null< Scalar< DataVector > * > linearized_hamiltonian_constraint, const Scalar< DataVector > &longitudinal_shift_minus_dt_conformal_metric_over_lapse_square, const Scalar< DataVector > &conformal_factor_minus_one, const Scalar< DataVector > &conformal_factor_correction)
 The linearization of add_distortion_hamiltonian_sources More...
 
void add_curved_hamiltonian_or_lapse_sources (gsl::not_null< Scalar< DataVector > * > hamiltonian_or_lapse_equation, const Scalar< DataVector > &conformal_ricci_scalar, const Scalar< DataVector > &field, double add_constant=0.)
 Add the contributions from a curved background geometry to the Hamiltonian constraint or lapse equation. More...
 
template<int ConformalMatterScale>
void add_lapse_sources (gsl::not_null< Scalar< DataVector > * > lapse_equation, const Scalar< DataVector > &conformal_energy_density, const Scalar< DataVector > &conformal_stress_trace, const Scalar< DataVector > &extrinsic_curvature_trace, const Scalar< DataVector > &dt_extrinsic_curvature_trace, const Scalar< DataVector > &shift_dot_deriv_extrinsic_curvature_trace, const Scalar< DataVector > &conformal_factor_minus_one, const Scalar< DataVector > &lapse_times_conformal_factor_minus_one)
 Add the nonlinear source to the lapse equation on a flat conformal background in Cartesian coordinates and with u¯ij=0=βi. More...
 
template<int ConformalMatterScale>
void add_linearized_lapse_sources (gsl::not_null< Scalar< DataVector > * > linearized_lapse_equation, const Scalar< DataVector > &conformal_energy_density, const Scalar< DataVector > &conformal_stress_trace, const Scalar< DataVector > &extrinsic_curvature_trace, const Scalar< DataVector > &dt_extrinsic_curvature_trace, const Scalar< DataVector > &shift_dot_deriv_extrinsic_curvature_trace, const Scalar< DataVector > &conformal_factor_minus_one, const Scalar< DataVector > &lapse_times_conformal_factor_minus_one, const Scalar< DataVector > &conformal_factor_correction, const Scalar< DataVector > &lapse_times_conformal_factor_correction)
 The linearization of add_lapse_sources More...
 
void add_distortion_hamiltonian_and_lapse_sources (gsl::not_null< Scalar< DataVector > * > hamiltonian_constraint, gsl::not_null< Scalar< DataVector > * > lapse_equation, const Scalar< DataVector > &longitudinal_shift_minus_dt_conformal_metric_square, const Scalar< DataVector > &conformal_factor_minus_one, const Scalar< DataVector > &lapse_times_conformal_factor_minus_one)
 Add the "distortion" source term to the Hamiltonian constraint and the lapse equation. More...
 
void add_linearized_distortion_hamiltonian_and_lapse_sources (gsl::not_null< Scalar< DataVector > * > hamiltonian_constraint, gsl::not_null< Scalar< DataVector > * > lapse_equation, const Scalar< DataVector > &longitudinal_shift_minus_dt_conformal_metric_square, const Scalar< DataVector > &conformal_factor_minus_one, const Scalar< DataVector > &lapse_times_conformal_factor_minus_one, const Scalar< DataVector > &conformal_factor_correction, const Scalar< DataVector > &lapse_times_conformal_factor_correction)
 The linearization of add_distortion_hamiltonian_and_lapse_sources More...
 
template<int ConformalMatterScale>
void add_curved_linearized_momentum_sources (gsl::not_null< Scalar< DataVector > * > linearized_hamiltonian_constraint, gsl::not_null< Scalar< DataVector > * > linearized_lapse_equation, gsl::not_null< tnsr::I< DataVector, 3 > * > linearized_momentum_constraint, const tnsr::I< DataVector, 3 > &conformal_momentum_density, const tnsr::i< DataVector, 3 > &extrinsic_curvature_trace_gradient, const tnsr::ii< DataVector, 3 > &conformal_metric, const tnsr::II< DataVector, 3 > &inv_conformal_metric, const Scalar< DataVector > &conformal_factor_minus_one, const Scalar< DataVector > &lapse_times_conformal_factor_minus_one, const tnsr::I< DataVector, 3 > &conformal_factor_flux, const tnsr::I< DataVector, 3 > &lapse_times_conformal_factor_flux, const tnsr::II< DataVector, 3 > &longitudinal_shift_minus_dt_conformal_metric, const Scalar< DataVector > &conformal_factor_correction, const Scalar< DataVector > &lapse_times_conformal_factor_correction, const tnsr::I< DataVector, 3 > &shift_correction, const tnsr::I< DataVector, 3 > &conformal_factor_flux_correction, const tnsr::I< DataVector, 3 > &lapse_times_conformal_factor_flux_correction, const tnsr::II< DataVector, 3 > &longitudinal_shift_correction)
 The linearization of add_curved_momentum_sources
 
template<int ConformalMatterScale>
void add_flat_cartesian_linearized_momentum_sources (gsl::not_null< Scalar< DataVector > * > linearized_hamiltonian_constraint, gsl::not_null< Scalar< DataVector > * > linearized_lapse_equation, gsl::not_null< tnsr::I< DataVector, 3 > * > linearized_momentum_constraint, const tnsr::I< DataVector, 3 > &conformal_momentum_density, const tnsr::i< DataVector, 3 > &extrinsic_curvature_trace_gradient, const Scalar< DataVector > &conformal_factor_minus_one, const Scalar< DataVector > &lapse_times_conformal_factor_minus_one, const tnsr::I< DataVector, 3 > &conformal_factor_flux, const tnsr::I< DataVector, 3 > &lapse_times_conformal_factor_flux, const tnsr::II< DataVector, 3 > &longitudinal_shift_minus_dt_conformal_metric, const Scalar< DataVector > &conformal_factor_correction, const Scalar< DataVector > &lapse_times_conformal_factor_correction, const tnsr::I< DataVector, 3 > &shift_correction, const tnsr::I< DataVector, 3 > &conformal_factor_flux_correction, const tnsr::I< DataVector, 3 > &lapse_times_conformal_factor_flux_correction, const tnsr::II< DataVector, 3 > &longitudinal_shift_correction)
 The linearization of add_flat_cartesian_momentum_sources
 
template<int ConformalMatterScale>
void add_curved_momentum_sources (gsl::not_null< Scalar< DataVector > * > hamiltonian_constraint, gsl::not_null< Scalar< DataVector > * > lapse_equation, gsl::not_null< tnsr::I< DataVector, 3 > * > momentum_constraint, const tnsr::I< DataVector, 3 > &conformal_momentum_density, const tnsr::i< DataVector, 3 > &extrinsic_curvature_trace_gradient, const tnsr::ii< DataVector, 3 > &conformal_metric, const tnsr::II< DataVector, 3 > &inv_conformal_metric, const tnsr::I< DataVector, 3 > &minus_div_dt_conformal_metric, const Scalar< DataVector > &conformal_factor_minus_one, const Scalar< DataVector > &lapse_times_conformal_factor_minus_one, const tnsr::I< DataVector, 3 > &conformal_factor_flux, const tnsr::I< DataVector, 3 > &lapse_times_conformal_factor_flux, const tnsr::II< DataVector, 3 > &longitudinal_shift_minus_dt_conformal_metric)
 Add the nonlinear source to the momentum constraint and add the "distortion" source term to the Hamiltonian constraint and lapse equation. More...
 
template<int ConformalMatterScale>
void add_flat_cartesian_momentum_sources (gsl::not_null< Scalar< DataVector > * > hamiltonian_constraint, gsl::not_null< Scalar< DataVector > * > lapse_equation, gsl::not_null< tnsr::I< DataVector, 3 > * > momentum_constraint, const tnsr::I< DataVector, 3 > &conformal_momentum_density, const tnsr::i< DataVector, 3 > &extrinsic_curvature_trace_gradient, const tnsr::I< DataVector, 3 > &minus_div_dt_conformal_metric, const Scalar< DataVector > &conformal_factor_minus_one, const Scalar< DataVector > &lapse_times_conformal_factor_minus_one, const tnsr::I< DataVector, 3 > &conformal_factor_flux, const tnsr::I< DataVector, 3 > &lapse_times_conformal_factor_flux, const tnsr::II< DataVector, 3 > &longitudinal_shift_minus_dt_conformal_metric)
 Add the nonlinear source to the momentum constraint and add the "distortion" source term to the Hamiltonian constraint and lapse equation. More...
 
void adm_linear_momentum_surface_integrand (gsl::not_null< tnsr::II< DataVector, 3 > * > result, const Scalar< DataVector > &conformal_factor, const tnsr::II< DataVector, 3 > &inv_spatial_metric, const tnsr::II< DataVector, 3 > &inv_extrinsic_curvature, const Scalar< DataVector > &trace_extrinsic_curvature)
 Surface integrand for the ADM linear momentum calculation. More...
 
tnsr::II< DataVector, 3 > adm_linear_momentum_surface_integrand (const Scalar< DataVector > &conformal_factor, const tnsr::II< DataVector, 3 > &inv_spatial_metric, const tnsr::II< DataVector, 3 > &inv_extrinsic_curvature, const Scalar< DataVector > &trace_extrinsic_curvature)
 Return-by-value overload.
 
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 > &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)
 Volume integrand for ADM linear momentum calculation defined as (see Eq. 20 in ): More...
 
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 > &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)
 Return-by-value overload.
 
void adm_mass_surface_integrand (gsl::not_null< tnsr::I< DataVector, 3 > * > result, const tnsr::i< DataVector, 3 > &deriv_conformal_factor, 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)
 Surface integrand for the ADM mass calculation. More...
 
tnsr::I< DataVector, 3 > adm_mass_surface_integrand (const tnsr::i< DataVector, 3 > &deriv_conformal_factor, 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)
 Return-by-value overload.
 
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 Scalar< DataVector > &energy_density, 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::i< DataVector, 3 > &conformal_christoffel_contracted, const tnsr::iJkk< DataVector, 3 > &deriv_conformal_christoffel_second_kind)
 Volume integrand for the ADM mass calculation. More...
 
Scalar< DataVectoradm_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 Scalar< DataVector > &energy_density, 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::i< DataVector, 3 > &conformal_christoffel_contracted, const tnsr::iJkk< DataVector, 3 > &deriv_conformal_christoffel_second_kind)
 Return-by-value overload.
 
void center_of_mass_surface_integrand (gsl::not_null< tnsr::I< DataVector, 3 > * > result, const Scalar< DataVector > &conformal_factor, const tnsr::I< DataVector, 3 > &coords)
 Surface integrand for the center of mass calculation. More...
 
tnsr::I< DataVector, 3 > center_of_mass_surface_integrand (const Scalar< DataVector > &conformal_factor, const tnsr::I< DataVector, 3 > &coords)
 Return-by-value overload.
 
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 > &coords)
 Volume integrand for the center of mass calculation. More...
 
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 > &coords)
 Return-by-value overload.
 
template<typename DataType >
void extrinsic_curvature (const gsl::not_null< tnsr::ii< DataType, 3 > * > result, const Scalar< DataType > &conformal_factor, const Scalar< DataType > &lapse, const tnsr::ii< DataType, 3 > &conformal_metric, const tnsr::II< DataType, 3 > &longitudinal_shift_minus_dt_conformal_metric, const Scalar< DataType > &trace_extrinsic_curvature)
 Extrinsic curvature computed from the conformal decomposition used in the XCTS system. More...
 
template<typename DataType >
tnsr::ii< DataType, 3 > extrinsic_curvature (const Scalar< DataType > &conformal_factor, const Scalar< DataType > &lapse, const tnsr::ii< DataType, 3 > &conformal_metric, const tnsr::II< DataType, 3 > &longitudinal_shift_minus_dt_conformal_metric, const Scalar< DataType > &trace_extrinsic_curvature)
 Return-by-value overload.
 
template<typename DataType >
void longitudinal_operator (gsl::not_null< tnsr::II< DataType, 3 > * > result, const tnsr::ii< DataType, 3 > &strain, const tnsr::II< DataType, 3 > &inv_metric)
 The longitudinal operator, or vector gradient, (Lβ)ij. More...
 
template<typename DataType >
void longitudinal_operator (gsl::not_null< tnsr::II< DataType, 3 > * > result, const tnsr::I< DataType, 3 > &shift, const tnsr::iJ< DataType, 3 > &deriv_shift, const tnsr::II< DataType, 3 > &inv_metric, const tnsr::Ijj< DataType, 3 > &christoffel_second_kind)
 The longitudinal operator, or vector gradient, (Lβ)ij. More...
 
template<typename DataType >
void longitudinal_operator_flat_cartesian (gsl::not_null< tnsr::II< DataType, 3 > * > result, const tnsr::ii< DataType, 3 > &strain)
 The conformal longitudinal operator (Lβ)ij on a flat conformal metric in Cartesian coordinates γij=δij. More...
 

Detailed Description

Items related to solving the Extended Conformal Thin Sandwich (XCTS) decomposition of the Einstein constraint equations.

The XCTS equations

(1)D¯2ψ18ψR¯112ψ5K2+18ψ7A¯ijA¯ij=2πψ5ρ(2)D¯i(L¯β)ij(L¯β)ijD¯iln(α¯)=α¯D¯i(α¯1u¯ij)+43α¯ψ6D¯jK+16πα¯ψ10Sj(3)D¯2(αψ)=αψ(78ψ8A¯ijA¯ij+512ψ4K2+18R¯+2πψ4(ρ+2S))ψ5tK+ψ5βiD¯iK(4)withA¯=12α¯((L¯β)iju¯ij)(5)andα¯=αψ6

are a set of nonlinear elliptic equations that the spacetime metric in general relativity must satisfy at all times. For an introduction see e.g. , in particular Box 3.3 which is largely mirrored here. We solve the XCTS equations for the conformal factor ψ, the product of lapse times conformal factor αψ and the shift vector βj. The remaining quantities in the equations, i.e. the conformal metric γ¯ij, the trace of the extrinsic curvature K, their respective time derivatives u¯ij and tK, the energy density ρ, the stress-energy trace S and the momentum density Si, are freely specifyable fields that define the physical scenario at hand. Of particular importance is the conformal metric, which defines the background geometry, the covariant derivative D¯, the Ricci scalar R¯ and the longitudinal operator

(6)(L¯β)ij=D¯iβj+D¯jβi23γ¯ijD¯kβk.

Note that the XCTS equations are essentially two Poisson equations and one Elasticity equation with nonlinear sources on a curved geometry. In this analogy, the longitudinal operator plays the role of the elastic constitutive relation that connects the symmetric "shift strain" D¯(iβj) with the "stress" (L¯β)ij of which we take the divergence in the momentum constraint. This particular constitutive relation is equivalent to an isotropic and homogeneous material with bulk modulus K=0 (not to be confused with the extrinsic curvature trace K in this context) and shear modulus μ=1 (see Elasticity::ConstitutiveRelations::IsotropicHomogeneous).

Once the XCTS equations are solved we can construct the spatial metric and extrinsic curvature as

(7)γij=ψ4γ¯ij(8)Kij=ψ2A¯ij+13γijK

from which we can compose the full spacetime metric.

Enumeration Type Documentation

◆ Equations

enum class Xcts::Equations
strong

Indicates a subset of the XCTS equations.

Enumerator
Hamiltonian 

Only the Hamiltonian constraint, solved for ψ.

HamiltonianAndLapse 

Both the Hamiltonian constraint and the lapse equation, solved for ψ and αψ.

HamiltonianLapseAndShift 

The full XCTS equations, solved for ψ, αψ and βexcess.

◆ Geometry

enum class Xcts::Geometry
strong

Types of conformal geometries for the XCTS equations.

Enumerator
FlatCartesian 

Euclidean (flat) manifold with Cartesian coordinates, i.e. the conformal metric has components γ¯ij=δij in these coordinates and thus all Christoffel symbols vanish: Γ¯jki=0.

Curved 

The conformal geometry is either curved or employs curved coordinates, so non-vanishing Christoffel symbols must be taken into account.

Function Documentation

◆ add_curved_hamiltonian_or_lapse_sources()

void Xcts::add_curved_hamiltonian_or_lapse_sources ( gsl::not_null< Scalar< DataVector > * >  hamiltonian_or_lapse_equation,
const Scalar< DataVector > &  conformal_ricci_scalar,
const Scalar< DataVector > &  field,
double  add_constant = 0. 
)

Add the contributions from a curved background geometry to the Hamiltonian constraint or lapse equation.

Adds 18ψR¯. This term appears both in the Hamiltonian constraint and the lapse equation (where in the latter ψ is replaced by αψ).

This term is linear.

See also
Xcts

◆ add_curved_momentum_sources()

template<int ConformalMatterScale>
void Xcts::add_curved_momentum_sources ( gsl::not_null< Scalar< DataVector > * >  hamiltonian_constraint,
gsl::not_null< Scalar< DataVector > * >  lapse_equation,
gsl::not_null< tnsr::I< DataVector, 3 > * >  momentum_constraint,
const tnsr::I< DataVector, 3 > &  conformal_momentum_density,
const tnsr::i< DataVector, 3 > &  extrinsic_curvature_trace_gradient,
const tnsr::ii< DataVector, 3 > &  conformal_metric,
const tnsr::II< DataVector, 3 > &  inv_conformal_metric,
const tnsr::I< DataVector, 3 > &  minus_div_dt_conformal_metric,
const Scalar< DataVector > &  conformal_factor_minus_one,
const Scalar< DataVector > &  lapse_times_conformal_factor_minus_one,
const tnsr::I< DataVector, 3 > &  conformal_factor_flux,
const tnsr::I< DataVector, 3 > &  lapse_times_conformal_factor_flux,
const tnsr::II< DataVector, 3 > &  longitudinal_shift_minus_dt_conformal_metric 
)

Add the nonlinear source to the momentum constraint and add the "distortion" source term to the Hamiltonian constraint and lapse equation.

Adds ((L¯β)iju¯ij)(j(αψ)αψ7jψψ)+ju¯ij+43αψψγ¯ijjK+16π(αψ)ψ3nS¯i to the momentum constraint, where n is the ConformalMatterScale and S¯i=ψnSi is the conformally-scaled momentum density.

Note that the ju¯ij term is not the full covariant divergence, but only the partial-derivatives part of it. The curved contribution to this term can be added together with the curved contribution to the flux divergence of the dynamic shift variable with the Elasticity::add_curved_sources function.

Also adds 18ψ7A¯ijA¯ij to the Hamiltonian constraint and 78αψ7A¯ijA¯ij to the lapse equation.

See also
Xcts

◆ add_distortion_hamiltonian_and_lapse_sources()

void Xcts::add_distortion_hamiltonian_and_lapse_sources ( gsl::not_null< Scalar< DataVector > * >  hamiltonian_constraint,
gsl::not_null< Scalar< DataVector > * >  lapse_equation,
const Scalar< DataVector > &  longitudinal_shift_minus_dt_conformal_metric_square,
const Scalar< DataVector > &  conformal_factor_minus_one,
const Scalar< DataVector > &  lapse_times_conformal_factor_minus_one 
)

Add the "distortion" source term to the Hamiltonian constraint and the lapse equation.

Adds 18ψ7A¯ijA¯ij to the Hamiltonian constraint and 78αψ7A¯ijA¯ij to the lapse equation.

See also
Xcts
Xcts::Tags::LongitudinalShiftMinusDtConformalMetricSquare

◆ add_distortion_hamiltonian_sources()

void Xcts::add_distortion_hamiltonian_sources ( gsl::not_null< Scalar< DataVector > * >  hamiltonian_constraint,
const Scalar< DataVector > &  longitudinal_shift_minus_dt_conformal_metric_over_lapse_square,
const Scalar< DataVector > &  conformal_factor_minus_one 
)

Add the "distortion" source term to the Hamiltonian constraint.

Adds 18ψ7A¯ijA¯ij.

See also
Xcts
Xcts::Tags::LongitudinalShiftMinusDtConformalMetricOverLapseSquare

◆ add_flat_cartesian_momentum_sources()

template<int ConformalMatterScale>
void Xcts::add_flat_cartesian_momentum_sources ( gsl::not_null< Scalar< DataVector > * >  hamiltonian_constraint,
gsl::not_null< Scalar< DataVector > * >  lapse_equation,
gsl::not_null< tnsr::I< DataVector, 3 > * >  momentum_constraint,
const tnsr::I< DataVector, 3 > &  conformal_momentum_density,
const tnsr::i< DataVector, 3 > &  extrinsic_curvature_trace_gradient,
const tnsr::I< DataVector, 3 > &  minus_div_dt_conformal_metric,
const Scalar< DataVector > &  conformal_factor_minus_one,
const Scalar< DataVector > &  lapse_times_conformal_factor_minus_one,
const tnsr::I< DataVector, 3 > &  conformal_factor_flux,
const tnsr::I< DataVector, 3 > &  lapse_times_conformal_factor_flux,
const tnsr::II< DataVector, 3 > &  longitudinal_shift_minus_dt_conformal_metric 
)

Add the nonlinear source to the momentum constraint and add the "distortion" source term to the Hamiltonian constraint and lapse equation.

Adds ((L¯β)iju¯ij)(j(αψ)αψ7jψψ)+ju¯ij+43αψψγ¯ijjK+16π(αψ)ψ3nS¯i to the momentum constraint, where n is the ConformalMatterScale and S¯i=ψnSi is the conformally-scaled momentum density.

Note that the ju¯ij term is not the full covariant divergence, but only the partial-derivatives part of it. The curved contribution to this term can be added together with the curved contribution to the flux divergence of the dynamic shift variable with the Elasticity::add_curved_sources function.

Also adds 18ψ7A¯ijA¯ij to the Hamiltonian constraint and 78αψ7A¯ijA¯ij to the lapse equation.

See also
Xcts

◆ add_hamiltonian_sources()

template<int ConformalMatterScale>
void Xcts::add_hamiltonian_sources ( gsl::not_null< Scalar< DataVector > * >  hamiltonian_constraint,
const Scalar< DataVector > &  conformal_energy_density,
const Scalar< DataVector > &  extrinsic_curvature_trace,
const Scalar< DataVector > &  conformal_factor_minus_one 
)

Add the nonlinear source to the Hamiltonian constraint on a flat conformal background in Cartesian coordinates and with u¯ij=0=βi.

Adds 112ψ5K22πψ5nρ¯ where n is the ConformalMatterScale and ρ¯=ψnρ is the conformally-scaled energy density. Additional sources can be added with add_distortion_hamiltonian_sources and add_curved_hamiltonian_or_lapse_sources.

See also
Xcts

◆ add_lapse_sources()

template<int ConformalMatterScale>
void Xcts::add_lapse_sources ( gsl::not_null< Scalar< DataVector > * >  lapse_equation,
const Scalar< DataVector > &  conformal_energy_density,
const Scalar< DataVector > &  conformal_stress_trace,
const Scalar< DataVector > &  extrinsic_curvature_trace,
const Scalar< DataVector > &  dt_extrinsic_curvature_trace,
const Scalar< DataVector > &  shift_dot_deriv_extrinsic_curvature_trace,
const Scalar< DataVector > &  conformal_factor_minus_one,
const Scalar< DataVector > &  lapse_times_conformal_factor_minus_one 
)

Add the nonlinear source to the lapse equation on a flat conformal background in Cartesian coordinates and with u¯ij=0=βi.

Adds (αψ)(512ψ4K2+2πψ4n(ρ¯+2S¯))+ψ5(βiiKtK) where n is the ConformalMatterScale and matter quantities are conformally-scaled. Additional sources can be added with add_distortion_hamiltonian_and_lapse_sources and add_curved_hamiltonian_or_lapse_sources.

See also
Xcts

◆ add_linearized_distortion_hamiltonian_and_lapse_sources()

void Xcts::add_linearized_distortion_hamiltonian_and_lapse_sources ( gsl::not_null< Scalar< DataVector > * >  hamiltonian_constraint,
gsl::not_null< Scalar< DataVector > * >  lapse_equation,
const Scalar< DataVector > &  longitudinal_shift_minus_dt_conformal_metric_square,
const Scalar< DataVector > &  conformal_factor_minus_one,
const Scalar< DataVector > &  lapse_times_conformal_factor_minus_one,
const Scalar< DataVector > &  conformal_factor_correction,
const Scalar< DataVector > &  lapse_times_conformal_factor_correction 
)

The linearization of add_distortion_hamiltonian_and_lapse_sources

Note that this linearization is only w.r.t. ψ and αψ.

◆ add_linearized_distortion_hamiltonian_sources()

void Xcts::add_linearized_distortion_hamiltonian_sources ( gsl::not_null< Scalar< DataVector > * >  linearized_hamiltonian_constraint,
const Scalar< DataVector > &  longitudinal_shift_minus_dt_conformal_metric_over_lapse_square,
const Scalar< DataVector > &  conformal_factor_minus_one,
const Scalar< DataVector > &  conformal_factor_correction 
)

The linearization of add_distortion_hamiltonian_sources

Note that this linearization is only w.r.t. ψ.

◆ add_linearized_lapse_sources()

template<int ConformalMatterScale>
void Xcts::add_linearized_lapse_sources ( gsl::not_null< Scalar< DataVector > * >  linearized_lapse_equation,
const Scalar< DataVector > &  conformal_energy_density,
const Scalar< DataVector > &  conformal_stress_trace,
const Scalar< DataVector > &  extrinsic_curvature_trace,
const Scalar< DataVector > &  dt_extrinsic_curvature_trace,
const Scalar< DataVector > &  shift_dot_deriv_extrinsic_curvature_trace,
const Scalar< DataVector > &  conformal_factor_minus_one,
const Scalar< DataVector > &  lapse_times_conformal_factor_minus_one,
const Scalar< DataVector > &  conformal_factor_correction,
const Scalar< DataVector > &  lapse_times_conformal_factor_correction 
)

The linearization of add_lapse_sources

Note that this linearization is only w.r.t. ψ and αψ. The linearization w.r.t. βi is added in add_curved_linearized_momentum_sources or add_flat_cartesian_linearized_momentum_sources.

◆ adm_linear_momentum_surface_integrand()

void Xcts::adm_linear_momentum_surface_integrand ( gsl::not_null< tnsr::II< DataVector, 3 > * >  result,
const Scalar< DataVector > &  conformal_factor,
const tnsr::II< DataVector, 3 > &  inv_spatial_metric,
const tnsr::II< DataVector, 3 > &  inv_extrinsic_curvature,
const Scalar< DataVector > &  trace_extrinsic_curvature 
)

Surface integrand for the ADM linear momentum calculation.

We define the ADM linear momentum integral as (see Eqs. 19-20 in ):

(9)PADMi=18πSψ10(KijKγij)dSj.

Note
For consistency with adm_linear_momentum_volume_integrand, this integrand needs to be contracted with the Euclidean face normal and integrated with the Euclidean area element.
Parameters
resultoutput pointer
conformal_factorthe conformal factor ψ
inv_spatial_metricthe inverse spatial metric γij
inv_extrinsic_curvaturethe inverse extrinsic curvature Kij
trace_extrinsic_curvaturethe trace of the extrinsic curvature K

◆ adm_linear_momentum_volume_integrand()

void Xcts::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 > &  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 
)

Volume integrand for ADM linear momentum calculation defined as (see Eq. 20 in ):

(10)PADMi=18πV(Γ¯jkiPjk+Γ¯jkjPjk2γ¯jkPjkγ¯ill(lnψ))dV,

where 1/(8π)Pjk is the result from adm_linear_momentum_surface_integrand.

Note
For consistency with adm_linear_momentum_surface_integrand, this integrand needs to be integrated with the Euclidean volume element.
Parameters
resultoutput pointer
surface_integrandthe quantity 1/(8π)Pij (result of adm_linear_momentum_surface_integrand)
conformal_factorthe conformal factor ψ
deriv_conformal_factorthe gradient of the conformal factor iψ
conformal_metricthe conformal metric γ¯ij
inv_conformal_metricthe inverse conformal metric γ¯ij
conformal_christoffel_second_kindthe conformal christoffel symbol Γ¯jki
conformal_christoffel_contractedthe contracted conformal christoffel symbol Γ¯i

◆ adm_mass_surface_integrand()

void Xcts::adm_mass_surface_integrand ( gsl::not_null< tnsr::I< DataVector, 3 > * >  result,
const tnsr::i< DataVector, 3 > &  deriv_conformal_factor,
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 
)

Surface integrand for the ADM mass calculation.

We define the ADM mass integral as (see Eq. 3.139 in ):

(11)MADM=116πS(γ¯jkΓ¯jkiγ¯ijΓ¯j8γ¯ijjψ)dS¯i.

Note
We don't use the other versions presented in of this integral because they make assumptions like γ¯=1, Γ¯iji=0 and faster fall-off of the conformal metric.
For consistency with adm_mass_volume_integrand, this integrand needs to be contracted with the conformal face normal and integrated with the conformal area element.
Parameters
resultoutput pointer
deriv_conformal_factorthe gradient of the conformal factor iψ
inv_conformal_metricthe inverse conformal metric γ¯ij
conformal_christoffel_second_kindthe conformal christoffel symbol Γ¯jki
conformal_christoffel_contractedthe conformal christoffel symbol contracted in its first two indices Γ¯i=Γ¯ijj

◆ adm_mass_volume_integrand()

void Xcts::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 Scalar< DataVector > &  energy_density,
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::i< DataVector, 3 > &  conformal_christoffel_contracted,
const tnsr::iJkk< DataVector, 3 > &  deriv_conformal_christoffel_second_kind 
)

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:

(12)MADM=116πV(iγ¯jkΓ¯jki+γ¯jkiΓ¯jki+Γ¯lγ¯jkΓ¯jkliγ¯ijΓ¯jγ¯ijiΓ¯jΓ¯lγ¯ljΓ¯j8D¯2ψ)dV¯,

where we can use the Hamiltonian constraint (Eq. 3.37 in ) to replace 8D¯2ψ with

(13)8D¯2ψ=ψR¯+23ψ5K214ψ51α2[(L¯β)iju¯ij][(L¯β)iju¯ij]16πψ5ρ.

Note
This is similar to Eq. 3.149 in , except that here we don't assume γ¯=1.
For consistency with adm_mass_surface_integrand, this integrand needs to be integrated with the conformal volume element.
Parameters
resultoutput pointer
conformal_factorthe conformal factor
conformal_ricci_scalarthe conformal Ricci scalar R¯
trace_extrinsic_curvaturethe extrinsic curvature trace K
longitudinal_shift_minus_dt_conformal_metric_over_lapse_squarethe quantity computed in Xcts::Tags::LongitudinalShiftMinusDtConformalMetricOverLapseSquare
energy_densitythe energy density ρ
inv_conformal_metricthe inverse conformal metric γ¯ij
deriv_inv_conformal_metricthe gradient of the inverse conformal metric iγ¯jk
conformal_christoffel_second_kindthe conformal christoffel symbol Γ¯jki
conformal_christoffel_contractedthe conformal christoffel symbol contracted in its first two indices Γ¯i=Γ¯ijj
deriv_conformal_christoffel_second_kindthe gradient of the conformal christoffel symbol iΓ¯klj

◆ center_of_mass_surface_integrand()

void Xcts::center_of_mass_surface_integrand ( gsl::not_null< tnsr::I< DataVector, 3 > * >  result,
const Scalar< DataVector > &  conformal_factor,
const tnsr::I< DataVector, 3 > &  coords 
)

Surface integrand for the center of mass calculation.

We define the center of mass integral as

(14)CCoMi=38πMADMS(ψ41)nidA,

where ni=xi/r and r=x2+y2+z2.

Analytically, this is identical to the definition in Eq. (25) of because

(15)SnidA=0.

Numerically, we have found that subtracting 1 from ψ4 results in less round-off errors, leading to a more accurate center of mass.

One way to interpret this integral is that we are summing over the unit vectors ni, rescaled by ψ4, in all directions. If ψ(r) is constant, no rescaling happens and all the unit vectors cancel out. If ψ(r) is not constant, then CCoM will emerge from the difference of large numbers (sum of rescaled ni in each subdomain). With larger and larger numbers being involved in this cancellation (i.e. with increasing radius of S), we loose numerical accuracy. In other words, we are seeking the subdominant terms. Since ψ41 in conformal flatness, subtracting 1 from it in the integrand makes the numbers involved in this cancellation smaller, reducing this issue. We have also tried different variations of this integrand with the same motivation, but (ψ41) is the best one when taking simplicity and accuracy gain into consideration.

Note
We don't include the ADM mass MADM in this integrand. After integrating the result of this function, you have to divide by MADM.
For consistency with center_of_mass_volume_integrand, this integrand needs to be integrated with the Euclidean area element.
See also
Xcts::adm_mass_surface_integrand
Warning
This integral assumes that the conformal metric falls off to flatness faster than 1/r2. That means that it cannot be directly used with the Kerr-Schild metric, which falls off as 1/r. This is not a problem for XCTS with Superposed Kerr-Schild (SKS) because of the exponential fall-off terms.
Parameters
resultoutput pointer
conformal_factorthe conformal factor ψ
coordsthe inertial coordinates xi

◆ center_of_mass_volume_integrand()

void Xcts::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 > &  coords 
)

Volume integrand for the center of mass calculation.

We cast the center of mass as an infinite volume integral by applying Gauss' law on the surface integral defined in center_of_mass_surface_integrand:

(16)CCoMi=38πMADMV(4ψ3jψninj+2rψ4ni)dV=34πMADMV1r2(2ψ3jψxixj+ψ4xi)dV,

where ni=xi/r and r=x2+y2+z2.

Note
For consistency with center_of_mass_surface_integrand, this integrand needs to be integrated with the Euclidean volume element.
See also
center_of_mass_surface_integrand
Parameters
resultoutput pointer
conformal_factorthe conformal factor ψ
deriv_conformal_factorthe gradient of the conformal factor iψ
coordsthe inertial coordinates xi

◆ extrinsic_curvature()

template<typename DataType >
void Xcts::extrinsic_curvature ( const gsl::not_null< tnsr::ii< DataType, 3 > * >  result,
const Scalar< DataType > &  conformal_factor,
const Scalar< DataType > &  lapse,
const tnsr::ii< DataType, 3 > &  conformal_metric,
const tnsr::II< DataType, 3 > &  longitudinal_shift_minus_dt_conformal_metric,
const Scalar< DataType > &  trace_extrinsic_curvature 
)

Extrinsic curvature computed from the conformal decomposition used in the XCTS system.

The extrinsic curvature decomposition is (see Eq. 3.113 in ):

Undefined control sequence \lapse

Parameters
resultoutput buffer for the extrinsic curvature
conformal_factorthe conformal factor ψ
lapsethe lapse α
conformal_metricthe conformal metric γ¯ij
longitudinal_shift_minus_dt_conformal_metricthe term (L¯β)iju¯ij. Note that (L¯β)ij is the conformal longitudinal shift, and (L¯β)ij=ψ4(Lβ)ij (Eq. 3.98 in ). See also Xcts::longitudinal_operator.
trace_extrinsic_curvaturethe trace of the extrinsic curvature, K=γijKij. Note that it is a conformal invariant, K=K¯ (by choice).

◆ longitudinal_operator() [1/2]

template<typename DataType >
void Xcts::longitudinal_operator ( gsl::not_null< tnsr::II< DataType, 3 > * >  result,
const tnsr::I< DataType, 3 > &  shift,
const tnsr::iJ< DataType, 3 > &  deriv_shift,
const tnsr::II< DataType, 3 > &  inv_metric,
const tnsr::Ijj< DataType, 3 > &  christoffel_second_kind 
)

The longitudinal operator, or vector gradient, (Lβ)ij.

Computes the longitudinal operator

(17)(Lβ)ij=iβj+jβi23γijkβk

of a vector field βi, where denotes the covariant derivative w.r.t. the metric γ (see e.g. Eq. (3.50) in ). Note that in the XCTS equations the longitudinal operator is typically applied to conformal quantities and w.r.t. the conformal metric γ¯.

In terms of the symmetric "strain" quantity Bij=(iγj)kβk the longitudinal operator is:

(18)(Lβ)ij=2(γikγjl13γijγkl)Bkl

Note that the strain can be computed with Elasticity::strain.

◆ longitudinal_operator() [2/2]

template<typename DataType >
void Xcts::longitudinal_operator ( gsl::not_null< tnsr::II< DataType, 3 > * >  result,
const tnsr::ii< DataType, 3 > &  strain,
const tnsr::II< DataType, 3 > &  inv_metric 
)

The longitudinal operator, or vector gradient, (Lβ)ij.

Computes the longitudinal operator

(19)(Lβ)ij=iβj+jβi23γijkβk

of a vector field βi, where denotes the covariant derivative w.r.t. the metric γ (see e.g. Eq. (3.50) in ). Note that in the XCTS equations the longitudinal operator is typically applied to conformal quantities and w.r.t. the conformal metric γ¯.

In terms of the symmetric "strain" quantity Bij=(iγj)kβk the longitudinal operator is:

(20)(Lβ)ij=2(γikγjl13γijγkl)Bkl

Note that the strain can be computed with Elasticity::strain.

◆ longitudinal_operator_flat_cartesian()

template<typename DataType >
void Xcts::longitudinal_operator_flat_cartesian ( gsl::not_null< tnsr::II< DataType, 3 > * >  result,
const tnsr::ii< DataType, 3 > &  strain 
)

The conformal longitudinal operator (Lβ)ij on a flat conformal metric in Cartesian coordinates γij=δij.

See also
Xcts::longitudinal_operator