SpECTRE
v2024.09.29
|
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 | |
struct | LinearizedSources |
The linearization of the sources | |
struct | Sources |
The sources | |
struct | SpacetimeQuantitiesComputer |
CachedTempBuffer computer class for 3+1 quantities from XCTS variables. See Xcts::SpacetimeQuantities . More... | |
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 | |
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 | |
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 [150]): 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< 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 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, | |
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, | |
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 | |
Items related to solving the Extended Conformal Thin Sandwich (XCTS) decomposition of the Einstein constraint equations.
The XCTS equations
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. [13], in particular Box 3.3 which is largely mirrored here. We solve the XCTS equations for the conformal factor
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" Elasticity::ConstitutiveRelations::IsotropicHomogeneous
).
Once the XCTS equations are solved we can construct the spatial metric and extrinsic curvature as
from which we can compose the full spacetime metric.
|
strong |
|
strong |
Types of conformal geometries for the XCTS equations.
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
This term is linear.
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 ConformalMatterScale
and
Note that the Elasticity::add_curved_sources
function.
Also adds
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
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
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 ConformalMatterScale
and
Note that the Elasticity::add_curved_sources
function.
Also adds
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
Adds ConformalMatterScale
and add_distortion_hamiltonian_sources
and add_curved_hamiltonian_or_lapse_sources
.
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
Adds 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
.
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.
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.
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. add_curved_linearized_momentum_sources
or add_flat_cartesian_linearized_momentum_sources
.
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 [150]):
adm_linear_momentum_volume_integrand
, this integrand needs to be contracted with the Euclidean face normal and integrated with the Euclidean area element.result | output pointer |
conformal_factor | the conformal factor |
inv_spatial_metric | the inverse spatial metric |
inv_extrinsic_curvature | the inverse extrinsic curvature |
trace_extrinsic_curvature | the trace of the extrinsic curvature |
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 [150]):
where adm_linear_momentum_surface_integrand
.
adm_linear_momentum_surface_integrand
, this integrand needs to be integrated with the Euclidean volume element.result | output pointer |
surface_integrand | the quantity adm_linear_momentum_surface_integrand ) |
conformal_factor | the conformal factor |
deriv_conformal_factor | the gradient of the conformal factor |
conformal_metric | the conformal metric |
inv_conformal_metric | the inverse conformal metric |
conformal_christoffel_second_kind | the conformal christoffel symbol |
conformal_christoffel_contracted | the contracted conformal christoffel symbol |
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 [13]):
adm_mass_volume_integrand
, this integrand needs to be contracted with the conformal face normal and integrated with the conformal area element.result | output pointer |
deriv_conformal_factor | the gradient of the conformal factor |
inv_conformal_metric | the inverse conformal metric |
conformal_christoffel_second_kind | the conformal christoffel symbol |
conformal_christoffel_contracted | the conformal christoffel symbol contracted in its first two indices |
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
:
where we can use the Hamiltonian constraint (Eq. 3.37 in [13]) to replace
adm_mass_surface_integrand
, this integrand needs to be integrated with the conformal volume element.result | output pointer |
conformal_factor | the conformal factor |
conformal_ricci_scalar | the conformal Ricci scalar |
trace_extrinsic_curvature | the extrinsic curvature trace |
longitudinal_shift_minus_dt_conformal_metric_over_lapse_square | the quantity computed in Xcts::Tags::LongitudinalShiftMinusDtConformalMetricOverLapseSquare |
energy_density | the energy density |
inv_conformal_metric | the inverse conformal metric |
deriv_inv_conformal_metric | the gradient of the inverse conformal metric |
conformal_christoffel_second_kind | the conformal christoffel symbol |
conformal_christoffel_contracted | the conformal christoffel symbol contracted in its first two indices |
deriv_conformal_christoffel_second_kind | the gradient of the conformal christoffel symbol |
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
where
Analytically, this is identical to the definition in Eq. (25) of [150] because
Numerically, we have found that subtracting
One way to interpret this integral is that we are summing over the unit vectors
center_of_mass_volume_integrand
, this integrand needs to be integrated with the Euclidean area element.Xcts::adm_mass_surface_integrand
result | output pointer |
conformal_factor | the conformal factor |
coords | the inertial coordinates |
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
:
where
center_of_mass_surface_integrand
, this integrand needs to be integrated with the Euclidean volume element.center_of_mass_surface_integrand
result | output pointer |
conformal_factor | the conformal factor |
deriv_conformal_factor | the gradient of the conformal factor |
coords | the inertial coordinates |
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 [13]):
result | output buffer for the extrinsic curvature |
conformal_factor | the conformal factor |
lapse | the lapse |
conformal_metric | the conformal metric |
longitudinal_shift_minus_dt_conformal_metric | the term |
trace_extrinsic_curvature | the trace of the extrinsic curvature, |
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,
Computes the longitudinal operator
of a vector field
In terms of the symmetric "strain" quantity
Note that the strain can be computed with Elasticity::strain
.
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,
Computes the longitudinal operator
of a vector field
In terms of the symmetric "strain" quantity
Note that the strain can be computed with Elasticity::strain
.
void Xcts::longitudinal_operator_flat_cartesian | ( | gsl::not_null< tnsr::II< DataType, 3 > * > | result, |
const tnsr::ii< DataType, 3 > & | strain | ||
) |
The conformal longitudinal operator
Xcts::longitudinal_operator