Lorene::Excised_slice Class Reference

Class to compute single black hole spacetime excised slices. More...

#include <excised_slice.h>

List of all members.

Public Member Functions

 Excised_slice (const Map &mp_i, int hor_type, int scheme_type)
 Standard constructor.
 Excised_slice (const Excised_slice &)
 Copy constructor.
 Excised_slice (const Map &mp_i, int hor_type, int scheme_type, FILE *fich)
 Constructor from a file (see sauve(FILE* )).
virtual ~Excised_slice ()
 Destructor.
void operator= (const Excised_slice &)
 Assignment to another Excised_slice.
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.
void secmembre_kerr (Sym_tensor &source_hh)
 If hor_type=1, computes the rhs of hyperbolic equation for conformal metric assuming stationarity; WARNING; up to now, we are only able to handle void spacetimes.
const Mapget_mp () const
 Returns the mapping.
int get_type_hor () const
 Returns the type of horizon that performs the excision.
int get_field_set () const
 Returns the field set chosen for the data.
const Scalarget_lapse () const
 Returns the lapse function N.
const Scalarget_conf_fact () const
 Returns the conformal factor.
const Vectorget_shift () const
 Returns the shift vector $\beta^i$.
const Sym_tensorget_hij () const
 Returns the deviation tensor $ h^{ij} $.
const Sym_tensorget_hatA () const
 Returns the rescaled tracefree extrinsic curvature $\hat{A}^{ij}$ (or its TT part if scheme_type=2).
const Vectorget_Xx () const
 Return the longitudinal part of Einstein Equations if the scheme is appropriate.
int set_type_hor ()
 Sets the horizon type.
Scalarset_lapse ()
 Sets the lapse.
Scalarset_conf_fact ()
 Sets the conformal factor.
Vectorset_shift ()
 Sets the shift vector.
Sym_tensorset_hij ()
 Sets the deviation tensor.
Sym_tensorset_hatA ()
 Sets the rescaled tracefree extrinsic curvature $\hat{A}^{ij}$ (or its TT part if scheme_type=2).
Vectorset_Xx ()
 Sets the longitudinal part of Einstein Equations if the scheme is apropriate.
virtual void sauve (FILE *) const
 Save in a file.
void Einstein_errors ()
 Prints out errors in Einstein equations for the data obtained.
Spheroid hor ()
 Spheroid at the excised boundary associated with the black hole MOTS on the slice.
double adm_mass ()
 Computation of the ADM mass of the BH spacetime.
double komar_angmom ()
 Computation of the Komar angular momentum w.r.t.
double virial_residue ()
 Computation of the Virial residual, as difference at infinity between the ADM mass and the Komar integral associated to the mass.

Protected Member Functions

virtual void del_deriv () const
 Deletes all the derived quantities.
virtual void set_der_0x0 () const
 Sets to 0x0 all the pointers on derived quantities.

Protected Attributes

const Mapmp
 Mapping associated with the slice.
int type_hor
 Chosen horizon type.
int field_set
 Chose field set type.
Scalar lapse
 Lapse function.
Scalar conf_fact
Vector shift
 Shift vector.
Sym_tensor hij
 Deviation tensor( non-flat part of the conformal 3-metric on the slice; see Bonazzola et al.
Sym_tensor hatA
 Rescaled tracefree extrinsic curvature tensor: rescaled same way as Cordero et al.
Vector Xx
 Longitudinal part of the extrinsic curvature.
Spheroidp_hor
 Computation of the spheroid associated with the black hole horizon.
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 * p_virial_residue
 Computation of the Virial residual, as difference at infinity between the ADM mass and the Komar integral associated to the mass.

Detailed Description

Class to compute single black hole spacetime excised slices.

It takes as arguments:

Definition at line 41 of file excised_slice.h.


Constructor & Destructor Documentation

Lorene::Excised_slice::Excised_slice ( const Map mp_i,
int  hor_type,
int  scheme_type 
)

Standard constructor.

Parameters:
mp_i Mapping on which the black hole slice will be defined
hor_type Type of horizon used for excision:

  • 0. General
  • 1. Isolated horizon-like
  • 2. Dynamical horizon-like
  • 3. Inner trapped surface
scheme_type Type of scheme use for metric field computation

  • 1. Original FCF scheme (Bonazzola et al., 2004)
  • 2. New FCF scheme (Cordero et al., 2009)

Definition at line 58 of file excised_slice.C.

References hatA, hij, lapse, set_der_0x0(), Lorene::Tensor::set_etat_zero(), shift, and Xx.

Lorene::Excised_slice::Excised_slice ( const Excised_slice ih  ) 

Copy constructor.

Definition at line 84 of file excised_slice.C.

References set_der_0x0().

Lorene::Excised_slice::Excised_slice ( const Map mp_i,
int  hor_type,
int  scheme_type,
FILE *  fich 
)

Constructor from a file (see sauve(FILE* )).

Parameters:
mp_i Mapping on which the black hole slice will be defined
hor_type Type of horizon used for excision:

  • 0. General
  • 1. Isolated horizon-like
  • 2. Dynamical horizon-like
  • 3. Inner trapped surface
scheme_type Type of scheme use for metric field computation

  • 1. Original FCF scheme (Bonazzola et al., 2004)
  • 2. New FCF scheme (Cordero et al., 2009)
fich input file (must have been created by the function sauve)

Definition at line 101 of file excised_slice.C.

References set_der_0x0().

Lorene::Excised_slice::~Excised_slice (  )  [virtual]

Destructor.

Definition at line 124 of file excised_slice.C.

References del_deriv().


Member Function Documentation

double Lorene::Excised_slice::adm_mass (  ) 
void Lorene::Excised_slice::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.

Parameters:
precis [input] threshold in the relative difference between the radial shift for two consecutive steps to stop the iterative procedure (default value: 1.e-11)
Omega_i,: Rotation parameter for the isolated horizon slice. The principal Killing vector associated to this quantity will be by default the $e_{}$ vector.
NorKappa_i,: boolean to tell whether boundary condition is put on the surface gravity (true) or the lapse (false)
NoK_i,: The corresponding Dirichlet value for this quantity.
isCF,: Tells whether the calculation is using the conformal flatness approximation (true) or not (false).
relax relaxation coefficient for all metric fields (0: no relaxation)
mer_max,: Maximum number of iterations allowed in the main loop
mer_max,: Maximum number of iterations allowed in the possible secondary loop for conformal geometry computation.
isvoid,: indicates whether or not the spacetime is vacuum. This is useless for now since matter terms are not handled in the class.

Definition at line 38 of file Excised_slice/isol_hole_compute_metric.C.

References Lorene::Scalar::allocate_all(), Lorene::Tensor::annule(), Lorene::Tensor::annule_domain(), Lorene::Scalar::annule_hard(), Lorene::Valeur::c_cf, Lorene::Valeur::coef(), Lorene::Metric::con(), Lorene::Metric_flat::con(), Lorene::contract(), Lorene::Metric::cov(), Lorene::Scalar::dec_dzpuis(), Lorene::Tensor::dec_dzpuis(), Lorene::Scalar::derive_con(), Lorene::Tensor::derive_con(), Lorene::Tensor_sym::derive_con(), Lorene::Tensor::derive_cov(), Lorene::Scalar::derive_cov(), Lorene::Tensor_sym::derive_cov(), Lorene::Sym_tensor::derive_lie(), Lorene::Metric::determinant(), Lorene::Sym_tensor::divergence(), Lorene::Vector::divergence(), Lorene::Scalar::dsdr(), Lorene::Vector::eta(), Lorene::Scalar::get_etat(), Lorene::Scalar::get_spectral_va(), hatA, hij, Lorene::Tensor::inc_dzpuis(), Lorene::Scalar::inc_dzpuis(), Lorene::Scalar::laplacian(), lapse, Lorene::log(), Lorene::max(), Lorene::maxabs(), mp, Lorene::Vector::mu(), Lorene::Scalar::mult_rsint(), secmembre_kerr(), Lorene::Vector::set(), Lorene::Tensor::set(), Lorene::Scalar::set_dzpuis(), Lorene::Scalar::set_etat_qcq(), Lorene::Scalar::set_grid_point(), Lorene::Scalar::set_spectral_va(), shift, Lorene::Scalar::sol_elliptic_boundary(), Lorene::sqrt(), Lorene::Vector::std_spectral_base(), Lorene::Scalar::std_spectral_base(), Lorene::Tensor::trace(), type_hor, Lorene::Tensor::up(), Lorene::Tensor::up_down(), Lorene::Scalar::val_grid_point(), and Lorene::Valeur::ylm().

void Lorene::Excised_slice::del_deriv (  )  const [protected, virtual]

Deletes all the derived quantities.

Definition at line 135 of file excised_slice.C.

References p_adm_mass, p_hor, p_komar_angmom, p_virial_residue, and set_der_0x0().

void Lorene::Excised_slice::Einstein_errors (  ) 
const Scalar& Lorene::Excised_slice::get_conf_fact (  )  const [inline]

Returns the conformal factor.

Definition at line 210 of file excised_slice.h.

int Lorene::Excised_slice::get_field_set (  )  const [inline]

Returns the field set chosen for the data.

Definition at line 204 of file excised_slice.h.

References field_set.

const Sym_tensor& Lorene::Excised_slice::get_hatA (  )  const [inline]

Returns the rescaled tracefree extrinsic curvature $\hat{A}^{ij}$ (or its TT part if scheme_type=2).

Definition at line 219 of file excised_slice.h.

References hatA.

const Sym_tensor& Lorene::Excised_slice::get_hij (  )  const [inline]

Returns the deviation tensor $ h^{ij} $.

Definition at line 216 of file excised_slice.h.

References hij.

const Scalar& Lorene::Excised_slice::get_lapse (  )  const [inline]

Returns the lapse function N.

Definition at line 207 of file excised_slice.h.

References lapse.

const Map& Lorene::Excised_slice::get_mp (  )  const [inline]

Returns the mapping.

Definition at line 197 of file excised_slice.h.

References mp.

const Vector& Lorene::Excised_slice::get_shift (  )  const [inline]

Returns the shift vector $\beta^i$.

Definition at line 213 of file excised_slice.h.

References shift.

int Lorene::Excised_slice::get_type_hor (  )  const [inline]

Returns the type of horizon that performs the excision.

Definition at line 200 of file excised_slice.h.

References type_hor.

const Vector& Lorene::Excised_slice::get_Xx (  )  const [inline]

Return the longitudinal part of Einstein Equations if the scheme is appropriate.

Definition at line 222 of file excised_slice.h.

References field_set, and Xx.

Spheroid Lorene::Excised_slice::hor (  ) 

Spheroid at the excised boundary associated with the black hole MOTS on the slice.

Set by default at the position $ r=1 $.

Definition at line 246 of file excised_slice.C.

References Lorene::Metric_flat::con(), hatA, hij, mp, Lorene::pow(), Lorene::Tensor::std_spectral_base(), Lorene::Scalar::std_spectral_base(), Lorene::Tensor::up_down(), and Lorene::Scalar::val_grid_point().

double Lorene::Excised_slice::komar_angmom (  ) 
void Lorene::Excised_slice::operator= ( const Excised_slice ih  ) 

Assignment to another Excised_slice.

Definition at line 159 of file excised_slice.C.

References del_deriv(), field_set, hatA, hij, lapse, mp, shift, type_hor, and Xx.

void Lorene::Excised_slice::sauve ( FILE *  fich  )  const [virtual]

Save in a file.

Definition at line 181 of file excised_slice.C.

References hatA, hij, lapse, Lorene::Tensor_sym::sauve(), Lorene::Tensor::sauve(), Lorene::Scalar::sauve(), shift, and Xx.

void Lorene::Excised_slice::secmembre_kerr ( Sym_tensor source_hh  ) 
Scalar& Lorene::Excised_slice::set_conf_fact (  )  [inline]

Sets the conformal factor.

Definition at line 234 of file excised_slice.h.

References del_deriv().

void Lorene::Excised_slice::set_der_0x0 (  )  const [protected, virtual]

Sets to 0x0 all the pointers on derived quantities.

Definition at line 144 of file excised_slice.C.

References p_adm_mass, p_hor, p_komar_angmom, and p_virial_residue.

Sym_tensor& Lorene::Excised_slice::set_hatA (  )  [inline]

Sets the rescaled tracefree extrinsic curvature $\hat{A}^{ij}$ (or its TT part if scheme_type=2).

Definition at line 243 of file excised_slice.h.

References del_deriv(), and hatA.

Sym_tensor& Lorene::Excised_slice::set_hij (  )  [inline]

Sets the deviation tensor.

Definition at line 240 of file excised_slice.h.

References del_deriv(), and hij.

Scalar& Lorene::Excised_slice::set_lapse (  )  [inline]

Sets the lapse.

Definition at line 231 of file excised_slice.h.

References del_deriv(), and lapse.

Vector& Lorene::Excised_slice::set_shift (  )  [inline]

Sets the shift vector.

Definition at line 237 of file excised_slice.h.

References del_deriv(), and shift.

int Lorene::Excised_slice::set_type_hor (  )  [inline]

Sets the horizon type.

Definition at line 228 of file excised_slice.h.

References del_deriv(), and type_hor.

Vector& Lorene::Excised_slice::set_Xx (  )  [inline]

Sets the longitudinal part of Einstein Equations if the scheme is apropriate.

Definition at line 246 of file excised_slice.h.

References del_deriv(), and Xx.

double Lorene::Excised_slice::virial_residue (  ) 

Computation of the Virial residual, as difference at infinity between the ADM mass and the Komar integral associated to the mass.

Definition at line 333 of file excised_slice.C.

References Lorene::Scalar::asymptot(), Lorene::Map::get_mg(), Lorene::Mg3d::get_nzone(), lapse, and mp.


Member Data Documentation

Chose field set type.

Definition at line 55 of file excised_slice.h.

Rescaled tracefree extrinsic curvature tensor: rescaled same way as Cordero et al.

( 2009). In the New FCF scheme, this member will represent the TT part of the extrinsic curvature only.

Definition at line 74 of file excised_slice.h.

Deviation tensor( non-flat part of the conformal 3-metric on the slice; see Bonazzola et al.

(2004)).

Definition at line 69 of file excised_slice.h.

Lapse function.

Definition at line 58 of file excised_slice.h.

const Map& Lorene::Excised_slice::mp [protected]

Mapping associated with the slice.

Definition at line 46 of file excised_slice.h.

double* Lorene::Excised_slice::p_adm_mass [mutable, protected]

Computation of the ADM mass of the BH spacetime.

Definition at line 91 of file excised_slice.h.

Spheroid* Lorene::Excised_slice::p_hor [mutable, protected]

Computation of the spheroid associated with the black hole horizon.

Definition at line 89 of file excised_slice.h.

double* Lorene::Excised_slice::p_komar_angmom [mutable, protected]

Computation of the Komar angular momentum w.r.t.

assumed rotational symmetry

Definition at line 96 of file excised_slice.h.

double* Lorene::Excised_slice::p_virial_residue [mutable, protected]

Computation of the Virial residual, as difference at infinity between the ADM mass and the Komar integral associated to the mass.

Definition at line 101 of file excised_slice.h.

Shift vector.

Definition at line 64 of file excised_slice.h.

Chosen horizon type.

Definition at line 52 of file excised_slice.h.

Longitudinal part of the extrinsic curvature.

Set to zero if we are in the original FCF scheme.

Definition at line 79 of file excised_slice.h.


The documentation for this class was generated from the following files:

Generated on 7 Dec 2019 for LORENE by  doxygen 1.6.1