76 #include "bin_bhns_extr.h" 78 #include "utilitaires.h" 89 bool relat,
bool kerrs,
bool multi)
90 : ref_triad(0.,
"Absolute frame Cartesian basis"),
91 star(mp, nzet, relat, eos, irrot, ref_triad, kerrs, multi)
106 : ref_triad(0.,
"Absolute frame Cartesian basis"),
110 mass_bh(bibi.mass_bh)
121 : ref_triad(0.,
"Absolute frame Cartesian basis"),
122 star(mp, eos, ref_triad, fich)
218 ost <<
"Binary BH-NS system" << endl ;
219 ost <<
"===================" << endl ;
222 ost <<
"Kerr-Schild background metric" << endl ;
223 ost <<
"-----------------------------" << endl ;
226 ost <<
"Conformally flat background metric" << endl ;
227 ost <<
"----------------------------------" << endl ;
230 ost <<
"Multipole falloff boundary condition" << endl ;
231 ost <<
"------------------------------------" << endl ;
234 ost <<
"1/r falloff boundary condition" << endl ;
235 ost <<
"------------------------------" << endl ;
238 <<
"Orbital angular velocity : " 239 <<
omega * f_unit <<
" rad/s" << endl ;
240 ost <<
"Coordinate separation between BH and NS : " 241 <<
separ / km <<
" km" << endl ;
243 ost << endl <<
"Neutron star : " ;
244 ost << endl <<
"============ " << endl ;
247 ost <<
"Relativistic star" << endl ;
248 ost <<
"-----------------" << endl ;
251 ost <<
"WARNING : BH-NS binary should be relativistic !!!" << endl ;
252 ost <<
"-------------------------------------------------" << endl ;
254 ost <<
"Number of domains occupied by the star : " <<
star.
get_nzet()
256 ost <<
"Equation of state : " << endl ;
260 <<
"Enthalpy at the coordinate origin : " 262 ost <<
"Proper baryon density at the coordinate origin : " 263 <<
star.
get_nbar()()(0,0,0,0) <<
" x 0.1 fm^-3" << endl ;
264 ost <<
"Proper energy density at the coordinate origin : " 265 <<
star.
get_ener()()(0,0,0,0) <<
" rho_nuc c^2" << endl ;
266 ost <<
"Pressure at the coordinate origin : " 269 <<
"Lapse N at the coordinate origin : " 271 ost <<
"Conformal factor A^2 at the coordinate origin : " 275 <<
"Equatorial radius (to BH) a_to : " 277 ost <<
"Equatorial radius (opp. to BH) a_opp : " 281 <<
"Baryon mass in isolation : " <<
star.
mass_b() / msol
283 <<
"Gravitational mass in isolation : " <<
star.
mass_g() / msol
285 <<
"Baryon mass in a binary system : " <<
mass_b_extr() / msol
289 ost <<
"Star in a binary system" << endl ;
290 ost <<
"-----------------------" << endl ;
293 ost <<
"irrotational configuration" << endl ;
296 ost <<
"corotating configuration" << endl ;
299 ost <<
"Absolute abscidia of the stellar center: " 301 ost <<
"Orientation with respect to the absolute frame : " 308 double d_tilde = d_ns / r0 ;
310 ost << endl <<
"Comparison with those by Baumgarte et al. :" << endl ;
311 ost <<
" Radius r0 : " << r0 / km <<
" km " << endl ;
312 ost <<
" Separation d : " << d_ns / km <<
" km " << endl ;
313 ost <<
" Normalized sep. (d/r0) : " << d_tilde << endl ;
315 ost << endl <<
"Black hole : " ;
316 ost << endl <<
"========== " << endl ;
317 ost <<
"Gravitational mass of BH : " 318 <<
mass_bh / msol <<
" M_sol" << endl ;
333 if (p_eos_poly != 0x0) {
335 double kappa = p_eos_poly->
get_kap() ;
336 double gamma = p_eos_poly->
get_gam() ;
337 double kap_ns2 =
pow( kappa, 0.5 /(gamma-1) ) ;
340 double r_poly = kap_ns2 /
sqrt(ggrav) ;
343 double t_poly = r_poly ;
346 double m_poly = r_poly / ggrav ;
357 ost << endl <<
"Quantities in polytropic units : " ;
358 ost << endl <<
"==============================" << endl ;
359 ost <<
" ( r_poly = " << r_poly / km <<
" km )" << endl ;
360 ost <<
" d_e_max : " <<
separ / r_poly << endl ;
363 ost <<
" d_bh/M_bh : " << d_ns/r_poly / (
mass_bh/m_poly) << endl ;
364 ost <<
" Omega : " <<
omega * t_poly << endl ;
365 ost <<
" Omega M_bh : " <<
omega * t_poly *
mass_bh / m_poly
367 ost <<
" M_bar(NS) : " <<
mass_b_extr() / m_poly << endl ;
368 ost <<
" M_bar(NS_0) : " <<
star.
mass_b() / m_poly << endl ;
371 ost <<
" M_grv(BH) : " <<
mass_bh / m_poly << endl ;
bool with_multipole() const
Returns true for the multipole falloff boundary condition, false for the one.
double * p_xa_barycenter_extr
Absolute coordinate X of the barycenter of the baryon density in the Kerr-Schild background metric or...
const Map & get_mp() const
Returns the mapping.
const Tenseur & get_a_car() const
Returns the total conformal factor .
void operator=(const Bin_bhns_extr &)
Assignment to another Bin_bhns_extr.
Cmp sqrt(const Cmp &)
Square root.
double xa_barycenter_extr() const
Absolute coordinate X of the barycenter of the baryon density in the Kerr-Schild background metric or...
Standard units of space, time and mass.
Equation of state base class.
double ray_eq() const
Coordinate radius at , [r_unit].
Base class for coordinate mappings.
double get_ori_x() const
Returns the x coordinate of the origin.
int get_nzet() const
Returns the number of domains occupied by the star.
ostream & operator>>(ostream &) const
Operator >> (function called by the operator <<)
bool is_irrotational() const
Returns true for an irrotational motion, false for a corotating one.
virtual double mass_g() const
Gravitational mass.
double get_rot_phi() const
Returns the angle between the x –axis and X –axis.
Et_bin_bhns_extr star
Neutron star.
virtual void sauve(FILE *) const
Save in a file.
double get_gam() const
Returns the adiabatic index (cf. Eq. (3))
const Tenseur & get_press() const
Returns the fluid pressure.
bool in_kerrschild() const
Returns true for the Kerr-Schild background metric, false for the conformally flat one...
void del_deriv() const
Deletes all the derived quantities.
void display_poly(ostream &) const
Display in polytropic units.
double ya_barycenter_extr() const
in the Kerr-Schild background metric
void sauve(FILE *) const
Save in a file.
double get_kap() const
Returns the pressure coefficient (cf.
const Base_vect_cart ref_triad
Cartesian triad of the absolute reference frame.
Polytropic equation of state (relativistic case).
const Tenseur & get_nnn() const
Returns the total lapse function N.
int fwrite_be(const int *aa, int size, int nb, FILE *fich)
Writes integer(s) into a binary file according to the big endian convention.
Cmp pow(const Cmp &, int)
Power .
int fread_be(int *aa, int size, int nb, FILE *fich)
Reads integer(s) from a binary file according to the big endian convention.
double ray_eq_pi() const
Coordinate radius at , [r_unit].
const Tenseur & get_ener() const
Returns the proper total energy density.
Bin_bhns_extr(Map &mp, int nzet, const Eos &eos, bool irrot, bool relat, bool kerrs, bool multi)
Standard constructor.
void set_der_0x0() const
Sets to 0x0 all the pointers on derived quantities.
~Bin_bhns_extr()
Destructor.
double separ
Absolute orbital separation between two centers of BH and NS.
double omega
Angular velocity with respect to an asymptotically inertial observer.
double * p_mass_b_extr
Baryon mass of the neutron star in the Kerr-Schild background metric or in the conformally flat one...
virtual double mass_b() const
Baryon mass.
const Tenseur & get_nbar() const
Returns the proper baryon density.
bool is_relativistic() const
Returns true for a relativistic star, false for a Newtonian one.
double mass_bh
Gravitational mass of BH.
const Eos & get_eos() const
Returns the equation of state.
const Tenseur & get_ent() const
Returns the enthalpy field.
double * p_ya_barycenter_extr
Absolute coordinate Y of the barycenter of the baryon density in the Kerr-Schild background metric...
double mass_b_extr() const
Baryon mass of the neutron star in the Kerr-Schild background metric or in the conformally flat...
Class for computing a Black hole - Neutron star binary system with an extreme mass ratio...