Class for stars in binary system in eXtended Conformal Thin Sandwich formulation. More...
#include <star.h>
Public Member Functions | |
Star_bin_xcts (Map &mp_i, int nzet_i, const Eos &eos_i, bool irrot) | |
Standard constructor. | |
Star_bin_xcts (const Star_bin_xcts &) | |
Copy constructor. | |
Star_bin_xcts (Map &mp_i, const Eos &eos_i, FILE *fich) | |
Constructor from a file (see sauve(FILE* ) ). | |
virtual | ~Star_bin_xcts () |
Destructor. | |
void | operator= (const Star_bin_xcts &) |
Assignment to another Star_bin_xcts . | |
Scalar & | set_pot_centri () |
Read/write the centrifugal potential. | |
Scalar & | set_Psi_auto () |
Read/write the conformal factor ![]() | |
Scalar & | set_Psi_comp () |
Read/write the conformal factor ![]() | |
Scalar & | set_chi_auto () |
Read/write the conformal factor ![]() | |
Scalar & | set_chi_comp () |
Read/write the conformal factor ![]() | |
Vector & | set_beta_auto () |
Read/write of ![]() | |
Vector & | set_beta () |
Read/write of ![]() | |
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 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_Psi_auto () const |
Returns the scalar field ![]() | |
const Scalar & | get_Psi_comp () const |
Returns the scalar field ![]() | |
const Scalar & | get_chi_auto () const |
Returns the scalar field ![]() | |
const Scalar & | get_chi_comp () const |
Returns the scalar field ![]() | |
const Vector & | get_dcov_chi () const |
Returns the covariant derivative of ![]() | |
const Vector & | get_dcov_Psi () const |
Returns the covariant derivative of the conformal factor ![]() | |
const Scalar & | get_Psi () const |
Return the conformal factor ![]() | |
const Scalar & | get_chi () const |
Return the function ![]() | |
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 Sym_tensor & | get_haij_auto () const |
Returns the part of the extrinsic curvature tensor ![]() beta_auto . | |
const Sym_tensor & | get_haij_comp () const |
Returns the part of the extrinsic curvature tensor ![]() beta_comp . | |
const Scalar & | get_hacar_auto () const |
Returns the part of ![]() beta_auto . | |
const Scalar & | get_hacar_comp () const |
Returns the part of ![]() beta_comp . | |
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_xcts &comp) |
Computes metric coefficients from known potentials, when the companion is another star. | |
void | update_metric (const Star_bin_xcts &comp, const Star_bin_xcts &star_prev, double relax) |
Same as update_metric (const Star_bin_xcts& ) but with relaxation. | |
void | update_metric_der_comp (const Star_bin_xcts &comp) |
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 () |
Computes haij_auto and hacar_auto from beta_auto , nn and Psi . | |
void | equilibrium (double ent_c, int mermax, int mermax_potvit, int mermax_poisson, double relax_poisson, double relax_potvit, double thres_adapt, const Tbl &fact, const Tbl *pent_limit, Tbl &diff) |
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_xcts &star_prev, double relax_ent, double relax_met, int mer, int fmer_met) |
Performs a relaxation on ent , Psi_auto , chi_auto and beta_auto . | |
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 | chi_auto |
Scalar field ![]() | |
Scalar | chi_comp |
Scalar field ![]() | |
Scalar | Psi_auto |
Scalar field ![]() | |
Scalar | Psi_comp |
Scalar field ![]() | |
Scalar | Psi |
Total conformal factor ![]() | |
Scalar | chi |
Total function ![]() | |
Scalar | psi4 |
Conformal factor ![]() | |
Vector | w_beta |
Solution for the vector part of the vector Poisson equation for ![]() | |
Scalar | khi |
Solution for the scalar part of the vector Poisson equation for ![]() | |
Vector | dcov_Psi |
Covariant derivative of the conformal factor ![]() | |
Vector | dcov_chi |
Covariant derivative of the function ![]() | |
Metric_flat | flat |
Flat metric defined on the mapping (Spherical components with respect to the mapping of the star) . | |
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 | haij_auto |
Part of the extrinsic curvature tensor ![]() beta_auto . | |
Sym_tensor | haij_comp |
Part of the extrinsic curvature tensor ![]() beta_comp . | |
Scalar | hacar_auto |
Part of the scalar ![]() beta_auto , i.e. | |
Scalar | hacar_comp |
Part of the scalar ![]() beta_auto and beta_comp , i.e. | |
Scalar | ssjm1_chi |
Effective source at the previous step for the resolution of the Poisson equation for . | |
Scalar | ssjm1_psi |
Effective source at the previous step for the resolution of the Poisson equation for . | |
Scalar | ssjm1_khi |
Effective source at the previous step for the resolution of the Poisson equation for khi . | |
Vector | ssjm1_wbeta |
Effective source at the previous step for wbeta of the vector Poisson equation for wbeta (needed for the solution of the vector Poisson equation for the shift ![]() | |
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_xcts |
ostream & | operator<< (ostream &, const Star &) |
Display. |
Class for stars in binary system in eXtended Conformal Thin Sandwich formulation.
*** UNDER DEEVELOPMENT *** ()
A Star_bin_xcts
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 1087 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 |
Definition at line 77 of file star_bin_xcts.C.
References Star::beta, beta_auto, beta_comp, bsn, Tensor::change_triad(), Vector::change_triad(), chi, chi_auto, chi_comp, d_psi, dcov_chi, dcov_Psi, Map::flat_met_cart(), Star::gamma, Map::get_bvect_cart(), hacar_auto, hacar_comp, haij_auto, haij_comp, khi, loggam, Star::mp, pot_centri, Psi, psi0, Psi_auto, Psi_comp, set_der_0x0(), Tensor::set_etat_zero(), ssjm1_chi, ssjm1_khi, ssjm1_psi, ssjm1_wbeta, Star::stress_euler, Star::u_euler, w_beta, and wit_w.
Star_bin_xcts::Star_bin_xcts | ( | const Star_bin_xcts & | 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 199 of file star_bin_xcts.C.
References Star::beta, beta_comp, bsn, Tensor::change_triad(), Vector::change_triad(), chi_comp, d_psi, dcov_chi, dcov_Psi, Map::flat_met_cart(), Star::gam_euler, Star::gamma, Map::get_bvect_cart(), Map::get_mg(), hacar_auto, hacar_comp, haij_auto, haij_comp, irrotational, khi, loggam, Star::mp, pot_centri, psi0, Psi_comp, set_der_0x0(), Tensor::set_etat_zero(), Star::stress_euler, Star::u_euler, w_beta, and wit_w.
Star_bin_xcts::~Star_bin_xcts | ( | ) | [virtual] |
void Star_bin_xcts::del_deriv | ( | ) | const [protected, virtual] |
Deletes all the derived quantities.
Reimplemented from Star.
Definition at line 300 of file star_bin_xcts.C.
References p_xa_barycenter, and set_der_0x0().
void Star_bin_xcts::del_hydro_euler | ( | ) | [protected, virtual] |
Sets to ETATNONDEF
(undefined state) the hydrodynamical quantities relative to the Eulerian observer.
Reimplemented from Star.
Definition at line 317 of file star_bin_xcts.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_xcts::equilibrium | ( | double | ent_c, | |
int | mermax, | |||
int | mermax_potvit, | |||
int | mermax_poisson, | |||
double | relax_poisson, | |||
double | relax_potvit, | |||
double | thres_adapt, | |||
const Tbl & | fact, | |||
const Tbl * | pent_limit, | |||
Tbl & | diff | |||
) |
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 | |
fact | [input] 1-D Tbl for the input of some factors:
| |
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. | |
diff | [output] 1-D Tbl for the storage of some error indicators |
Definition at line 87 of file star_bin_equilibrium_xcts.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(), beta_auto, chi, chi_auto, chi_comp, contract(), dcov_chi, dcov_Psi, Tensor::dec_dzpuis(), Scalar::derive_con(), Tensor::derive_con(), diffrel(), Vector::divergence(), Scalar::dsdr(), Star::ener_euler, Star::ent, Star::equation_of_state(), flat, Map::get_bvect_cart(), Map::get_mg(), Mg3d::get_np(), Mg3d::get_nr(), Mg3d::get_nt(), Mg3d::get_nzone(), hacar_auto, hacar_comp, haij_auto, Map_et::homothetie(), hydro_euler(), irrotational, khi, Scalar::laplacian(), log(), loggam, max(), Star::mp, Star::nn, norme(), Star::nzet, Scalar::poisson(), pot_centri, pow(), Star::press, Psi, psi4, Psi_auto, Psi_comp, Star::ray_eq(), Star::ray_pole(), Map::reevaluate_symy(), Map::resize(), Star::s_euler, Vector::set(), Tenseur::set(), Tbl::set(), Tenseur::set_etat_qcq(), Tbl::set_etat_qcq(), Scalar::set_spectral_va(), Tenseur::set_std_base(), Valeur::smooth(), sqrt(), ssjm1_chi, ssjm1_khi, ssjm1_psi, ssjm1_wbeta, Scalar::std_spectral_base(), Star::u_euler, Scalar::val_grid_point(), Scalar::val_point(), Map::val_r(), velocity_potential(), and w_beta.
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_xcts::extrinsic_curvature | ( | ) |
Computes haij_auto
and hacar_auto
from beta_auto
, nn
and Psi
.
Definition at line 51 of file star_bin_extr_curv_xcts.C.
References beta_auto, chi, Metric_flat::con(), contract(), del_deriv(), Tensor::derive_con(), Vector::divergence(), flat, hacar_auto, haij_auto, pow(), Psi, Tensor::set(), Tensor::std_spectral_base(), and Tensor::up_down().
void Star_bin_xcts::fait_d_psi | ( | ) |
Computes the gradient of the total velocity potential .
Definition at line 627 of file star_bin_xcts.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_xcts::get_beta_auto | ( | ) | const [inline] |
const Vector& Star_bin_xcts::get_beta_comp | ( | ) | const [inline] |
const Vector& Star_bin_xcts::get_bsn | ( | ) | const [inline] |
const Scalar& Star_bin_xcts::get_chi | ( | ) | const [inline] |
const Scalar& Star_bin_xcts::get_chi_auto | ( | ) | const [inline] |
const Scalar& Star_bin_xcts::get_chi_comp | ( | ) | const [inline] |
const Vector& Star_bin_xcts::get_d_psi | ( | ) | const [inline] |
const Vector& Star_bin_xcts::get_dcov_chi | ( | ) | const [inline] |
const Vector& Star_bin_xcts::get_dcov_Psi | ( | ) | 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_xcts::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 Scalar& Star_bin_xcts::get_hacar_auto | ( | ) | const [inline] |
Returns the part of generated by
beta_auto
.
Definition at line 1413 of file star.h.
References hacar_auto.
const Scalar& Star_bin_xcts::get_hacar_comp | ( | ) | const [inline] |
Returns the part of generated by
beta_comp
.
Definition at line 1418 of file star.h.
References hacar_comp.
const Sym_tensor& Star_bin_xcts::get_haij_auto | ( | ) | const [inline] |
const Sym_tensor& Star_bin_xcts::get_haij_comp | ( | ) | const [inline] |
const Scalar& Star_bin_xcts::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 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_xcts::get_pot_centri | ( | ) | const [inline] |
const Scalar& Star::get_press | ( | ) | const [inline, inherited] |
const Scalar& Star_bin_xcts::get_Psi | ( | ) | const [inline] |
const Scalar& Star_bin_xcts::get_psi0 | ( | ) | const [inline] |
const Scalar& Star_bin_xcts::get_psi4 | ( | ) | const [inline] |
const Scalar& Star_bin_xcts::get_Psi_auto | ( | ) | const [inline] |
const Scalar& Star_bin_xcts::get_Psi_comp | ( | ) | 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 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_xcts::get_wit_w | ( | ) | const [inline] |
void Star_bin_xcts::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 45 of file star_bin_hydro_xcts.C.
References Tensor::annule(), 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, Star::nzet, 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_xcts::is_irrotational | ( | ) | const [inline] |
Returns true
for an irrotational motion, false
for a corotating one.
Definition at line 1313 of file star.h.
References irrotational.
void Star_bin_xcts::kinematics | ( | double | omega, | |
double | x_axe | |||
) |
Computes the quantities bsn
and pot_centri
.
The calculation is performed starting from the quantities nn
, beta
, Psi
, 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 55 of file star_bin_kinema_xcts.C.
References Tensor::annule_domain(), 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_xcts::mass_b | ( | ) | const [virtual] |
Baryon mass.
Implements Star.
Definition at line 46 of file star_bin_global_xcts.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_xcts::mass_g | ( | ) | const [virtual] |
Gravitational mass.
Implements Star.
Definition at line 70 of file star_bin_global_xcts.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_xcts::operator= | ( | const Star_bin_xcts & | star | ) |
Assignment to another Star_bin_xcts
.
Reimplemented from Star.
Definition at line 332 of file star_bin_xcts.C.
References beta_auto, beta_comp, bsn, chi, chi_auto, chi_comp, d_psi, dcov_chi, dcov_Psi, del_deriv(), flat, hacar_auto, hacar_comp, haij_auto, haij_comp, irrotational, khi, loggam, pot_centri, Psi, psi0, psi4, Psi_auto, Psi_comp, ssjm1_chi, ssjm1_khi, ssjm1_psi, ssjm1_wbeta, w_beta, and wit_w.
ostream & Star_bin_xcts::operator>> | ( | ostream & | ost | ) | const [protected, virtual] |
Operator >> (virtual function called by the operator <<).
Reimplemented from Star.
Definition at line 459 of file star_bin_xcts.C.
References Star::beta, beta_auto, bsn, d_psi, Scalar::dsdr(), Star::ener, Star::ent, Star::eos, Star::gam_euler, Map::get_ori_x(), hacar_auto, hacar_comp, haij_auto, haij_comp, irrotational, loggam, mass_b(), mass_g(), max(), min(), Star::mp, Star::nbar, Star::nn, Star::nzet, pow(), Star::press, Psi_auto, Psi_comp, Star::ray_eq(), Star::ray_eq_pi(), Star::ray_eq_pis2(), Star::ray_pole(), Scalar::std_spectral_base(), Star::u_euler, Scalar::val_grid_point(), Scalar::val_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_xcts::relaxation | ( | const Star_bin_xcts & | star_prev, | |
double | relax_ent, | |||
double | relax_met, | |||
int | mer, | |||
int | fmer_met | |||
) |
Performs a relaxation on ent
, Psi_auto
, chi_auto
and beta_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 685 of file star_bin_xcts.C.
References beta_auto, chi_auto, del_deriv(), Star::ent, Star::equation_of_state(), and Psi_auto.
void Star_bin_xcts::sauve | ( | FILE * | fich | ) | const [virtual] |
Save in a file.
Reimplemented from Star.
Definition at line 428 of file star_bin_xcts.C.
References beta_auto, chi_auto, Star::gam_euler, irrotational, khi, psi0, Psi_auto, Tensor::sauve(), Scalar::sauve(), ssjm1_chi, ssjm1_khi, ssjm1_psi, ssjm1_wbeta, and w_beta.
Vector & Star_bin_xcts::set_beta | ( | ) |
Read/write of .
Definition at line 414 of file star_bin_xcts.C.
References Star::beta, and del_deriv().
Vector & Star_bin_xcts::set_beta_auto | ( | ) |
Read/write of .
Definition at line 407 of file star_bin_xcts.C.
References beta_auto, and del_deriv().
Scalar & Star_bin_xcts::set_chi_auto | ( | ) |
Read/write the conformal factor .
Definition at line 400 of file star_bin_xcts.C.
References chi_auto, and del_deriv().
Scalar & Star_bin_xcts::set_chi_comp | ( | ) |
Read/write the conformal factor .
Definition at line 393 of file star_bin_xcts.C.
References chi_comp, and del_deriv().
void Star_bin_xcts::set_der_0x0 | ( | ) | const [protected, virtual] |
Sets to 0x0
all the pointers on derived quantities.
Reimplemented from Star.
Definition at line 309 of file star_bin_xcts.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().
Map& Star::set_mp | ( | ) | [inline, inherited] |
Scalar & Star_bin_xcts::set_pot_centri | ( | ) |
Read/write the centrifugal potential.
Definition at line 372 of file star_bin_xcts.C.
References del_deriv(), and pot_centri.
Scalar & Star_bin_xcts::set_Psi_auto | ( | ) |
Read/write the conformal factor .
Definition at line 386 of file star_bin_xcts.C.
References del_deriv(), and Psi_auto.
Scalar & Star_bin_xcts::set_Psi_comp | ( | ) |
Read/write the conformal factor .
Definition at line 379 of file star_bin_xcts.C.
References del_deriv(), and Psi_comp.
void Star_bin_xcts::update_metric | ( | const Star_bin_xcts & | comp, | |
const Star_bin_xcts & | star_prev, | |||
double | relax | |||
) |
Same as update_metric
(const Star_bin_xcts& ) 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 148 of file star_bin_upmetr_xcts.C.
References Star::beta, beta_auto, beta_comp, Vector::change_triad(), chi, chi_auto, chi_comp, Metric_flat::con(), del_deriv(), extrinsic_curvature(), flat, Star::gamma, Map::get_bvect_cart(), Tensor::get_triad(), Scalar::import(), log(), Star::logn, Star::mp, Star::nn, pow(), Psi, psi4, Psi_auto, Psi_comp, Vector::set(), Tensor::set_etat_qcq(), Scalar::set_etat_qcq(), Scalar::set_etat_zero(), Tensor::set_triad(), Vector::std_spectral_base(), and Scalar::std_spectral_base().
void Star_bin_xcts::update_metric | ( | const Star_bin_xcts & | comp | ) |
Computes metric coefficients from known potentials, when the companion is another star.
The calculation is performed starting from the quantities Psi_auto
, chi_auto
, beta_auto
, comp.Psi_auto
, comp.chi_auto
, comp.beta_auto
which are supposed to be up to date. From these, the following fields are updated: Psi_comp
, chi_comp
, beta_comp
, nn
, psi4
, beta
.
comp | companion star. | |
omega | angular velocity with respect to an asymptotically inertial observer |
Definition at line 63 of file star_bin_upmetr_xcts.C.
References Star::beta, beta_auto, beta_comp, Vector::change_triad(), chi, chi_auto, chi_comp, Metric_flat::con(), del_deriv(), extrinsic_curvature(), flat, Star::gamma, Map::get_bvect_cart(), Tensor::get_triad(), Scalar::import(), log(), Star::logn, Star::mp, Star::nn, pow(), Psi, psi4, Psi_auto, Psi_comp, Vector::set(), Tensor::set_etat_qcq(), Scalar::set_etat_qcq(), Scalar::set_etat_zero(), Tensor::set_triad(), Vector::std_spectral_base(), and Scalar::std_spectral_base().
void Star_bin_xcts::update_metric_der_comp | ( | const Star_bin_xcts & | comp | ) |
Computes the derivative of metric functions related to the companion star.
Definition at line 58 of file star_bin_upmetr_der_xcts.C.
References beta_auto, beta_comp, Tensor::change_triad(), Vector::change_triad(), chi, chi_auto, Metric_flat::con(), contract(), dcov_chi, dcov_Psi, Tensor::dec_dzpuis(), del_deriv(), Scalar::derive_cov(), Vector::divergence(), flat, Map::get_bvect_cart(), Tensor::get_triad(), hacar_comp, haij_auto, haij_comp, Scalar::import(), Tensor::inc_dzpuis(), Star::mp, pow(), Psi, Psi_auto, Tensor::set(), Vector::set(), Valeur::set_base(), Scalar::set_spectral_va(), Tensor::std_spectral_base(), and Tensor::up_down().
double Star_bin_xcts::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 66 of file star_bin_vel_pot_xcts.C.
References Param::add_double(), Param::add_int(), Param::add_int_mod(), Tensor::annule(), Scalar::annule(), bsn, Vector::change_triad(), chi, contract(), d_psi, Eos::der_nbar_ent(), Scalar::derive_con(), Scalar::derive_cov(), diffrel(), 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(), Scalar::laplacian(), log(), Star::mp, norme(), Star::nzet, Map::poisson_compact(), Psi, psi0, psi4, Vector::set(), Cmp::set(), Scalar::set_domain(), Tensor::set_etat_qcq(), Scalar::set_spectral_va(), Tensor::set_triad(), Scalar::std_spectral_base(), Cmp::va, Valeur::ylm(), and Valeur::ylm_i().
double Star_bin_xcts::xa_barycenter | ( | ) | const [virtual] |
Absolute coordinate X of the barycenter of the baryon density,.
Definition at line 95 of file star_bin_global_xcts.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_xcts::beta_auto [protected] |
Vector Star_bin_xcts::beta_comp [protected] |
Vector Star_bin_xcts::bsn [protected] |
Scalar Star_bin_xcts::chi [protected] |
Scalar Star_bin_xcts::chi_auto [protected] |
Scalar Star_bin_xcts::chi_comp [protected] |
Vector Star_bin_xcts::d_psi [protected] |
Vector Star_bin_xcts::dcov_chi [protected] |
Vector Star_bin_xcts::dcov_Psi [protected] |
Scalar Star::ener [protected, inherited] |
Scalar Star::ener_euler [protected, inherited] |
Metric_flat Star_bin_xcts::flat [protected] |
Scalar Star::gam_euler [protected, inherited] |
Metric Star::gamma [protected, inherited] |
Scalar Star_bin_xcts::hacar_auto [protected] |
Scalar Star_bin_xcts::hacar_comp [protected] |
Sym_tensor Star_bin_xcts::haij_auto [protected] |
Sym_tensor Star_bin_xcts::haij_comp [protected] |
bool Star_bin_xcts::irrotational [protected] |
Scalar Star_bin_xcts::khi [protected] |
Scalar Star_bin_xcts::loggam [protected] |
Scalar Star::logn [protected, inherited] |
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_xcts::p_xa_barycenter [mutable, protected] |
Tbl* Star::p_xi_surf [mutable, protected, inherited] |
Scalar Star_bin_xcts::pot_centri [protected] |
Scalar Star::press [protected, inherited] |
Scalar Star_bin_xcts::Psi [protected] |
Scalar Star_bin_xcts::psi0 [protected] |
Scalar Star_bin_xcts::psi4 [protected] |
Scalar Star_bin_xcts::Psi_auto [protected] |
Scalar Star_bin_xcts::Psi_comp [protected] |
Scalar Star::s_euler [protected, inherited] |
Scalar Star_bin_xcts::ssjm1_chi [protected] |
Scalar Star_bin_xcts::ssjm1_khi [protected] |
Scalar Star_bin_xcts::ssjm1_psi [protected] |
Vector Star_bin_xcts::ssjm1_wbeta [protected] |
Sym_tensor Star::stress_euler [protected, inherited] |
Vector Star::u_euler [protected, inherited] |
Vector Star_bin_xcts::w_beta [protected] |
Vector Star_bin_xcts::wit_w [protected] |