LORENE
Stars and black holes

Classes

class  Lorene::Bhole
 Black hole. More...
 
class  Lorene::Bhole_binaire
 Binary black holes system. More...
 
class  Lorene::Bin_bhns
 Class for computing a black hole - neutron star binary system with comparable mass () More...
 
class  Lorene::Bin_bhns_extr
 Class for computing a Black hole - Neutron star binary system with an extreme mass ratio. More...
 
class  Lorene::Binary
 Binary systems. More...
 
class  Lorene::Binary_xcts
 Binary systems in eXtended Conformal Thin Sandwich formulation. More...
 
class  Lorene::Black_hole
 Base class for black holes. More...
 
class  Lorene::Et_bin_bhns_extr
 Class for a neutron star in black hole - neutron star binary systems. More...
 
class  Lorene::Et_rot_bifluid
 Class for two-fluid rotating relativistic stars. More...
 
class  Lorene::Et_rot_diff
 Class for differentially rotating stars. More...
 
class  Lorene::Et_rot_mag
 Class for magnetized (isolator or perfect conductor), rigidly rotating stars. More...
 
class  Lorene::Etoile
 
Base class for stars *** DEPRECATED : use class Star instead ***. More...
 
class  Lorene::Etoile_bin
 Class for stars in binary system. More...
 
class  Lorene::Etoile_rot
 Class for isolated rotating stars *** DEPRECATED : use class Star_rot instead ***. More...
 
class  Lorene::Excision_hor
 Surface where boundary conditions for quantities in the bulk will be calculated It relies on geometrical properties of the associated Spheroid() (*** WARNING! under development***) More...
 
class  Lorene::Excision_surf
 Surface where boundary conditions for quantities in the bulk will be calculated It relies on geometrical properties of the associated Spheroid() (*** WARNING! under development***) More...
 
class  Lorene::Gravastar
 Class for perfect fluid rotating gravastar. More...
 
class  Lorene::Hole_bhns
 Class for black holes in black hole-neutron star binaries. More...
 
class  Lorene::Single_hor
 Binary black holes system. More...
 
class  Lorene::Spheroid
 Spheroidal 2-surfaces embedded in a time-slice of the 3+1 formalism. More...
 
class  Lorene::Star
 Base class for stars. More...
 
class  Lorene::Star_bin
 Class for stars in binary system. More...
 
class  Lorene::Star_bin_xcts
 Class for stars in binary system in eXtended Conformal Thin Sandwich formulation. More...
 
class  Lorene::Star_bhns
 Class for stars in black hole-neutron star binaries. More...
 
class  Lorene::Star_rot
 Class for isolated rotating stars. More...
 
class  Lorene::Star_rot_CFC
 Class for relativistic rotating stars in Conformal Flatness Condition and maximal slicing. More...
 
class  Lorene::Star_rot_Dirac
 Class for relativistic rotating stars in Dirac gauge and maximal slicing. More...
 
class  Lorene::Star_rot_Dirac_diff
 Class for relativistic differentially rotating stars in Dirac gauge and maximal slicing. More...
 

Functions

bool Lorene::ah_finder (const Metric &gamma, const Sym_tensor &k_dd_in, Valeur &h, Scalar &ex_fcn, double a_axis, double b_axis, double c_axis, bool verbose=true, bool print=false, double precis=1.e-8, double precis_exp=1.e-6, int it_max=200, int it_relax=200, double relax_fac=1.)
 Apparent horizon linked functions. More...
 

Detailed Description

Function Documentation

◆ ah_finder()

bool Lorene::ah_finder ( const Metric gamma,
const Sym_tensor k_dd_in,
Valeur h,
Scalar ex_fcn,
double  a_axis,
double  b_axis,
double  c_axis,
bool  verbose = true,
bool  print = false,
double  precis = 1.e-8,
double  precis_exp = 1.e-6,
int  it_max = 200,
int  it_relax = 200,
double  relax_fac = 1. 
)

Apparent horizon linked functions.

Apparent horizon finder.()

Find the apparent horizon (AH) on a spatial slice for given 3-metric $\gamma_{ij}$ and extrinsic curvature $K_{ij}$. Method: We solve the apparent-horizon equation $\Theta=0$ (where $\Theta$ is the expansion function of the outgoing null geodesic) as a generalized angular Poisson equation ( $\Delta_{\theta\phi}u+\lambda u=\sigma$) in the form:

\[ \Delta_{\theta\phi}h - 2h = \Psi^4 |DF| h^2\Theta + \Delta_{\theta\phi}h - 2h \]

where $h=h(\theta,\phi)$ (h) is the 2-surface of the AH, $\Psi$ is the conformal factor, and $|DF|:= (\gamma^{ij}D_iFD_jF)^{1/2} $ (where $F$ is the level set function $ F:=r-h(\theta,\phi) $ and $D_i$ is the covariant derivative w.r.t. $\gamma_{ij}$).
We solve the above equation iteratively with the right hand side of the equation as a source term.

Parameters
gamma: [input] the 3-metric $\gamma_{ij}$ w.r.t. which the AH is to be found.
k_dd_in: [input] the extrinsic curvature $K_{ij}$.
Valeurh : [output] the 2-surface of the apparent horizon
Scalarex_fcn : [output] the expansion function defined from the level set function $F:= r - h(\theta,\phi)$
a_axis: [input] the initial guess for $h$ is a triaxial ellipsoidal surface with a_axis the axis along the x-axis of
the Cartesian grid
b_axis: [input] axis along the y-axis (cf a_axis)
c_axis: [input] axis along the z-axis (cf a_axis)
boolprint : [input] screen printout during iterations? (default: false)
precis: [input] threshold in the relative difference between the 2-surface $h$ of two consecutive steps to stop the iterative procedure (default value: 1.e-8)
precis_exp: [input] maximum error of the expansion function evaluated on the 2-surface $h$ (should be zero by definition) after the iteration is stopped (default value: 1.e-6)
it_max: [input] maximum number of steps (default value: 200)
it_relax: [input] step at which relaxation is switched on (default value: 200)
relax_fac: [input] relaxation factor (default value: 1 for no relaxation)
Returns
bool ah_flag : a flag to indicate whether an apparent horizon is found

Definition at line 97 of file app_hor_finder.C.

References Lorene::abs(), Lorene::Scalar::allocate_all(), Lorene::Valeur::annule_hard(), Lorene::Valeur::c, Lorene::Valeur::c_cf, Lorene::Valeur::coef_i(), Lorene::Metric::con(), Lorene::Metric_flat::con(), Lorene::Metric::connect(), Lorene::contract(), Lorene::Map::cosp, Lorene::Map::cost, Lorene::Scalar::derive_con(), Lorene::Tensor::derive_cov(), Lorene::Scalar::derive_cov(), Lorene::Metric::determinant(), Lorene::Metric_flat::determinant(), Lorene::diffrel(), Lorene::Vector::divergence(), Lorene::Map::flat_met_spher(), Lorene::Mg3d::get_angu(), Lorene::Map::get_bvect_spher(), Lorene::Connection::get_delta(), Lorene::Map::get_mg(), Lorene::Metric::get_mp(), Lorene::Mg3d::get_np(), Lorene::Mg3d::get_nr(), Lorene::Mg3d::get_nt(), Lorene::Mg3d::get_nzone(), Lorene::max(), Lorene::min(), Lorene::Map::phi, Lorene::Mtbl_cf::poisson_angu(), Lorene::pow(), Lorene::Map::r, Lorene::Valeur::set(), Lorene::Scalar::set_etat_zero(), Lorene::Scalar::set_grid_point(), Lorene::Map::sinp, Lorene::Map::sint, Lorene::sqrt(), Lorene::Valeur::std_base_scal(), Lorene::Scalar::std_spectral_base(), Lorene::Map::tet, Lorene::Tensor::trace(), Lorene::Scalar::val_point(), Lorene::Valeur::ylm(), and Lorene::Valeur::ylm_i().