LORENE
Lorene::Binaire Class Reference

Binary systems. More...

#include <binaire.h>

Public Member Functions

 Binaire (Map &mp1, int nzet1, const Eos &eos1, int irrot1, Map &mp2, int nzet2, const Eos &eos2, int irrot2, int relat)
 Standard constructor. More...
 
 Binaire (const Binaire &)
 
 Binaire (Map &mp1, const Eos &eos1, Map &mp2, const Eos &eos2, FILE *fich)
 Copy constructor. More...
 
void operator= (const Binaire &)
 Assignment to another { Binaire}. More...
 
Etoile_binset (int i)
 Read/write of the star no. i. More...
 
double & set_omega ()
 Sets the orbital angular velocity [{ f_unit}]. More...
 
double & set_x_axe ()
 Sets the absolute coordinate X of the rotation axis [{ r_unit}]. More...
 
const Etoile_binoperator() (int i) const
 Returns a reference to the star no. i. More...
 
double get_omega () const
 Returns the orbital angular velocity [{ f_unit}]. More...
 
double get_x_axe () const
 Returns the absolute coordinate X of the rotation axis [{ r_unit}]. More...
 
double separation () const
 Returns the coordinate separation of the two stellar centers [{ r_unit}]. More...
 
void sauve (FILE *) const
 
void display_poly (ostream &) const
 Display in polytropic units. More...
 
void write_global (ostream &) const
 Write global quantities in a formatted file. More...
 
double mass_adm () const
 Total ADM mass. More...
 
double mass_kom () const
 Total Komar mass. More...
 
const Tblangu_mom () const
 Total angular momentum. More...
 
double total_ener () const
 Total energy (excluding the rest mass energy). More...
 
double virial () const
 Estimates the relative error on the virial theorem (for a relativistic one, it returns $|1 - M_{ Komar} / M_{ ADM}|$) More...
 
double virial_gb () const
 Estimates the relative error on the virial theorem calculated by E.Gourgoulhon and S.Bonazzola (Class. More...
 
double virial_fus () const
 Estimates the relative error on the virial theorem calculated by J.L.Friedman, K.Uryu, and M.Shibata (PRD accepted, gr-qc/0108070) More...
 
double ham_constr () const
 Estimates the relative error on the Hamiltonian constraint equation by comparing $ A$ with {equation} -4 A^2 E - {A^2 4} K_{ij} K^{ij} - {1 2} {}_i A {}^i A {equation}. More...
 
const Tblmom_constr () const
 Estimates the relative error on the momentum constraint equation by comparing ${}_j K^{ij}$ with {equation} 8 J^i - 5 K^{ij} {}_j A {equation}. More...
 
void orbit (double fact_omeg_min, double fact_omeg_max, double &xgg1, double &xgg2)
 Computes the orbital angular velocity { omega} and the position of the rotation axis { x_axe}. More...
 
void orbit_eqmass (double fact_omeg_min, double fact_omeg_max, double mass1, double mass2, double &xgg1, double &xgg2)
 Computes the orbital angular velocity { omega} and the position of the rotation axis { x_axe}. More...
 
void analytical_omega ()
 Sets the orbital angular velocity to some 2-PN analytical value (Keplerian value in the Newtonian case) More...
 
void analytical_shift ()
 Sets some analytical template for the shift vector (via the members { w_shift} and { khi_shift} of the two { Etoile_bin}. More...
 

Private Member Functions

void del_deriv () const
 Destructor. More...
 
void set_der_0x0 () const
 Sets to { 0x0} all the pointers on derived quantities. More...
 
ostream & operator>> (ostream &) const
 Operator >> (function called by the operator <<). More...
 

Private Attributes

const Base_vect_cart ref_triad
 Cartesian triad of the absolute reference frame. More...
 
Etoile_bin star1
 First star of the system. More...
 
Etoile_bin star2
 Second star of the system. More...
 
Etoile_binet [2]
 Array of the two stars (to perform loops on the stars): { et[0]} contains the address of { star1} and { et[1]} that of { star2}. More...
 
double omega
 Angular velocity with respect to an asymptotically inertial observer. More...
 
double x_axe
 Absolute X coordinate of the rotation axis. More...
 
double * p_mass_adm
 Total ADM mass of the system. More...
 
double * p_mass_kom
 Total Komar mass of the system. More...
 
Tblp_angu_mom
 Total angular momentum of the system. More...
 
double * p_total_ener
 Total energy of the system. More...
 
double * p_virial
 Virial theorem error. More...
 
double * p_virial_gb
 Virial theorem error by E.Gourgoulhon and S.Bonazzola. More...
 
double * p_virial_fus
 Virial theorem error by J.L.Friedman, K.Uryu, and M.Shibata. More...
 
double * p_ham_constr
 Relative error on the Hamiltonian constraint. More...
 
Tblp_mom_constr
 Relative error on the momentum constraint. More...
 

Friends

ostream & operator<< (ostream &, const Binaire &)
 Save in a file. More...
 

Detailed Description

Binary systems.

Version
#$Id: binaire.h,v 1.7 2017/02/24 15:34:59 j_novak Exp $#

Definition at line 107 of file binaire.h.

Constructor & Destructor Documentation

◆ Binaire() [1/2]

Lorene::Binaire::Binaire ( Map mp1,
int  nzet1,
const Eos eos1,
int  irrot1,
Map mp2,
int  nzet2,
const Eos eos2,
int  irrot2,
int  relat 
)

Standard constructor.

Parameters
mp1Mapping on which { star1} will be defined
nzet1Number of domains occupied by { star1}
eos1Equation of state of { star1}
irrot1should be { true} if { star1} is irrotational, { false} if { star1} is corotating
mp2Mapping on which { star2} will be defined
nzet2Number of domains occupied by { star2}
eos2Equation of state of { star2}
irrot2should be { true} if { star2} is irrotational, { false} if { star2} is corotating
relatshould be { true} for a relativistic configuration, { false} for a Newtonian one

Definition at line 112 of file binaire.C.

References et, omega, set_der_0x0(), star1, star2, and x_axe.

◆ Binaire() [2/2]

Lorene::Binaire::Binaire ( Map mp1,
const Eos eos1,
Map mp2,
const Eos eos2,
FILE *  fich 
)

Copy constructor.

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

Parameters
mp1Mapping on which { star1} will be defined
eos1Equation of state of { star1}
mp2Mapping on which { star2} will be defined
eos2Equation of state of { star2}
fichinput file (must have been created by the function { sauve})

Definition at line 148 of file binaire.C.

References et, Lorene::fread_be(), omega, set_der_0x0(), star1, star2, and x_axe.

Member Function Documentation

◆ analytical_omega()

void Lorene::Binaire::analytical_omega ( )

Sets the orbital angular velocity to some 2-PN analytical value (Keplerian value in the Newtonian case)

Definition at line 66 of file binaire_omega_ana.C.

References del_deriv(), Lorene::Etoile_bin::is_irrotational(), Lorene::Etoile::is_relativistic(), Lorene::Etoile_bin::mass_g(), omega, Lorene::pow(), Lorene::Etoile::ray_eq(), separation(), Lorene::sqrt(), star1, and star2.

◆ analytical_shift()

◆ angu_mom()

◆ del_deriv()

void Lorene::Binaire::del_deriv ( ) const
private

Destructor.

Deletes all the derived quantities

Definition at line 180 of file binaire.C.

References p_angu_mom, p_ham_constr, p_mass_adm, p_mass_kom, p_mom_constr, p_total_ener, p_virial, p_virial_fus, p_virial_gb, and set_der_0x0().

◆ display_poly()

◆ get_omega()

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

Returns the orbital angular velocity [{ f_unit}].

Definition at line 247 of file binaire.h.

References omega.

◆ get_x_axe()

double Lorene::Binaire::get_x_axe ( ) const
inline

Returns the absolute coordinate X of the rotation axis [{ r_unit}].

Definition at line 250 of file binaire.h.

References x_axe.

◆ ham_constr()

◆ mass_adm()

◆ mass_kom()

◆ mom_constr()

◆ operator()()

const Etoile_bin& Lorene::Binaire::operator() ( int  i) const
inline

Returns a reference to the star no. i.

Definition at line 242 of file binaire.h.

References et.

◆ operator=()

void Lorene::Binaire::operator= ( const Binaire bibi)

Assignment to another { Binaire}.

Definition at line 220 of file binaire.C.

References del_deriv(), omega, ref_triad, star1, star2, and x_axe.

◆ operator>>()

ostream & Lorene::Binaire::operator>> ( ostream &  ost) const
private

Operator >> (function called by the operator <<).

Definition at line 260 of file binaire.C.

References omega, separation(), star1, star2, and x_axe.

◆ orbit()

void Lorene::Binaire::orbit ( double  fact_omeg_min,
double  fact_omeg_max,
double &  xgg1,
double &  xgg2 
)

Computes the orbital angular velocity { omega} and the position of the rotation axis { x_axe}.

Parameters
fact_omeg_min[input] : determines the lower bound of the interval { [omega_min, omega_max]} in which { omega} is searched by { omega_min = fact_omeg_min * omega}, where { omega} is the previous value of the angular velocity (typical value : { fact_omeg_min = 0.5})
fact_omeg_max[input] : determines the higher bound of the interval { [omega_min, omega_max]} in which { omega} is searched by { omega_max = fact_omeg_max * omega}, where { omega} is the previous value of the angular velocity. (typical value : { fact_omeg_max = 1.5})
xgg1[output] : x coordinate (relative to star 1 mapping) of the `‘center of mass’' of star 1
xgg2[output] : x coordinate (relative to star 2 mapping) of the `‘center of mass’' of star 2

Definition at line 125 of file binaire_orbite.C.

References Lorene::Param::add_double(), Lorene::Param::add_int(), Lorene::Tenseur::change_triad(), Lorene::Cmp::dsdx(), et, Lorene::flat_scalar_prod(), Lorene::Etoile::get_a_car(), Lorene::Map::get_bvect_cart(), Lorene::Etoile::get_d_logn_auto_div(), Lorene::Etoile::get_ent(), Lorene::Etoile::get_eos(), Lorene::Etoile_bin::get_loggam(), Lorene::Etoile::get_logn_auto_regu(), Lorene::Etoile_bin::get_logn_comp(), Lorene::Etoile::get_mp(), Lorene::Etoile::get_nnn(), Lorene::Map::get_rot_phi(), Lorene::Etoile::get_shift(), Lorene::Tbl::get_taille(), Lorene::Tenseur::get_triad(), Lorene::Etoile::is_relativistic(), omega, ref_triad, save_profile(), x_axe, Lorene::Etoile_bin::xa_barycenter(), Lorene::zero_list(), Lorene::zerosec(), and Lorene::zerosec_b().

◆ orbit_eqmass()

void Lorene::Binaire::orbit_eqmass ( double  fact_omeg_min,
double  fact_omeg_max,
double  mass1,
double  mass2,
double &  xgg1,
double &  xgg2 
)

Computes the orbital angular velocity { omega} and the position of the rotation axis { x_axe}.

Parameters
fact_omeg_min[input] : determines the lower bound of the interval { [omega_min, omega_max]} in which { omega} is searched by { omega_min = fact_omeg_min * omega}, where { omega} is the previous value of the angular velocity (typical value : { fact_omeg_min = 0.5})
fact_omeg_max[input] : determines the higher bound of the interval { [omega_min, omega_max]} in which { omega} is searched by { omega_max = fact_omeg_max * omega}, where { omega} is the previous value of the angular velocity. (typical value : { fact_omeg_max = 1.5})
mass1[input] : baryon rest mass of NS1
mass2[input] : baryon rest mass of NS2
xgg1[output] : x coordinate (relative to star 1 mapping) of the `‘center of mass’' of star 1
xgg2[output] : x coordinate (relative to star 2 mapping) of the `‘center of mass’' of star 2

Definition at line 467 of file binaire_orbite.C.

References Lorene::Param::add_double(), Lorene::Param::add_int(), Lorene::Tenseur::change_triad(), Lorene::Cmp::dsdx(), et, Lorene::flat_scalar_prod(), Lorene::Etoile::get_a_car(), Lorene::Map::get_bvect_cart(), Lorene::Etoile::get_d_logn_auto_div(), Lorene::Etoile::get_eos(), Lorene::Etoile_bin::get_loggam(), Lorene::Etoile::get_logn_auto_regu(), Lorene::Etoile_bin::get_logn_comp(), Lorene::Etoile::get_mp(), Lorene::Etoile::get_nnn(), Lorene::Map::get_rot_phi(), Lorene::Etoile::get_shift(), Lorene::Tbl::get_taille(), Lorene::Tenseur::get_triad(), Lorene::Etoile::is_relativistic(), omega, ref_triad, x_axe, Lorene::Etoile_bin::xa_barycenter(), Lorene::zero_list(), Lorene::zerosec(), and Lorene::zerosec_b().

◆ separation()

double Lorene::Binaire::separation ( ) const

Returns the coordinate separation of the two stellar centers [{ r_unit}].

Definition at line 460 of file binaire.C.

References Lorene::Etoile::get_mp(), Lorene::Map::get_ori_x(), Lorene::Map::get_ori_y(), Lorene::Map::get_ori_z(), Lorene::sqrt(), star1, and star2.

◆ set()

Etoile_bin& Lorene::Binaire::set ( int  i)
inline

Read/write of the star no. i.

Definition at line 226 of file binaire.h.

References del_deriv(), and et.

◆ set_der_0x0()

void Lorene::Binaire::set_der_0x0 ( ) const
private

Sets to { 0x0} all the pointers on derived quantities.

Definition at line 198 of file binaire.C.

References p_angu_mom, p_ham_constr, p_mass_adm, p_mass_kom, p_mom_constr, p_total_ener, p_virial, p_virial_fus, and p_virial_gb.

◆ set_omega()

double& Lorene::Binaire::set_omega ( )
inline

Sets the orbital angular velocity [{ f_unit}].

Definition at line 232 of file binaire.h.

References omega.

◆ set_x_axe()

double& Lorene::Binaire::set_x_axe ( )
inline

Sets the absolute coordinate X of the rotation axis [{ r_unit}].

Definition at line 235 of file binaire.h.

References x_axe.

◆ total_ener()

double Lorene::Binaire::total_ener ( ) const

Total energy (excluding the rest mass energy).

In the Newtonian case, it is defined as the sum of kinetic, internal, and gravitational potential energies.

In the relativistic case, it is defined as $M_{ ADM} - M_{ bar,1} - M_{ bar,2}$.

Definition at line 292 of file binaire_global.C.

References et, Lorene::flat_scalar_prod(), Lorene::Etoile::get_ener(), Lorene::Etoile::get_logn_auto(), Lorene::Etoile_bin::get_logn_comp(), Lorene::Etoile::get_nbar(), Lorene::Etoile::get_u_euler(), Lorene::Cmp::integrale(), Lorene::Etoile::is_relativistic(), mass_adm(), Lorene::Etoile_bin::mass_b(), p_total_ener, star1, and star2.

◆ virial()

double Lorene::Binaire::virial ( ) const

◆ virial_fus()

double Lorene::Binaire::virial_fus ( ) const

◆ virial_gb()

◆ write_global()

Friends And Related Function Documentation

◆ operator<<

ostream& operator<< ( ostream &  ost,
const Binaire bibi 
)
friend

Save in a file.

Definition at line 254 of file binaire.C.

Member Data Documentation

◆ et

Etoile_bin* Lorene::Binaire::et[2]
private

Array of the two stars (to perform loops on the stars): { et[0]} contains the address of { star1} and { et[1]} that of { star2}.

Definition at line 127 of file binaire.h.

◆ omega

double Lorene::Binaire::omega
private

Angular velocity with respect to an asymptotically inertial observer.

Definition at line 132 of file binaire.h.

◆ p_angu_mom

Tbl* Lorene::Binaire::p_angu_mom
mutableprivate

Total angular momentum of the system.

Definition at line 148 of file binaire.h.

◆ p_ham_constr

double* Lorene::Binaire::p_ham_constr
mutableprivate

Relative error on the Hamiltonian constraint.

Definition at line 163 of file binaire.h.

◆ p_mass_adm

double* Lorene::Binaire::p_mass_adm
mutableprivate

Total ADM mass of the system.

Definition at line 142 of file binaire.h.

◆ p_mass_kom

double* Lorene::Binaire::p_mass_kom
mutableprivate

Total Komar mass of the system.

Definition at line 145 of file binaire.h.

◆ p_mom_constr

Tbl* Lorene::Binaire::p_mom_constr
mutableprivate

Relative error on the momentum constraint.

Definition at line 166 of file binaire.h.

◆ p_total_ener

double* Lorene::Binaire::p_total_ener
mutableprivate

Total energy of the system.

Definition at line 151 of file binaire.h.

◆ p_virial

double* Lorene::Binaire::p_virial
mutableprivate

Virial theorem error.

Definition at line 154 of file binaire.h.

◆ p_virial_fus

double* Lorene::Binaire::p_virial_fus
mutableprivate

Virial theorem error by J.L.Friedman, K.Uryu, and M.Shibata.

Definition at line 160 of file binaire.h.

◆ p_virial_gb

double* Lorene::Binaire::p_virial_gb
mutableprivate

Virial theorem error by E.Gourgoulhon and S.Bonazzola.

Definition at line 157 of file binaire.h.

◆ ref_triad

const Base_vect_cart Lorene::Binaire::ref_triad
private

Cartesian triad of the absolute reference frame.

Definition at line 115 of file binaire.h.

◆ star1

Etoile_bin Lorene::Binaire::star1
private

First star of the system.

Definition at line 118 of file binaire.h.

◆ star2

Etoile_bin Lorene::Binaire::star2
private

Second star of the system.

Definition at line 121 of file binaire.h.

◆ x_axe

double Lorene::Binaire::x_axe
private

Absolute X coordinate of the rotation axis.

Definition at line 136 of file binaire.h.


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