11#ifndef __EXCISEDSLICE_H_
12#define __EXCISEDSLICE_H_
183 void compute_stat_metric(
double precis,
double Omega_i,
bool NorKappa_i,
Scalar NoK_i,
bool isCF_i=
true,
double relax=0.2,
int mer_max = 2000,
int mer_max2=200,
bool isvoid =
true) ;
224 else if (
field_set==1) cout <<
"Error: the scheme used is the original FCF; this variable is irrelevant" << endl;
225 else cout <<
"error in the scheme definition; please check the class consistency" << endl;}
252 virtual void sauve(FILE* )
const ;
Class to compute single black hole spacetime excised slices.
Scalar & set_conf_fact()
Sets the conformal factor.
const Scalar & get_conf_fact() const
Returns the conformal factor.
Sym_tensor hatA
Rescaled tracefree extrinsic curvature tensor: rescaled same way as Cordero et al.
double * p_adm_mass
Computation of the ADM mass of the BH spacetime.
double komar_angmom()
Computation of the Komar angular momentum w.r.t.
void compute_stat_metric(double precis, double Omega_i, bool NorKappa_i, Scalar NoK_i, bool isCF_i=true, double relax=0.2, int mer_max=2000, int mer_max2=200, bool isvoid=true)
If hor_type=1, computes a quasi-stationary 3-slice from the chosen parameters.
int type_hor
Chosen horizon type.
void Einstein_errors()
Prints out errors in Einstein equations for the data obtained.
double * p_komar_angmom
Computation of the Komar angular momentum w.r.t.
virtual void del_deriv() const
Deletes all the derived quantities.
Sym_tensor & set_hij()
Sets the deviation tensor.
Vector & set_Xx()
Sets the longitudinal part of Einstein Equations if the scheme is apropriate.
Spheroid * p_hor
Computation of the spheroid associated with the black hole horizon.
const Scalar & get_lapse() const
Returns the lapse function N.
Scalar & set_lapse()
Sets the lapse.
const Vector & get_shift() const
Returns the shift vector .
double virial_residue()
Computation of the Virial residual, as difference at infinity between the ADM mass and the Komar inte...
void secmembre_kerr(Sym_tensor &source_hh)
If hor_type=1, computes the rhs of hyperbolic equation for conformal metric assuming stationarity; WA...
const Sym_tensor & get_hatA() const
Returns the rescaled tracefree extrinsic curvature (or its TT part if scheme_type=2)
int get_field_set() const
Returns the field set chosen for the data.
Sym_tensor hij
Deviation tensor( non-flat part of the conformal 3-metric on the slice; see Bonazzola et al.
const Map & get_mp() const
Returns the mapping.
virtual void set_der_0x0() const
Sets to 0x0 all the pointers on derived quantities.
int set_type_hor()
Sets the horizon type.
Sym_tensor & set_hatA()
Sets the rescaled tracefree extrinsic curvature (or its TT part if scheme_type=2)
Scalar lapse
Lapse function.
double * p_virial_residue
Computation of the Virial residual, as difference at infinity between the ADM mass and the Komar inte...
const Sym_tensor & get_hij() const
Returns the deviation tensor .
const Vector & get_Xx() const
Return the longitudinal part of Einstein Equations if the scheme is appropriate.
void operator=(const Excised_slice &)
Assignment to another Excised_slice.
Vector & set_shift()
Sets the shift vector.
virtual ~Excised_slice()
Destructor.
const Map & mp
Mapping associated with the slice.
int get_type_hor() const
Returns the type of horizon that performs the excision.
int field_set
Chose field set type.
Vector shift
Shift vector.
virtual void sauve(FILE *) const
Save in a file.
double adm_mass()
Computation of the ADM mass of the BH spacetime.
Vector Xx
Longitudinal part of the extrinsic curvature.
Spheroid hor()
Spheroid at the excised boundary associated with the black hole MOTS on the slice.
Base class for coordinate mappings.
Tensor field of valence 0 (or component of a tensorial field).
Spheroidal 2-surfaces embedded in a time-slice of the 3+1 formalism.
Class intended to describe valence-2 symmetric tensors.
Tensor field of valence 1.