LORENE
Lorene::Bin_hor Class Reference

Public Member Functions

 Bin_hor (Map_af &mp1, Map_af &mp2)
 Standard constructor. More...
 
 Bin_hor (const Bin_hor &)
 Copy constructor. More...
 
 Bin_hor (Map_af &mp1, Map_af &mp2, FILE *fich)
 Constructor from a binary file. More...
 
virtual ~Bin_hor ()
 Destructor. More...
 
void sauve (FILE *fich) const
 Total or partial saves in a binary file. More...
 
void write_global (ostream &, double lim_nn, int bound_nn, int bound_psi, int bound_beta, double alpha) const
 Write global quantities in a formatted file. More...
 
void operator= (const Bin_hor &)
 Affectation operator. More...
 
Single_horset (int i)
 Read/write of a component of the system. More...
 
void set_omega (double ome)
 Sets the orbital velocity to ome. More...
 
const Single_horoperator() (int i) const
 Read only of a component of the system. More...
 
double get_omega () const
 Returns the angular velocity. More...
 
void init_bin_hor ()
 Initialisation of the system. More...
 
double viriel () const
 Computes the viriel error, that is the difference between the ADM and the Komar masses, calculated by the asymptotic behaviours of respectively $\Psi$ and N . More...
 
void extrinsic_curvature ()
 Calculation of the extrinsic curvature tensor. More...
 
void decouple ()
 Calculates decouple which is used to obtain tkij_auto and tkij_comp. More...
 
void set_statiques (double precis, double relax, int bound_nn, double lim_nn, int bound_psi)
 Initialize the systeme to Misner Lindquist solution, that is solving for N and $\Psi$ in the case $\Omega = 0$. More...
 
double coal (double ang_vel, double relax, int nb_om, int nb_it, int bound_nn, double lim_nn, int bound_psi, int bound_beta, double omega_eff, double alpha, ostream &fich_iteration, ostream &fich_correction, ostream &fich_viriel, ostream &fich_kss, int step, int search_mass, double mass_irr, const int sortie=0)
 Solves the equation for a particular angular velocity, the systeme being initialized to Misner-Lindquist solution. More...
 
void solve_lapse (double precis, double relax, int bound_nn, double lim_nn)
 Solves the equation for the lapse : The fields are the total values except those with subscript $_a$, which are the fields generated by each holes (a = 1, 2). More...
 
void solve_psi (double precis, double relax, int bound_psi)
 Solves the equation for the conformal factor : The fields are the total values excpet those with subscript $_a$, which are the fields generated by each holes (a = 1, 2). More...
 
void solve_shift (double precis, double relax, int bound_beta, double omega_eff)
 Solves the equation for the shift, using the Oohara-Nakarmure scheme : The fields are the total values excpet those with subscript $_a$, which are the fields generated by each holes (a = 1, 2). More...
 
void import_bh (const Bin_hor &bin)
 Function to initialize a Bin_hor from a solution computed with a smaller number of colocation points. More...
 
double adm_mass () const
 Calculates the ADM mass of the system. More...
 
double komar_mass () const
 Calculates the Komar mass of the system using : $M = \frac{1}{4 \pi} \oint_{\infty} D^i N {\mathrm d} S_i$. More...
 
double ang_mom_hor () const
 Calculates the angular momentum of the black hole using the formula at the horizon. More...
 
double ang_mom_adm () const
 Calculates the angular momentum of the black hole. More...
 
double proper_distance (const int nr=65) const
 Calculation of the proper distance between the two spheres of inversion, along the x axis. More...
 
Sym_tensor hh_Samaya_hole1 ()
 Calculation of the hole1 part of the Post-Newtonian correction to $h^{ij}$. More...
 
Sym_tensor hh_Samaya_hole2 ()
 Calculation of the hole2 part of the Post-Newtonian correction to $h^{ij}$. More...
 
void set_hh_Samaya ()
 Calculation of the Post-Newtonian correction to $h^{ij}$. More...
 

Private Attributes

Single_hor hole1
 Black hole one. More...
 
Single_hor hole2
 Black hole two. More...
 
Single_horholes [2]
 Array on the black holes. More...
 
double omega
 Angular velocity. More...
 

Detailed Description

Definition at line 1337 of file isol_hor.h.

Constructor & Destructor Documentation

◆ Bin_hor() [1/3]

Lorene::Bin_hor::Bin_hor ( Map_af mp1,
Map_af mp2 
)

Standard constructor.

Parameters
mp1affine mapping for the first black hole
mp2affine mapping for the second black hole
depth_innumber of stored time slices; this parameter is used to set the scheme_order member with scheme_order = depth_in - 1. scheme_order can be changed afterwards by the method set_scheme_order(int).

Definition at line 100 of file bin_hor.C.

References hole1, hole2, and holes.

◆ Bin_hor() [2/3]

Lorene::Bin_hor::Bin_hor ( const Bin_hor source)

Copy constructor.

Definition at line 110 of file bin_hor.C.

References hole1, hole2, and holes.

◆ Bin_hor() [3/3]

Lorene::Bin_hor::Bin_hor ( Map_af mp1,
Map_af mp2,
FILE *  fich 
)

Constructor from a binary file.

Parameters
mp1affine mapping of the first black hole
mp2affine mapping of the second black hole
fichfile containing the saved Bin_hor
partial_readindicates whether the full object must be read in file or whether the final construction is devoted to a constructor of a derived class
depth_innumber of stored time slices; this parameter is used to set the scheme_order member with scheme_order = depth_in - 1. scheme_order can be changed afterwards by the method set_scheme_order(int).

Definition at line 120 of file bin_hor.C.

References Lorene::fread_be(), hole1, hole2, holes, and omega.

◆ ~Bin_hor()

Lorene::Bin_hor::~Bin_hor ( )
virtual

Destructor.

Definition at line 135 of file bin_hor.C.

Member Function Documentation

◆ adm_mass()

◆ ang_mom_adm()

◆ ang_mom_hor()

◆ coal()

double Lorene::Bin_hor::coal ( double  ang_vel,
double  relax,
int  nb_om,
int  nb_it,
int  bound_nn,
double  lim_nn,
int  bound_psi,
int  bound_beta,
double  omega_eff,
double  alpha,
ostream &  fich_iteration,
ostream &  fich_correction,
ostream &  fich_viriel,
ostream &  fich_kss,
int  step,
int  search_mass,
double  mass_irr,
const int  sortie = 0 
)

Solves the equation for a particular angular velocity, the systeme being initialized to Misner-Lindquist solution.

Parameters
angu[input] : angular velocity used for the boundary condition on $\vec{\beta}$.
relax[input] : relaxation parameter.
nb_om[input] : number of intermediates velocities to go from 0 to omega , typically 10.
nb_it[input] : number of iteration when omega is fixed
bound_nn[input] : type of the boundary condition for the lapse.
lim_nn[input] : value (double) of the coefficient for the boundary condition. Only used for boundary_nn_Dir(double), boundary_nn_Neu_eff(double) and boundary_nn_Dir_eff(double).
bound_psi[input] : type of the boundary condition for psi.
bound_nn[input] : type of the boundary condition for the shift.
stepcurrent step of the iteration
sortie[input] : flag for the output on files (0 no output files).
Returns
: the virial error.

Definition at line 136 of file binhor_coal.C.

References ang_mom_hor(), Lorene::Single_hor::area_hor(), Lorene::Single_hor::beta_auto, Lorene::contract(), Lorene::diffrelmax(), extrinsic_curvature(), Lorene::Single_hor::ff, Lorene::Map::flat_met_spher(), Lorene::Single_hor::get_gam(), Lorene::Single_hor::get_k_dd(), Lorene::Map::get_mg(), Lorene::Mg3d::get_np(), Lorene::Mg3d::get_nt(), Lorene::Mg3d::get_nzone(), get_omega(), hole1, hole2, Lorene::Map_af::homothetie_interne(), Lorene::log10(), Lorene::Single_hor::mp, Lorene::Single_hor::omega_hor(), Lorene::pow(), Lorene::Metric::radial_vect(), Lorene::Single_hor::radius, Lorene::Single_hor::regul, set_omega(), solve_lapse(), solve_psi(), solve_shift(), Lorene::sqrt(), Lorene::Single_hor::tgam, Lorene::Scalar::val_grid_point(), and viriel().

◆ decouple()

◆ extrinsic_curvature()

◆ get_omega()

double Lorene::Bin_hor::get_omega ( ) const
inline

Returns the angular velocity.

Definition at line 1422 of file isol_hor.h.

References omega.

◆ hh_Samaya_hole1()

◆ hh_Samaya_hole2()

◆ import_bh()

void Lorene::Bin_hor::import_bh ( const Bin_hor bin)

Function to initialize a Bin_hor from a solution computed with a smaller number of colocation points.

◆ init_bin_hor()

void Lorene::Bin_hor::init_bin_hor ( )

Initialisation of the system.

Each hole is set close to a Schwarzschild one and the parts of the fields generated by the companion are calculated.

The angular velocity is set to zero.

Definition at line 164 of file bin_hor.C.

References decouple(), extrinsic_curvature(), hole1, hole2, Lorene::Single_hor::init_bhole(), Lorene::Single_hor::n_comp_import(), Lorene::Single_hor::psi_comp_import(), and set_omega().

◆ komar_mass()

◆ operator()()

const Single_hor& Lorene::Bin_hor::operator() ( int  i) const
inline

Read only of a component of the system.

i must be equal to 1 or 2.

Definition at line 1417 of file isol_hor.h.

◆ operator=()

void Lorene::Bin_hor::operator= ( const Bin_hor source)

Affectation operator.

Definition at line 142 of file bin_hor.C.

References hole1, hole2, and omega.

◆ proper_distance()

double Lorene::Bin_hor::proper_distance ( const int  nr = 65) const

Calculation of the proper distance between the two spheres of inversion, along the x axis.

Parameters
nr[input] : number of points used for the calculation.

◆ sauve()

void Lorene::Bin_hor::sauve ( FILE *  fich) const

Total or partial saves in a binary file.

Parameters
fichbinary file
partial_saveindicates whether the whole object must be saved.

Definition at line 154 of file bin_hor.C.

References Lorene::fwrite_be(), hole1, hole2, omega, and Lorene::Single_hor::sauve().

◆ set()

Single_hor& Lorene::Bin_hor::set ( int  i)
inline

Read/write of a component of the system.

i must be equal to 1 or 2.

Definition at line 1403 of file isol_hor.h.

References holes.

◆ set_hh_Samaya()

◆ set_omega()

void Lorene::Bin_hor::set_omega ( double  ome)
inline

Sets the orbital velocity to ome.

Definition at line 1409 of file isol_hor.h.

References hole1, hole2, omega, and Lorene::Single_hor::set_omega().

◆ set_statiques()

void Lorene::Bin_hor::set_statiques ( double  precis,
double  relax,
int  bound_nn,
double  lim_nn,
int  bound_psi 
)

Initialize the systeme to Misner Lindquist solution, that is solving for N and $\Psi$ in the case $\Omega = 0$.

Parameters
precis[input] : precision for the convergence (on N ).
relax[input] : relaxation parameter.
bound_nn[input] : type of the boundary condition for the lapse.
lim_nn[input] : value (double) of the coefficient for the boundary condition. Only used for boundary_nn_Dir(double), boundary_nn_Neu_eff(double) and boundary_nn_Dir_eff(double).
bound_psi[input] : type of the boundary condition for psi.

Definition at line 99 of file binhor_coal.C.

References Lorene::diffrelmax(), extrinsic_curvature(), Lorene::Map::get_mg(), Lorene::Mg3d::get_nzone(), hole1, hole2, init_bin_hor(), Lorene::Single_hor::init_met_trK(), Lorene::Single_hor::mp, Lorene::Single_hor::n_auto, set_omega(), solve_lapse(), and solve_psi().

◆ solve_lapse()

void Lorene::Bin_hor::solve_lapse ( double  precis,
double  relax,
int  bound_nn,
double  lim_nn 
)

Solves the equation for the lapse : The fields are the total values except those with subscript $_a$, which are the fields generated by each holes (a = 1, 2).

Parameters
precis[input] : precision, for the boudary conditions are obtained by iteration.
relax[input] : relaxation parameter.
bound_nn[input] : type of the boundary condition at the horizon.
lim_nn[input] : value (double) of the coefficient for the boundary condition. Only used for boundary_nn_Dir(double), boundary_nn_Neu_eff(double) and boundary_nn_Dir_eff(double).

Definition at line 127 of file binhor_equations.C.

References Lorene::Single_hor::aa, Lorene::Single_hor::aa_auto, Lorene::Single_hor::beta_auto, Lorene::contract(), Lorene::Single_hor::decouple, Lorene::Scalar::derive_con(), Lorene::Tensor::derive_cov(), Lorene::Scalar::derive_cov(), Lorene::diffrel(), Lorene::Single_hor::dn, Lorene::Single_hor::dpsi, Lorene::Single_hor::ff, Lorene::Map::get_mg(), Lorene::Mg3d::get_nzone(), Lorene::Single_hor::get_psi4(), Lorene::Single_hor::hh, hole1, hole2, Lorene::Scalar::inc_dzpuis(), Lorene::Single_hor::mp, Lorene::Single_hor::n_auto, Lorene::Single_hor::n_comp_import(), Lorene::Single_hor::nn, Lorene::norme(), Lorene::Single_hor::psi, Lorene::Single_hor::psi_auto, Lorene::Scalar::raccord(), Lorene::Single_hor::tgam, Lorene::Single_hor::trK, Lorene::Single_hor::trK_point, and Lorene::Tensor::up_down().

◆ solve_psi()

void Lorene::Bin_hor::solve_psi ( double  precis,
double  relax,
int  bound_psi 
)

Solves the equation for the conformal factor : The fields are the total values excpet those with subscript $_a$, which are the fields generated by each holes (a = 1, 2).

Parameters
precis[input] : precision, for the boudary conditions are being obtained by iteration.
relax[input] : relaxation parameter.
bound_psi[input] : type of the boundary condition at the horizon.

Definition at line 304 of file binhor_equations.C.

References Lorene::Single_hor::aa, Lorene::Single_hor::aa_auto, Lorene::Single_hor::boundary_psi_app_hor(), Lorene::contract(), Lorene::Single_hor::decouple, Lorene::Tensor::derive_cov(), Lorene::Scalar::derive_cov(), Lorene::diffrel(), Lorene::Single_hor::ff, Lorene::Map::get_mg(), Lorene::Mg3d::get_nzone(), Lorene::Single_hor::get_psi4(), Lorene::Single_hor::hh, hole1, hole2, Lorene::Scalar::inc_dzpuis(), Lorene::Single_hor::mp, Lorene::norme(), Lorene::Single_hor::psi, Lorene::Single_hor::psi_auto, Lorene::Single_hor::psi_comp_import(), Lorene::Scalar::raccord(), Lorene::Single_hor::set_der_0x0(), Lorene::Scalar::std_spectral_base(), Lorene::Single_hor::tgam, Lorene::Single_hor::trK, and Lorene::Tensor::up_down().

◆ solve_shift()

void Lorene::Bin_hor::solve_shift ( double  precis,
double  relax,
int  bound_beta,
double  omega_eff 
)

Solves the equation for the shift, using the Oohara-Nakarmure scheme : The fields are the total values excpet those with subscript $_a$, which are the fields generated by each holes (a = 1, 2).

Parameters
precis[input] : precision for the solver, the boundary conditions and the inversion of the operator being obtained by iteration.
relax[input] : relaxation parameter.
bound_beta[input] : type of the boundary condition at the horizon.

Definition at line 456 of file binhor_equations.C.

References Lorene::Single_hor::aa, Lorene::Single_hor::aa_auto, Lorene::Scalar::annule_hard(), Lorene::Single_hor::beta_auto, Lorene::Single_hor::beta_comp_import(), Lorene::Vector::change_triad(), Lorene::Metric::connect(), Lorene::contract(), Lorene::Tensor::dec_dzpuis(), Lorene::Single_hor::decouple, Lorene::Tensor::derive_con(), Lorene::Scalar::derive_con(), Lorene::Tensor::derive_cov(), Lorene::Scalar::derive_cov(), Lorene::Tensor::derive_lie(), Lorene::diffrel(), Lorene::Sym_tensor::divergence(), Lorene::Vector::divergence(), Lorene::Single_hor::dn, Lorene::Single_hor::dpsi, Lorene::Single_hor::ff, Lorene::Single_hor::gamt_point, Lorene::Map::get_bvect_cart(), Lorene::Map::get_bvect_spher(), Lorene::Connection::get_delta(), Lorene::Map::get_mg(), Lorene::Mg3d::get_np(), Lorene::Mg3d::get_nt(), Lorene::Mg3d::get_nzone(), Lorene::Map::get_rot_phi(), Lorene::Single_hor::hh, hole1, hole2, Lorene::Tensor::inc_dzpuis(), Lorene::Scalar::inc_dzpuis(), Lorene::Single_hor::mp, Lorene::Single_hor::n_auto, Lorene::Single_hor::n_comp, Lorene::Single_hor::nn, Lorene::norme(), omega, Lorene::Single_hor::psi, Lorene::Scalar::raccord(), Lorene::Single_hor::regul, Lorene::Single_hor::regularisation(), Lorene::Coord::set(), Lorene::Tbl::set(), Lorene::Vector::set(), Lorene::Tensor::set(), Lorene::Mg3d::std_base_vect_cart(), Lorene::Tensor::std_spectral_base(), Lorene::Single_hor::tgam, Lorene::Single_hor::trK, Lorene::Tensor::up_down(), Lorene::Map::xa, and Lorene::Map::ya.

◆ viriel()

double Lorene::Bin_hor::viriel ( ) const

Computes the viriel error, that is the difference between the ADM and the Komar masses, calculated by the asymptotic behaviours of respectively $\Psi$ and N .

Definition at line 74 of file binhor_viriel.C.

References Lorene::Scalar::asymptot(), Lorene::Map::get_mg(), Lorene::Mg3d::get_nzone(), hole1, hole2, Lorene::Single_hor::mp, Lorene::Single_hor::n_auto, and Lorene::Single_hor::psi_auto.

◆ write_global()

Member Data Documentation

◆ hole1

Single_hor Lorene::Bin_hor::hole1
private

Black hole one.

Definition at line 1342 of file isol_hor.h.

◆ hole2

Single_hor Lorene::Bin_hor::hole2
private

Black hole two.

Definition at line 1343 of file isol_hor.h.

◆ holes

Single_hor* Lorene::Bin_hor::holes[2]
private

Array on the black holes.

Definition at line 1346 of file isol_hor.h.

◆ omega

double Lorene::Bin_hor::omega
private

Angular velocity.

Definition at line 1348 of file isol_hor.h.


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