19#include "excision_surf.h"
20#include "excision_hor.h"
135 Isol_hole(
const Map& mp_i,
double Omega_i,
bool NorKappa_i,
Scalar NoK_i,
bool isCF_i =
false) ;
154 Isol_hole(
const Map& mp_i,
double Omega_i,
bool NorKappa_i,
Scalar NoK_i,
bool isCF_i, FILE* fich) ;
181 void compute_stat_metric(
double precis,
double relax,
int mer_max,
int mer_max2,
bool isvoid =
true) ;
206 else cout <<
"The boundary condition is imposed on the surface gravity!" << endl;
214 else cout <<
"The boundary condition is imposed on the lapse!" <<endl;
239 virtual void sauve(FILE* )
const ;
Class to compute quasistationary single black hole spacetimes in vacuum.
Spheroid hor()
Spheroid at the excised boundary associated with the black hole MOTS on the slice.
double * p_adm_mass
Computation of the ADM mass of the BH spacetime.
double * p_komar_angmom
Computation of the Komar angular momentum w.r.t.
double get_Omega() const
Returns the rotation rate.
const Map & get_mp() const
Returns the mapping.
virtual ~Isol_hole()
Destructor.
void Einstein_errors()
Prints out errors in Einstein equations for the data obtained.
Scalar boundNoK
Indicates if the boundary value for the lapse or the surface gravity is used.
virtual void del_deriv() const
Deletes all the derived quantities.
const Sym_tensor & get_hatA() const
Returns the rescaled tracefree extrinsic curvature .
Scalar lapse
Lapse function.
bool NorKappa
Rotation rate of the horizon in the azimuthal direction.
Sym_tensor hij
Deviation tensor( non-flat part of the conformal 3-metric on the slice; see Bonazzola et al.
const Sym_tensor & get_hij() const
Returns the deviation tensor .
double * p_virial_residue
Computation of the Virial residual, as difference at infinity between the ADM mass and the Komar inte...
void operator=(const Isol_hole &)
Assignment to another Isol_hole.
const Scalar & get_boundN() const
Returns the boundary value used for the lapse (if it is the one used)
void secmembre_kerr(Sym_tensor &source_hh)
Computes the rhs of hyperbolic equation for conformal metric assuming statioarity; WARNING; up to now...
Vector shift
Shift vector.
virtual void sauve(FILE *) const
Save in a file.
Sym_tensor hatA
Rescaled tracefree extrinsic curvature tensor: see Cordero et al.
const Map & mp
Mapping associated with the star.
virtual void set_der_0x0() const
Sets to 0x0 all the pointers on derived quantities.
const Scalar & get_lapse() const
Returns the lapse function N.
bool isCF
Stores the boundary value of the lapse or surface gravity.
const Scalar & get_conf_fact() const
Returns the conformal factor.
double adm_mass()
Computation of the ADM mass of the BH spacetime.
Spheroid * p_hor
Computation of the spheroid associated with the black hole horizon.
void compute_stat_metric(double precis, double relax, int mer_max, int mer_max2, bool isvoid=true)
Computes a quasi-stationary 3-slice from the chosen parameters.
double komar_angmom()
Computation of the Komar angular momentum w.r.t.
const Vector & get_shift() const
Returns the shift vector .
const Scalar & get_Kappa() const
Returns the surface gravity value at the boundary (if it is the one used)
double virial_residue()
Computation of the Virial residual, as difference at infinity between the ADM mass and the Komar inte...
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.