Class for stars in binary system. More...
#include <star.h>
Public Member Functions | |
Star_bin (Map &mp_i, int nzet_i, const Eos &eos_i, bool irrot, bool conf_flat) | |
Standard constructor. | |
Star_bin (const Star_bin &) | |
Copy constructor. | |
Star_bin (Map &mp_i, const Eos &eos_i, FILE *fich) | |
Constructor from a file (see sauve(FILE* ) ). | |
virtual | ~Star_bin () |
Destructor. | |
void | operator= (const Star_bin &) |
Assignment to another Star_bin . | |
Scalar & | set_pot_centri () |
Read/write the centrifugal potential. | |
Scalar & | set_logn_comp () |
Read/write of the logarithm of the lapse generated principally by the companion. | |
void | set_logn_auto (const Scalar &logn_auto_new) |
Assignment of a new logn_auto. | |
void | set_lnq_auto (const Scalar &lnq_auto_new) |
Assignment of a new lnq_auto. | |
Vector & | set_beta_auto () |
Read/write of ![]() | |
Vector & | set_beta () |
Read/write of ![]() | |
void | set_conf_flat (bool confflat) |
Write if conformally flat. | |
bool | is_irrotational () const |
Returns true for an irrotational motion, false for a corotating one. | |
const Scalar & | get_psi0 () const |
Returns the non-translational part of the velocity potential. | |
const Vector & | get_d_psi () const |
Returns the covariant derivative of the velocity potential (Spherical components with respect to the mapping of the star). | |
const Vector & | get_wit_w () const |
Returns the spatial projection of the fluid 3-velocity with respect to the co-orbiting observer. | |
const Scalar & | get_loggam () const |
Returns the logarithm of the Lorentz factor between the fluid and the co-orbiting observer. | |
const Vector & | get_bsn () const |
Returns the shift vector, divided by N, of the rotating coordinates, ![]() | |
const Scalar & | get_pot_centri () const |
Returns the centrifugal potential. | |
const Scalar & | get_logn_auto () const |
Returns the part of the lapse logarithm (gravitational potential at the Newtonian limit) generated principally by the star. | |
const Scalar & | get_logn_comp () const |
Returns the part of the lapse logarithm (gravitational potential at the Newtonian limit) generated principally by the companion star. | |
const Vector & | get_beta_auto () const |
Returns the part of the shift vector ![]() | |
const Vector & | get_beta_comp () const |
Returns the part of the shift vector ![]() | |
const Scalar & | get_lnq_auto () const |
Returns the part of the vector field ![]() | |
const Scalar & | get_lnq_comp () const |
Returns the part of the vector field ![]() | |
const Vector & | get_dcov_logn () const |
Returns the covariant derivative of ![]() | |
const Vector & | get_dcon_logn () const |
Returns the contravariant derivative of ![]() | |
const Vector & | get_dcov_phi () const |
Returns the covariant derivative of ![]() | |
const Vector & | get_dcon_phi () const |
Returns the contravariant derivative of ![]() | |
const Scalar & | get_psi4 () const |
Return the conformal factor ![]() | |
const Metric & | get_flat () const |
Return the flat metric defined on the mapping (Spherical components with respect to the mapping of the star). | |
const Metric & | get_gtilde () const |
Return the conformal 3-metric ![]() | |
const Sym_tensor & | get_hij () const |
Return the total deviation of the inverse conformal metric ![]() | |
const Sym_tensor & | get_hij_auto () const |
Return the deviation of the inverse conformal metric ![]() | |
const Sym_tensor & | get_hij_comp () const |
Return the deviation of the inverse conformal metric ![]() | |
const Sym_tensor & | get_tkij_auto () const |
Returns the part of the extrinsic curvature tensor ![]() beta_auto . | |
const Sym_tensor & | get_tkij_comp () const |
Returns the part of the extrinsic curvature tensor ![]() beta_comp . | |
const Scalar & | get_kcar_auto () const |
Returns the part of ![]() beta_auto . | |
const Scalar & | get_kcar_comp () const |
Returns the part of ![]() beta_comp . | |
const Scalar | get_decouple () const |
Returns the function used to construct beta_auto from beta . | |
bool | is_conf_flat () const |
virtual void | sauve (FILE *) const |
Save in a file. | |
virtual double | mass_b () const |
Baryon mass. | |
virtual double | mass_g () const |
Gravitational mass. | |
virtual double | xa_barycenter () const |
Absolute coordinate X of the barycenter of the baryon density,. | |
virtual void | hydro_euler () |
Computes the hydrodynamical quantities relative to the Eulerian observer from those in the fluid frame, as well as wit_w and loggam . | |
void | update_metric (const Star_bin &comp, double omega) |
Computes metric coefficients from known potentials, when the companion is another star. | |
void | update_metric (const Star_bin &comp, const Star_bin &star_prev, double relax, double omega) |
Same as update_metric (const Star_bin& ) but with relaxation. | |
void | update_metric_der_comp (const Star_bin &comp, double omega) |
Computes the derivative of metric functions related to the companion star. | |
void | kinematics (double omega, double x_axe) |
Computes the quantities bsn and pot_centri . | |
void | fait_d_psi () |
Computes the gradient of the total velocity potential ![]() | |
void | extrinsic_curvature (double omega) |
Computes tkij_auto and akcar_auto from beta_auto , nn and Q . | |
void | equilibrium (double ent_c, int mermax, int mermax_potvit, int mermax_poisson, double relax_poisson, double relax_potvit, double thres_adapt, Tbl &diff, double om) |
Computes an equilibrium configuration. | |
double | velocity_potential (int mermax, double precis, double relax) |
Computes the non-translational part of the velocity scalar potential ![]() | |
void | relaxation (const Star_bin &star_prev, double relax_ent, double relax_met, int mer, int fmer_met) |
Performs a relaxation on ent , logn_auto , lnq_auto , beta_auto and hij_auto . | |
void | test_K_Hi () const |
Test if the gauge conditions we impose are well satisfied. | |
void | helical (double omega) const |
Test of the helical symmetry. | |
Map & | set_mp () |
Read/write of the mapping. | |
void | set_enthalpy (const Scalar &) |
Assignment of the enthalpy field. | |
void | equation_of_state () |
Computes the proper baryon and energy density, as well as pressure from the enthalpy. | |
virtual void | equilibrium_spher (double ent_c, double precis=1.e-14, const Tbl *pent_limit=0x0) |
Computes a spherical static configuration. | |
const Map & | get_mp () const |
Returns the mapping. | |
int | get_nzet () const |
Returns the number of domains occupied by the star. | |
const Eos & | get_eos () const |
Returns the equation of state. | |
const Scalar & | get_ent () const |
Returns the enthalpy field. | |
const Scalar & | get_nbar () const |
Returns the proper baryon density. | |
const Scalar & | get_ener () const |
Returns the proper total energy density. | |
const Scalar & | get_press () const |
Returns the fluid pressure. | |
const Scalar & | get_ener_euler () const |
Returns the total energy density with respect to the Eulerian observer. | |
const Scalar & | get_s_euler () const |
Returns the trace of the stress tensor in the Eulerian frame. | |
const Scalar & | get_gam_euler () const |
Returns the Lorentz factor between the fluid and Eulerian observers. | |
const Vector & | get_u_euler () const |
Returns the fluid 3-velocity with respect to the Eulerian observer. | |
const Tensor & | get_stress_euler () const |
Returns the spatial part of the stress-energy tensor with respect to the Eulerian observer. | |
const Scalar & | get_logn () const |
Returns the logarithm of the lapse N. | |
const Scalar & | get_nn () const |
Returns the lapse function N. | |
const Vector & | get_beta () const |
Returns the shift vector ![]() | |
const Scalar & | get_lnq () const |
const Metric & | get_gamma () const |
Returns the 3-metric ![]() | |
double | ray_eq () const |
Coordinate radius at ![]() ![]() | |
double | ray_eq_pis2 () const |
Coordinate radius at ![]() ![]() | |
double | ray_eq_pi () const |
Coordinate radius at ![]() ![]() | |
double | ray_eq_3pis2 () const |
Coordinate radius at ![]() ![]() | |
double | ray_pole () const |
Coordinate radius at ![]() | |
virtual const Itbl & | l_surf () const |
Description of the stellar surface: returns a 2-D Itbl containing the values of the domain index l on the surface at the collocation points in ![]() | |
const Tbl & | xi_surf () const |
Description of the stellar surface: returns a 2-D Tbl containing the values of the radial coordinate ![]() ![]() | |
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. | |
virtual void | del_hydro_euler () |
Sets to ETATNONDEF (undefined state) the hydrodynamical quantities relative to the Eulerian observer. | |
virtual ostream & | operator>> (ostream &) const |
Operator >> (virtual function called by the operator <<). | |
Protected Attributes | |
bool | irrotational |
true for an irrotational star, false for a corotating one | |
Scalar | psi0 |
Scalar potential ![]() | |
Vector | d_psi |
Gradient of ![]() | |
Vector | wit_w |
Spatial projection of the fluid 3-velocity with respect to the co-orbiting observer. | |
Scalar | loggam |
Logarithm of the Lorentz factor between the fluid and the co-orbiting observer. | |
Vector | bsn |
3-vector shift, divided by N, of the rotating coordinates, ![]() | |
Scalar | pot_centri |
Centrifugal potential. | |
Scalar | logn_auto |
Part of the lapse logarithm (gravitational potential at the Newtonian limit) generated principally by the star. | |
Scalar | logn_comp |
Part of the lapse logarithm (gravitational potential at the Newtonian limit) generated principally by the companion star. | |
Vector | dcov_logn |
Covariant derivative of the total logarithm of the lapse. | |
Vector | dcon_logn |
Contravariant derivative of the total logarithm of the lapse. | |
Scalar | lnq_auto |
Scalar field ![]() | |
Scalar | lnq_comp |
Scalar field ![]() | |
Scalar | psi4 |
Conformal factor ![]() | |
Vector | dcov_phi |
Covariant derivative of the logarithm of the conformal factor. | |
Vector | dcon_phi |
Contravariant derivative of the logarithm of the conformal factor. | |
Metric_flat | flat |
Flat metric defined on the mapping (Spherical components with respect to the mapping of the star) . | |
Metric | gtilde |
Conformal metric ![]() | |
Vector | beta_auto |
Part of the shift vector generated principally by the star (Spherical components with respect to the mapping of the star). | |
Vector | beta_comp |
Part of the shift vector generated principally by the star (Spherical components with respect to the mapping of the star). | |
Sym_tensor | hij |
Total deviation of the inverse conformal metric ![]() | |
Sym_tensor | hij_auto |
Deviation of the inverse conformal metric ![]() | |
Sym_tensor | hij_comp |
Deviation of the inverse conformal metric ![]() | |
Sym_tensor | tkij_auto |
Part of the extrinsic curvature tensor ![]() beta_auto . | |
Sym_tensor | tkij_comp |
Part of the extrinsic curvature tensor ![]() beta_comp . | |
Scalar | kcar_auto |
Part of the scalar ![]() beta_auto , i.e. | |
Scalar | kcar_comp |
Part of the scalar ![]() beta_auto and beta_comp , i.e. | |
Scalar | ssjm1_logn |
Effective source at the previous step for the resolution of the Poisson equation for logn_auto . | |
Scalar | ssjm1_lnq |
Effective source at the previous step for the resolution of the Poisson equation for lnq_auto . | |
Scalar | ssjm1_khi |
Effective source at the previous step for the resolution of the Poisson equation for khi . | |
Vector | ssjm1_wbeta |
Scalar | ssjm1_h11 |
Effective source at the previous step for the resolution of the Poisson equation for h00_auto . | |
Scalar | ssjm1_h21 |
Effective source at the previous step for the resolution of the Poisson equation for h10_auto . | |
Scalar | ssjm1_h31 |
Effective source at the previous step for the resolution of the Poisson equation for h20_auto . | |
Scalar | ssjm1_h22 |
Effective source at the previous step for the resolution of the Poisson equation for h11_auto . | |
Scalar | ssjm1_h32 |
Effective source at the previous step for the resolution of the Poisson equation for h21_auto . | |
Scalar | ssjm1_h33 |
Effective source at the previous step for the resolution of the Poisson equation for h22_auto . | |
Scalar | decouple |
Function used to construct the part ![]() ![]() | |
bool | conf_flat |
true if the 3-metric is conformally flat, false for a more general metric. | |
double * | p_xa_barycenter |
Absolute coordinate X of the barycenter of the baryon density. | |
Map & | mp |
Mapping associated with the star. | |
int | nzet |
Number of domains of *mp occupied by the star. | |
const Eos & | eos |
Equation of state of the stellar matter. | |
Scalar | ent |
Log-enthalpy. | |
Scalar | nbar |
Baryon density in the fluid frame. | |
Scalar | ener |
Total energy density in the fluid frame. | |
Scalar | press |
Fluid pressure. | |
Scalar | ener_euler |
Total energy density in the Eulerian frame. | |
Scalar | s_euler |
Trace of the stress scalar in the Eulerian frame. | |
Scalar | gam_euler |
Lorentz factor between the fluid and Eulerian observers. | |
Vector | u_euler |
Fluid 3-velocity with respect to the Eulerian observer. | |
Sym_tensor | stress_euler |
Spatial part of the stress-energy tensor with respect to the Eulerian observer. | |
Scalar | logn |
Logarithm of the lapse N . | |
Scalar | nn |
Lapse function N . | |
Vector | beta |
Shift vector. | |
Scalar | lnq |
Metric | gamma |
3-metric | |
double * | p_ray_eq |
Coordinate radius at ![]() ![]() | |
double * | p_ray_eq_pis2 |
Coordinate radius at ![]() ![]() | |
double * | p_ray_eq_pi |
Coordinate radius at ![]() ![]() | |
double * | p_ray_eq_3pis2 |
Coordinate radius at ![]() ![]() | |
double * | p_ray_pole |
Coordinate radius at ![]() | |
Itbl * | p_l_surf |
Description of the stellar surface: 2-D Itbl containing the values of the domain index l on the surface at the collocation points in ![]() | |
Tbl * | p_xi_surf |
Description of the stellar surface: 2-D Tbl containing the values of the radial coordinate ![]() ![]() | |
double * | p_mass_b |
Baryon mass. | |
double * | p_mass_g |
Gravitational mass. | |
Friends | |
class | Binary |
ostream & | operator<< (ostream &, const Star &) |
Display. |
Class for stars in binary system.
*** UNDER DEEVELOPMENT *** ()
A Star_bin
can be construted in two states, represented by the bool
member irrotational:
(i) irrotational (i.e. the fluid motion is irrotational) or (ii) rigidly corotating with respect to the orbital motion (synchronized binary).
Definition at line 479 of file star.h.
Standard constructor.
mp_i | Mapping on which the star will be defined | |
nzet_i | Number of domains occupied by the star | |
eos_i | Equation of state of the stellar matter | |
irrot | should be true for an irrotational star, false for a corotating one | |
conf_flat | should be true for a conformally flat metric false for a general one |
Definition at line 111 of file star_bin.C.
References Star::beta, beta_auto, beta_comp, bsn, Tensor::change_triad(), Vector::change_triad(), d_psi, dcon_logn, dcon_phi, dcov_logn, dcov_phi, Map::flat_met_cart(), Star::gamma, Map::get_bvect_cart(), hij, hij_auto, hij_comp, kcar_auto, kcar_comp, lnq_auto, lnq_comp, loggam, logn_auto, logn_comp, Star::mp, pot_centri, psi0, psi4, set_der_0x0(), Tensor::set_etat_zero(), ssjm1_h11, ssjm1_h21, ssjm1_h22, ssjm1_h31, ssjm1_h32, ssjm1_h33, ssjm1_khi, ssjm1_lnq, ssjm1_logn, Star::stress_euler, tkij_auto, tkij_comp, Star::u_euler, and wit_w.
Star_bin::Star_bin | ( | const Star_bin & | star | ) |
Constructor from a file (see sauve(FILE* )
).
mp_i | Mapping on which the star will be defined | |
eos_i | Equation of state of the stellar matter | |
fich | input file (must have been created by the function sauve ) |
Definition at line 254 of file star_bin.C.
References Star::beta, beta_comp, bsn, Tensor::change_triad(), Vector::change_triad(), conf_flat, d_psi, dcon_logn, dcon_phi, dcov_logn, dcov_phi, Map::flat_met_cart(), Star::gam_euler, Star::gamma, Map::get_bvect_cart(), Map::get_mg(), hij, hij_comp, irrotational, kcar_auto, kcar_comp, lnq_comp, loggam, logn_comp, Star::mp, pot_centri, psi0, psi4, set_der_0x0(), Tensor::set_etat_zero(), Star::stress_euler, tkij_auto, tkij_comp, Star::u_euler, and wit_w.
Star_bin::~Star_bin | ( | ) | [virtual] |
void Star_bin::del_deriv | ( | ) | const [protected, virtual] |
Deletes all the derived quantities.
Reimplemented from Star.
Definition at line 365 of file star_bin.C.
References p_xa_barycenter, and set_der_0x0().
void Star_bin::del_hydro_euler | ( | ) | [protected, virtual] |
Sets to ETATNONDEF
(undefined state) the hydrodynamical quantities relative to the Eulerian observer.
Reimplemented from Star.
Definition at line 385 of file star_bin.C.
References del_deriv().
void Star::equation_of_state | ( | ) | [inherited] |
Computes the proper baryon and energy density, as well as pressure from the enthalpy.
Reimplemented in Gravastar.
Definition at line 458 of file star.C.
References Param::add_int(), Scalar::allocate_all(), Star::del_deriv(), Star::ener, Eos::ener_ent(), Star::ent, Star::eos, Mg3d::get_grille3d(), Map::get_mg(), Mg3d::get_np(), Mg3d::get_nr(), Mg3d::get_nt(), Mg3d::get_nzone(), Star::mp, Star::nbar, Eos::nbar_ent(), Star::nzet, Star::press, Eos::press_ent(), Mtbl::set(), Scalar::set_domain(), Scalar::set_etat_qcq(), Tbl::set_etat_qcq(), Mtbl::set_etat_qcq(), Scalar::std_spectral_base(), Mtbl::t, and Grille3d::x.
void Star_bin::equilibrium | ( | double | ent_c, | |
int | mermax, | |||
int | mermax_potvit, | |||
int | mermax_poisson, | |||
double | relax_poisson, | |||
double | relax_potvit, | |||
double | thres_adapt, | |||
Tbl & | diff, | |||
double | om | |||
) |
Computes an equilibrium configuration.
ent_c | [input] Central enthalpy | |
mermax | [input] Maximum number of steps | |
mermax_poisson | [input] Maximum number of steps in poisson scalar | |
relax_poisson | [input] Relaxation factor in poisson scalar | |
mermax_potvit | [input] Maximum number of steps in Map_radial::poisson_compact | |
relax_potvit | [input] Relaxation factor in Map_radial::poisson_compact | |
thres_adapt | [input] Threshold on dH/dr for the adaptation of the mapping | |
diff | [output] 1-D Tbl for the storage of some error indicators |
Definition at line 139 of file star_bin_equilibrium.C.
References abs(), Map::adapt(), Param::add_cmp_mod(), Param::add_double(), Param::add_int(), Param::add_int_mod(), Param::add_tbl(), Param::add_tenseur_mod(), Scalar::annule_hard(), Cmp::annule_hard(), Star::beta, beta_auto, Metric_flat::con(), Metric::con(), conf_flat, contract(), Metric::cov(), dcon_logn, dcon_phi, dcov_logn, dcov_phi, Scalar::dec_dzpuis(), Tensor::dec_dzpuis(), Scalar::derive_con(), Tensor::derive_con(), Tensor_sym::derive_con(), Tensor_sym::derive_cov(), Tensor::derive_cov(), Scalar::derive_cov(), Sym_tensor::derive_lie(), diffrel(), Sym_tensor::divergence(), Vector::divergence(), Scalar::dsdr(), Star::ener_euler, Star::ent, Star::equation_of_state(), exp(), flat, Map::get_bvect_cart(), Connection::get_delta(), Map::get_mg(), Mg3d::get_np(), Mg3d::get_nr(), Mg3d::get_nt(), Mg3d::get_nzone(), Map::get_ori_x(), Map::get_rot_phi(), gtilde, hij, hij_auto, Map_et::homothetie(), hydro_euler(), Tenseur::inc_dzpuis(), Tensor::inc_dzpuis(), irrotational, kcar_auto, kcar_comp, Scalar::laplacian(), lnq_auto, loggam, Star::logn, logn_auto, logn_comp, max(), Star::mp, Star::nn, norme(), Star::nzet, Scalar::poisson(), pot_centri, pow(), Star::press, psi4, Map::r, Star::ray_eq(), Star::ray_eq_pi(), Star::ray_pole(), Map::reevaluate_symy(), Map::resize(), Star::s_euler, Coord::set(), Tensor::set(), Itbl::set(), Vector::set(), Tenseur::set(), Tbl::set(), Scalar::set_domain(), Tenseur::set_etat_qcq(), Tbl::set_etat_qcq(), Cmp::set_etat_qcq(), Scalar::set_outer_boundary(), Scalar::set_spectral_va(), Tenseur::set_std_base(), sqrt(), ssjm1_h11, ssjm1_h21, ssjm1_h22, ssjm1_h31, ssjm1_h32, ssjm1_h33, ssjm1_khi, ssjm1_lnq, ssjm1_logn, Scalar::std_spectral_base(), Star::stress_euler, tkij_auto, tkij_comp, Star::u_euler, Scalar::val_grid_point(), Scalar::val_point(), Map::val_r(), velocity_potential(), Map::xa, Map::ya, and Map::za.
void Star::equilibrium_spher | ( | double | ent_c, | |
double | precis = 1.e-14 , |
|||
const Tbl * | pent_limit = 0x0 | |||
) | [virtual, inherited] |
Computes a spherical static configuration.
ent_c | [input] central value of the enthalpy | |
precis | [input] threshold in the relative difference between the enthalpy fields of two consecutive steps to stop the iterative procedure (default value: 1.e-14) | |
ent_limit | [input] : array of enthalpy values to be set at the boundaries between the domains; if set to 0x0 (default), the initial values will be kept. |
Definition at line 91 of file star_equil_spher.C.
References Map_et::adapt(), Param::add_double(), Param::add_int(), Param::add_int_mod(), Param::add_tbl(), Scalar::annule(), diffrel(), Scalar::dsdr(), Map_af::dsdr(), Star::ener, Star::ener_euler, Star::ent, Star::equation_of_state(), exp(), Star::gam_euler, Star::gamma, Map_et::get_alpha(), Map_af::get_alpha(), Map_et::get_beta(), Map_af::get_beta(), Map::get_bvect_cart(), Map::get_mg(), Mg3d::get_nr(), Mg3d::get_nt(), Mg3d::get_nzone(), Map_af::homothetie(), Scalar::integrale(), Star::logn, Star::mass_b(), Star::mass_g(), Star::mp, Star::nn, norme(), Star::nzet, Map_af::poisson(), Star::press, Star::s_euler, Vector::set(), Map_af::set_alpha(), Map_af::set_beta(), Scalar::set_dzpuis(), Cmp::set_etat_qcq(), Scalar::set_etat_zero(), sqrt(), Scalar::std_spectral_base(), Star::u_euler, Scalar::val_grid_point(), and Map::val_r().
void Star_bin::extrinsic_curvature | ( | double | omega | ) |
Computes tkij_auto
and akcar_auto
from beta_auto
, nn
and Q
.
omega | angular velocity with respect to an asymptotically inertial observer |
Definition at line 71 of file star_bin_extr_curv.C.
References Scalar::annule_hard(), beta_auto, Metric::con(), contract(), del_deriv(), Tensor::derive_con(), Sym_tensor::derive_lie(), Vector::divergence(), Map::get_bvect_cart(), Map::get_mg(), Mg3d::get_nzone(), Map::get_rot_phi(), gtilde, hij_auto, kcar_auto, Star::mp, Star::nn, Tensor::set(), Coord::set(), Mg3d::std_base_vect_cart(), tkij_auto, Tensor::up_down(), Map::xa, and Map::ya.
void Star_bin::fait_d_psi | ( | ) |
Computes the gradient of the total velocity potential .
Definition at line 641 of file star_bin.C.
References Tensor::annule(), Scalar::annule_hard(), bsn, Vector::change_triad(), d_psi, Scalar::derive_cov(), Star::ent, exp(), flat, Star::gam_euler, Map::get_bvect_cart(), Map::get_mg(), Mg3d::get_nzone(), irrotational, Star::mp, Star::nzet, psi0, psi4, Vector::set(), Tensor::set(), Valeur::set_base(), Tensor::set_etat_nondef(), Vector::std_spectral_base(), and Cmp::va.
const Vector& Star::get_beta | ( | ) | const [inline, inherited] |
const Vector& Star_bin::get_beta_auto | ( | ) | const [inline] |
const Vector& Star_bin::get_beta_comp | ( | ) | const [inline] |
const Vector& Star_bin::get_bsn | ( | ) | const [inline] |
const Vector& Star_bin::get_d_psi | ( | ) | const [inline] |
const Vector& Star_bin::get_dcon_logn | ( | ) | const [inline] |
const Vector& Star_bin::get_dcon_phi | ( | ) | const [inline] |
const Vector& Star_bin::get_dcov_logn | ( | ) | const [inline] |
const Vector& Star_bin::get_dcov_phi | ( | ) | const [inline] |
const Scalar Star_bin::get_decouple | ( | ) | const [inline] |
const Scalar& Star::get_ener | ( | ) | const [inline, inherited] |
Returns the proper total energy density.
Definition at line 366 of file star.h.
References Star::ener.
const Scalar& Star::get_ener_euler | ( | ) | const [inline, inherited] |
Returns the total energy density with respect to the Eulerian observer.
Definition at line 372 of file star.h.
References Star::ener_euler.
const Scalar& Star::get_ent | ( | ) | const [inline, inherited] |
const Eos& Star::get_eos | ( | ) | const [inline, inherited] |
const Metric& Star_bin::get_flat | ( | ) | const [inline] |
const Scalar& Star::get_gam_euler | ( | ) | const [inline, inherited] |
Returns the Lorentz factor between the fluid and Eulerian observers.
Definition at line 378 of file star.h.
References Star::gam_euler.
const Metric& Star::get_gamma | ( | ) | const [inline, inherited] |
const Metric& Star_bin::get_gtilde | ( | ) | const [inline] |
const Sym_tensor& Star_bin::get_hij | ( | ) | const [inline] |
const Sym_tensor& Star_bin::get_hij_auto | ( | ) | const [inline] |
const Sym_tensor& Star_bin::get_hij_comp | ( | ) | const [inline] |
const Scalar& Star_bin::get_kcar_auto | ( | ) | const [inline] |
const Scalar& Star_bin::get_kcar_comp | ( | ) | const [inline] |
const Scalar& Star_bin::get_lnq_auto | ( | ) | const [inline] |
const Scalar& Star_bin::get_lnq_comp | ( | ) | const [inline] |
const Scalar& Star_bin::get_loggam | ( | ) | const [inline] |
const Scalar& Star::get_logn | ( | ) | const [inline, inherited] |
Returns the logarithm of the lapse N.
In the Newtonian case, this is the Newtonian gravitational potential (in units of ).
Definition at line 392 of file star.h.
References Star::logn.
const Scalar& Star_bin::get_logn_auto | ( | ) | const [inline] |
const Scalar& Star_bin::get_logn_comp | ( | ) | const [inline] |
const Map& Star::get_mp | ( | ) | const [inline, inherited] |
const Scalar& Star::get_nbar | ( | ) | const [inline, inherited] |
const Scalar& Star::get_nn | ( | ) | const [inline, inherited] |
int Star::get_nzet | ( | ) | const [inline, inherited] |
Returns the number of domains occupied by the star.
Definition at line 354 of file star.h.
References Star::nzet.
const Scalar& Star_bin::get_pot_centri | ( | ) | const [inline] |
const Scalar& Star::get_press | ( | ) | const [inline, inherited] |
const Scalar& Star_bin::get_psi0 | ( | ) | const [inline] |
const Scalar& Star_bin::get_psi4 | ( | ) | const [inline] |
const Scalar& Star::get_s_euler | ( | ) | const [inline, inherited] |
Returns the trace of the stress tensor in the Eulerian frame.
Definition at line 375 of file star.h.
References Star::s_euler.
const Tensor& Star::get_stress_euler | ( | ) | const [inline, inherited] |
Returns the spatial part of the stress-energy tensor with respect to the Eulerian observer.
Definition at line 386 of file star.h.
References Star::stress_euler.
const Sym_tensor& Star_bin::get_tkij_auto | ( | ) | const [inline] |
const Sym_tensor& Star_bin::get_tkij_comp | ( | ) | const [inline] |
const Vector& Star::get_u_euler | ( | ) | const [inline, inherited] |
Returns the fluid 3-velocity with respect to the Eulerian observer.
Definition at line 381 of file star.h.
References Star::u_euler.
const Vector& Star_bin::get_wit_w | ( | ) | const [inline] |
void Star_bin::helical | ( | double | omega | ) | const |
Test of the helical symmetry.
void Star_bin::hydro_euler | ( | ) | [virtual] |
Computes the hydrodynamical quantities relative to the Eulerian observer from those in the fluid frame, as well as wit_w
and loggam
.
The calculation is performed starting from the quantities ent
, ener
, press
and bsn
, which are supposed to be up to date. From these, the following fields are updated: gam_euler
, u_euler
, ener_euler
, s_euler
, stress_euler
, wit_w
and loggam
.
Reimplemented from Star.
Definition at line 62 of file star_bin_hydro.C.
References Tensor::annule_domain(), bsn, Metric::con(), contract(), Metric::cov(), d_psi, del_deriv(), Star::ener, Star::ener_euler, Star::ent, exp(), Star::gam_euler, Star::gamma, Map::get_bvect_cart(), Map::get_mg(), Mg3d::get_nzone(), irrotational, log(), loggam, Star::mp, norme(), Star::press, Star::s_euler, Tensor::set(), Scalar::set_dzpuis(), Tensor::set_etat_zero(), sqrt(), Vector::std_spectral_base(), Scalar::std_spectral_base(), Star::stress_euler, Star::u_euler, and wit_w.
bool Star_bin::is_irrotational | ( | ) | const [inline] |
Returns true
for an irrotational motion, false
for a corotating one.
Definition at line 769 of file star.h.
References irrotational.
void Star_bin::kinematics | ( | double | omega, | |
double | x_axe | |||
) |
Computes the quantities bsn
and pot_centri
.
The calculation is performed starting from the quantities nn
, beta
, Q
, which are supposed to be up to date.
omega | angular velocity with respect to an asymptotically inertial observer | |
x_axe | absolute X coordinate of the rotation axis |
Definition at line 71 of file star_bin_kinema.C.
References Scalar::annule(), Tensor::annule(), Star::beta, bsn, Vector::change_triad(), contract(), Metric::cov(), del_deriv(), Star::gamma, Map::get_bvect_cart(), Map::get_mg(), Mg3d::get_nzone(), Map::get_rot_phi(), log(), Star::mp, Star::nn, pot_centri, Vector::set(), sqrt(), Scalar::std_spectral_base(), Vector::std_spectral_base(), Map::xa, and Map::ya.
const Itbl & Star::l_surf | ( | ) | const [virtual, inherited] |
Description of the stellar surface: returns a 2-D Itbl
containing the values of the domain index l on the surface at the collocation points in .
The stellar surface is defined as the location where the enthalpy (member ent
) vanishes.
Reimplemented in Star_rot.
Definition at line 59 of file star_global.C.
References Star::ent, Map::get_mg(), Mg3d::get_np(), Mg3d::get_nt(), Scalar::get_spectral_va(), Star::mp, Star::nzet, Star::p_l_surf, and Star::p_xi_surf.
double Star_bin::mass_b | ( | ) | const [virtual] |
Baryon mass.
Implements Star.
Definition at line 63 of file star_bin_global.C.
References Metric::determinant(), Star::gam_euler, Star::gamma, Scalar::integrale(), Star::nbar, Star::p_mass_b, sqrt(), and Scalar::std_spectral_base().
double Star_bin::mass_g | ( | ) | const [virtual] |
Gravitational mass.
Implements Star.
Definition at line 89 of file star_bin_global.C.
References Metric::determinant(), Star::ener_euler, Star::gamma, Scalar::integrale(), Star::nn, Star::p_mass_g, Star::s_euler, sqrt(), and Scalar::std_spectral_base().
void Star_bin::operator= | ( | const Star_bin & | star | ) |
Assignment to another Star_bin
.
Reimplemented from Star.
Definition at line 400 of file star_bin.C.
References beta_auto, beta_comp, bsn, conf_flat, d_psi, dcon_logn, dcon_phi, dcov_logn, dcov_phi, decouple, del_deriv(), flat, gtilde, hij, hij_auto, hij_comp, irrotational, kcar_auto, kcar_comp, lnq_auto, lnq_comp, loggam, logn_auto, logn_comp, pot_centri, psi0, psi4, ssjm1_h11, ssjm1_h21, ssjm1_h22, ssjm1_h31, ssjm1_h32, ssjm1_h33, ssjm1_khi, ssjm1_lnq, ssjm1_logn, tkij_auto, tkij_comp, and wit_w.
ostream & Star_bin::operator>> | ( | ostream & | ost | ) | const [protected, virtual] |
Operator >> (virtual function called by the operator <<).
Reimplemented from Star.
Definition at line 518 of file star_bin.C.
References Star::beta, beta_auto, bsn, d_psi, Star::gam_euler, Map::get_ori_x(), irrotational, kcar_auto, kcar_comp, loggam, logn_auto, logn_comp, max(), min(), Star::mp, Star::ray_eq(), Star::ray_eq_pi(), tkij_auto, tkij_comp, Star::u_euler, Scalar::val_grid_point(), Map::val_r(), wit_w, and xa_barycenter().
double Star::ray_eq | ( | ) | const [inherited] |
Coordinate radius at ,
[r_unit].
Definition at line 104 of file star_global.C.
References Map::get_mg(), Mg3d::get_nt(), Mg3d::get_type_p(), Mg3d::get_type_t(), Star::l_surf(), Star::mp, Star::p_ray_eq, Map::val_r(), and Star::xi_surf().
double Star::ray_eq_3pis2 | ( | ) | const [inherited] |
Coordinate radius at ,
[r_unit].
Definition at line 229 of file star_global.C.
References Map::get_mg(), Mg3d::get_np(), Mg3d::get_nt(), Mg3d::get_type_p(), Mg3d::get_type_t(), Star::l_surf(), Star::mp, Star::p_ray_eq_3pis2, Star::ray_eq_pis2(), Map::val_r(), and Star::xi_surf().
double Star::ray_eq_pi | ( | ) | const [inherited] |
Coordinate radius at ,
[r_unit].
Definition at line 182 of file star_global.C.
References Map::get_mg(), Mg3d::get_np(), Mg3d::get_nt(), Mg3d::get_type_p(), Mg3d::get_type_t(), Star::l_surf(), Star::mp, Star::p_ray_eq_pi, Star::ray_eq(), Map::val_r(), and Star::xi_surf().
double Star::ray_eq_pis2 | ( | ) | const [inherited] |
Coordinate radius at ,
[r_unit].
Definition at line 134 of file star_global.C.
References Map::get_mg(), Mg3d::get_np(), Mg3d::get_nt(), Mg3d::get_type_p(), Mg3d::get_type_t(), Star::l_surf(), Star::mp, Star::p_ray_eq_pis2, Map::val_r(), and Star::xi_surf().
double Star::ray_pole | ( | ) | const [inherited] |
Coordinate radius at [r_unit].
Definition at line 274 of file star_global.C.
References Map::get_mg(), Mg3d::get_type_t(), Star::l_surf(), Star::mp, Star::p_ray_pole, Map::val_r(), and Star::xi_surf().
void Star_bin::relaxation | ( | const Star_bin & | star_prev, | |
double | relax_ent, | |||
double | relax_met, | |||
int | mer, | |||
int | fmer_met | |||
) |
Performs a relaxation on ent
, logn_auto
, lnq_auto
, beta_auto
and hij_auto
.
star_prev | [input] star at the previous step. | |
relax_ent | [input] Relaxation factor for ent | |
relax_met | [input] Relaxation factor for logn_auto , lnq_auto , beta_auto , only if (mer % fmer_met == 0). | |
mer | [input] Step number | |
fmer_met | [input] Step interval between metric updates |
Definition at line 697 of file star_bin.C.
References beta_auto, del_deriv(), Star::ent, Star::equation_of_state(), hij_auto, lnq_auto, and logn_auto.
void Star_bin::sauve | ( | FILE * | fich | ) | const [virtual] |
Save in a file.
Reimplemented from Star.
Definition at line 485 of file star_bin.C.
References beta_auto, conf_flat, Star::gam_euler, hij_auto, irrotational, lnq_auto, logn_auto, psi0, Tensor_sym::sauve(), Tensor::sauve(), Scalar::sauve(), ssjm1_h11, ssjm1_h21, ssjm1_h22, ssjm1_h31, ssjm1_h32, ssjm1_h33, ssjm1_khi, ssjm1_lnq, and ssjm1_logn.
Vector & Star_bin::set_beta | ( | ) |
Vector & Star_bin::set_beta_auto | ( | ) |
void Star_bin::set_conf_flat | ( | bool | confflat | ) | [inline] |
void Star_bin::set_der_0x0 | ( | ) | const [protected, virtual] |
Sets to 0x0
all the pointers on derived quantities.
Reimplemented from Star.
Definition at line 377 of file star_bin.C.
References p_xa_barycenter.
void Star::set_enthalpy | ( | const Scalar & | ent_i | ) | [inherited] |
Assignment of the enthalpy field.
Definition at line 375 of file star.C.
References Star::del_deriv(), Star::ent, and Star::equation_of_state().
void Star_bin::set_lnq_auto | ( | const Scalar & | lnq_auto_new | ) | [inline] |
void Star_bin::set_logn_auto | ( | const Scalar & | logn_auto_new | ) | [inline] |
Scalar & Star_bin::set_logn_comp | ( | ) |
Read/write of the logarithm of the lapse generated principally by the companion.
Definition at line 457 of file star_bin.C.
References del_deriv(), and logn_comp.
Map& Star::set_mp | ( | ) | [inline, inherited] |
Scalar & Star_bin::set_pot_centri | ( | ) |
Read/write the centrifugal potential.
Definition at line 450 of file star_bin.C.
References del_deriv(), and pot_centri.
void Star_bin::test_K_Hi | ( | ) | const |
Test if the gauge conditions we impose are well satisfied.
Definition at line 723 of file star_bin.C.
References Metric::con(), Tensor_sym::derive_cov(), Sym_tensor::divergence(), flat, Map::get_mg(), Mg3d::get_np(), Mg3d::get_nr(), Mg3d::get_nt(), gtilde, Star::mp, and norme().
void Star_bin::update_metric | ( | const Star_bin & | comp, | |
const Star_bin & | star_prev, | |||
double | relax, | |||
double | omega | |||
) |
Same as update_metric
(const Star_bin& ) but with relaxation.
comp | companion star. | |
star_prev | previous value of the star. | |
relax | relaxation parameter. | |
omega | angular velocity with respect to an asymptotically inertial observer |
Definition at line 272 of file star_bin_upmetr.C.
References Star::beta, beta_auto, beta_comp, Tensor::change_triad(), Vector::change_triad(), Metric_flat::con(), conf_flat, del_deriv(), Metric::determinant(), exp(), extrinsic_curvature(), flat, Star::gamma, Map::get_bvect_cart(), Map::get_mg(), Star::get_mp(), Mg3d::get_np(), Mg3d::get_nr(), Mg3d::get_nt(), Mg3d::get_nzone(), Tensor::get_triad(), gtilde, hij, hij_auto, hij_comp, Scalar::import(), lnq_auto, lnq_comp, Star::logn, logn_auto, logn_comp, Star::mp, Star::nn, psi4, Tensor::set(), Vector::set(), Tensor::set_etat_qcq(), Scalar::set_etat_qcq(), Tensor::set_etat_zero(), Scalar::set_etat_zero(), Tensor::set_triad(), Tensor::std_spectral_base(), Vector::std_spectral_base(), and Scalar::std_spectral_base().
void Star_bin::update_metric | ( | const Star_bin & | comp, | |
double | omega | |||
) |
Computes metric coefficients from known potentials, when the companion is another star.
The calculation is performed starting from the quantities logn_auto
, lnq_auto
, beta_auto
, hij_auto
, comp.logn_auto
, comp.lnq_auto
, comp.beta_auto
, comp.hij_auto
which are supposed to be up to date. From these, the following fields are updated: logn_comp
, lnq_comp
, beta_comp
, hij_comp
, nn
, psi4
, beta
,
comp | companion star. | |
omega | angular velocity with respect to an asymptotically inertial observer |
Definition at line 89 of file star_bin_upmetr.C.
References Star::beta, beta_auto, beta_comp, Tensor::change_triad(), Vector::change_triad(), Metric_flat::con(), conf_flat, del_deriv(), Metric::determinant(), exp(), extrinsic_curvature(), flat, Star::gamma, Map::get_bvect_cart(), Map::get_mg(), Star::get_mp(), Mg3d::get_np(), Mg3d::get_nr(), Mg3d::get_nt(), Mg3d::get_nzone(), Tensor::get_triad(), gtilde, hij, hij_auto, hij_comp, Scalar::import(), lnq_auto, lnq_comp, Star::logn, logn_auto, logn_comp, Star::mp, Star::nn, psi4, Tensor::set(), Vector::set(), Tensor::set_etat_qcq(), Scalar::set_etat_qcq(), Tensor::set_etat_zero(), Scalar::set_etat_zero(), Tensor::set_triad(), Tensor::std_spectral_base(), Vector::std_spectral_base(), and Scalar::std_spectral_base().
void Star_bin::update_metric_der_comp | ( | const Star_bin & | comp, | |
double | omega | |||
) |
Computes the derivative of metric functions related to the companion star.
omega | angular velocity with respect to an asymptotically inertial observer |
Definition at line 97 of file star_bin_upmetr_der.C.
References Scalar::annule_hard(), beta_comp, Vector::change_triad(), Metric::con(), contract(), dcon_logn, dcon_phi, dcov_logn, dcov_phi, Tensor::dec_dzpuis(), del_deriv(), Tensor::derive_con(), Scalar::derive_con(), Tensor::derive_cov(), Scalar::derive_cov(), Sym_tensor::derive_lie(), Vector::divergence(), flat, Map::get_bvect_cart(), get_flat(), Map::get_mg(), Tensor::get_mp(), Mg3d::get_nzone(), Map::get_rot_phi(), Tensor::get_triad(), gtilde, hij_comp, Scalar::import(), Tensor::inc_dzpuis(), kcar_comp, lnq_auto, logn_auto, Star::mp, Star::nn, Tensor::set(), Coord::set(), Vector::set(), Valeur::set_base(), Scalar::set_spectral_va(), Mg3d::std_base_vect_cart(), tkij_auto, tkij_comp, Tensor::up_down(), Map::xa, and Map::ya.
double Star_bin::velocity_potential | ( | int | mermax, | |
double | precis, | |||
double | relax | |||
) |
Computes the non-translational part of the velocity scalar potential by solving the continuity equation.
mermax | [input] Maximum number of steps in the iteration | |
precis | [input] Required precision: the iteration will be stopped when the relative difference on ![]() precis . | |
relax | [input] Relaxation factor. |
Definition at line 72 of file star_bin_vel_pot.C.
References Param::add_double(), Param::add_int(), Param::add_int_mod(), Tensor::annule(), Scalar::annule(), bsn, Vector::change_triad(), contract(), d_psi, Eos::der_nbar_ent(), Scalar::derive_con(), Tensor::derive_cov(), Scalar::derive_cov(), diffrel(), Tensor::down(), Star::ent, Star::eos, exp(), flat, Map::flat_met_spher(), Star::gam_euler, Map::get_bvect_cart(), Map::get_bvect_spher(), Scalar::get_etat(), Map::get_mg(), Mg3d::get_nzone(), Tensor::get_triad(), hij, Tensor::inc_dzpuis(), Scalar::laplacian(), Star::mp, norme(), Star::nzet, Map::poisson_compact(), psi0, psi4, Vector::set(), Cmp::set(), Scalar::set_domain(), Scalar::set_spectral_va(), Tensor::set_triad(), Scalar::std_spectral_base(), Cmp::va, Valeur::ylm(), and Valeur::ylm_i().
double Star_bin::xa_barycenter | ( | ) | const [virtual] |
Absolute coordinate X of the barycenter of the baryon density,.
Definition at line 114 of file star_bin_global.C.
References Tensor::annule_domain(), Metric::determinant(), Star::gam_euler, Star::gamma, Map::get_mg(), Mg3d::get_nzone(), Scalar::integrale(), mass_b(), Star::mp, Star::nbar, p_xa_barycenter, sqrt(), Scalar::std_spectral_base(), and Map::xa.
const Tbl & Star::xi_surf | ( | ) | const [inherited] |
Description of the stellar surface: returns a 2-D Tbl
containing the values of the radial coordinate on the surface at the collocation points in
.
The stellar surface is defined as the location where the enthalpy (member ent
) vanishes.
Definition at line 85 of file star_global.C.
References Star::l_surf(), Star::p_l_surf, and Star::p_xi_surf.
ostream& operator<< | ( | ostream & | , | |
const Star & | ||||
) | [friend, inherited] |
Display.
Vector Star::beta [protected, inherited] |
Vector Star_bin::beta_auto [protected] |
Vector Star_bin::beta_comp [protected] |
Vector Star_bin::bsn [protected] |
bool Star_bin::conf_flat [protected] |
Vector Star_bin::d_psi [protected] |
Vector Star_bin::dcon_logn [protected] |
Vector Star_bin::dcon_phi [protected] |
Vector Star_bin::dcov_logn [protected] |
Vector Star_bin::dcov_phi [protected] |
Scalar Star_bin::decouple [protected] |
Scalar Star::ener [protected, inherited] |
Scalar Star::ener_euler [protected, inherited] |
Metric_flat Star_bin::flat [protected] |
Scalar Star::gam_euler [protected, inherited] |
Metric Star::gamma [protected, inherited] |
Metric Star_bin::gtilde [protected] |
Sym_tensor Star_bin::hij [protected] |
Sym_tensor Star_bin::hij_auto [protected] |
Sym_tensor Star_bin::hij_comp [protected] |
bool Star_bin::irrotational [protected] |
Scalar Star_bin::kcar_auto [protected] |
Scalar Star_bin::kcar_comp [protected] |
Scalar Star_bin::lnq_auto [protected] |
Scalar Star_bin::lnq_comp [protected] |
Scalar Star_bin::loggam [protected] |
Scalar Star::logn [protected, inherited] |
Scalar Star_bin::logn_auto [protected] |
Scalar Star_bin::logn_comp [protected] |
Scalar Star::nbar [protected, inherited] |
int Star::nzet [protected, inherited] |
Itbl* Star::p_l_surf [mutable, protected, inherited] |
double* Star::p_mass_b [mutable, protected, inherited] |
double* Star::p_mass_g [mutable, protected, inherited] |
double* Star::p_ray_eq [mutable, protected, inherited] |
double* Star::p_ray_eq_3pis2 [mutable, protected, inherited] |
double* Star::p_ray_eq_pi [mutable, protected, inherited] |
double* Star::p_ray_eq_pis2 [mutable, protected, inherited] |
double* Star::p_ray_pole [mutable, protected, inherited] |
double* Star_bin::p_xa_barycenter [mutable, protected] |
Tbl* Star::p_xi_surf [mutable, protected, inherited] |
Scalar Star_bin::pot_centri [protected] |
Scalar Star::press [protected, inherited] |
Scalar Star_bin::psi0 [protected] |
Scalar Star_bin::psi4 [protected] |
Scalar Star::s_euler [protected, inherited] |
Scalar Star_bin::ssjm1_h11 [protected] |
Scalar Star_bin::ssjm1_h21 [protected] |
Scalar Star_bin::ssjm1_h22 [protected] |
Scalar Star_bin::ssjm1_h31 [protected] |
Scalar Star_bin::ssjm1_h32 [protected] |
Scalar Star_bin::ssjm1_h33 [protected] |
Scalar Star_bin::ssjm1_khi [protected] |
Scalar Star_bin::ssjm1_lnq [protected] |
Scalar Star_bin::ssjm1_logn [protected] |
Sym_tensor Star::stress_euler [protected, inherited] |
Sym_tensor Star_bin::tkij_auto [protected] |
Sym_tensor Star_bin::tkij_comp [protected] |
Vector Star::u_euler [protected, inherited] |
Vector Star_bin::wit_w [protected] |