LORENE
|
Base class for stars. More...
#include <star.h>
Public Member Functions | |
Star (Map &mp_i, int nzet_i, const Eos &eos_i) | |
Standard constructor. More... | |
Star (const Star &) | |
Copy constructor. More... | |
Star (Map &mp_i, const Eos &eos_i, FILE *fich) | |
Constructor from a file (see sauve(FILE* ) ). More... | |
virtual | ~Star () |
Destructor. More... | |
void | operator= (const Star &) |
Assignment to another Star . More... | |
Map & | set_mp () |
Read/write of the mapping. More... | |
void | set_enthalpy (const Scalar &) |
Assignment of the enthalpy field. More... | |
void | equation_of_state () |
Computes the proper baryon and energy density, as well as pressure from the enthalpy. More... | |
virtual void | hydro_euler () |
Computes the hydrodynamical quantities relative to the Eulerian observer from those in the fluid frame (nbar , ener and press ). More... | |
virtual void | equilibrium_spher (double ent_c, double precis=1.e-14, const Tbl *pent_limit=0x0) |
Computes a spherical static configuration. More... | |
const Map & | get_mp () const |
Returns the mapping. More... | |
int | get_nzet () const |
Returns the number of domains occupied by the star. More... | |
const Eos & | get_eos () const |
Returns the equation of state. More... | |
const Scalar & | get_ent () const |
Returns the enthalpy field. More... | |
const Scalar & | get_nbar () const |
Returns the proper baryon density. More... | |
const Scalar & | get_ener () const |
Returns the proper total energy density. More... | |
const Scalar & | get_press () const |
Returns the fluid pressure. More... | |
const Scalar & | get_ener_euler () const |
Returns the total energy density with respect to the Eulerian observer. More... | |
const Scalar & | get_s_euler () const |
Returns the trace of the stress tensor in the Eulerian frame. More... | |
const Scalar & | get_gam_euler () const |
Returns the Lorentz factor between the fluid and Eulerian observers. More... | |
const Vector & | get_u_euler () const |
Returns the fluid 3-velocity with respect to the Eulerian observer. More... | |
const Tensor & | get_stress_euler () const |
Returns the spatial part of the stress-energy tensor with respect to the Eulerian observer. More... | |
const Scalar & | get_logn () const |
Returns the logarithm of the lapse N. More... | |
const Scalar & | get_nn () const |
Returns the lapse function N. More... | |
const Vector & | get_beta () const |
Returns the shift vector . More... | |
const Scalar & | get_lnq () const |
const Metric & | get_gamma () const |
Returns the 3-metric . More... | |
virtual void | sauve (FILE *) const |
Save in a file. More... | |
double | ray_eq () const |
Coordinate radius at , [r_unit]. More... | |
double | ray_eq_pis2 () const |
Coordinate radius at , [r_unit]. More... | |
double | ray_eq_pi () const |
Coordinate radius at , [r_unit]. More... | |
double | ray_eq_3pis2 () const |
Coordinate radius at , [r_unit]. More... | |
double | ray_pole () const |
Coordinate radius at [r_unit]. More... | |
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 . More... | |
const Tbl & | xi_surf () const |
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 . More... | |
virtual double | mass_b () const =0 |
Baryon mass. More... | |
virtual double | mass_g () const =0 |
Gravitational mass. More... | |
Protected Member Functions | |
virtual void | del_deriv () const |
Deletes all the derived quantities. More... | |
virtual void | set_der_0x0 () const |
Sets to 0x0 all the pointers on derived quantities. More... | |
virtual void | del_hydro_euler () |
Sets to ETATNONDEF (undefined state) the hydrodynamical quantities relative to the Eulerian observer. More... | |
virtual ostream & | operator>> (ostream &) const |
Operator >> (virtual function called by the operator <<). More... | |
Protected Attributes | |
Map & | mp |
Mapping associated with the star. More... | |
int | nzet |
Number of domains of *mp occupied by the star. More... | |
const Eos & | eos |
Equation of state of the stellar matter. More... | |
Scalar | ent |
Log-enthalpy. More... | |
Scalar | nbar |
Baryon density in the fluid frame. More... | |
Scalar | ener |
Total energy density in the fluid frame. More... | |
Scalar | press |
Fluid pressure. More... | |
Scalar | ener_euler |
Total energy density in the Eulerian frame. More... | |
Scalar | s_euler |
Trace of the stress scalar in the Eulerian frame. More... | |
Scalar | gam_euler |
Lorentz factor between the fluid and Eulerian observers. More... | |
Vector | u_euler |
Fluid 3-velocity with respect to the Eulerian observer. More... | |
Sym_tensor | stress_euler |
Spatial part of the stress-energy tensor with respect to the Eulerian observer. More... | |
Scalar | logn |
Logarithm of the lapse N . More... | |
Scalar | nn |
Lapse function N . More... | |
Vector | beta |
Shift vector. More... | |
Scalar | lnq |
Metric | gamma |
3-metric More... | |
double * | p_ray_eq |
Coordinate radius at , . More... | |
double * | p_ray_eq_pis2 |
Coordinate radius at , . More... | |
double * | p_ray_eq_pi |
Coordinate radius at , . More... | |
double * | p_ray_eq_3pis2 |
Coordinate radius at , . More... | |
double * | p_ray_pole |
Coordinate radius at . More... | |
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 . More... | |
Tbl * | p_xi_surf |
Description of the stellar surface: 2-D Tbl containing the values of the radial coordinate on the surface at the collocation points in . More... | |
double * | p_mass_b |
Baryon mass. More... | |
double * | p_mass_g |
Gravitational mass. More... | |
Friends | |
ostream & | operator<< (ostream &, const Star &) |
Display. More... | |
Base class for stars.
()
A Star
is constructed upon (i) a mapping (derived class of Map
), the center of which defines the center of the star, and (ii) an equation of state (derived class of Eos
).
It contains tensor fields (class Tensor
) which describe the hydrodynamical quantities as well as the gravitational field (spacetime metric).
According to the 3+1 formalism, the spacetime metric is written
where is the 3-metric, described by a Lorene object of class Metric
.
The 3+1 formalism introduces two kinds of privileged observers: the fluid comoving observer and the Eulerian observer, whose 4-velocity is the unit future directed normal to the t = const hypersurfaces. The hydrodynamical quantities measured by the fluid observer correspond to the members ent
, nbar
, ener
, and press
. The hydrodynamical quantities measured by the Eulerian observer correspond to the members ener_euler
, s_euler
, gam_euler
, and u_euler
.
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 |
Definition at line 127 of file star.C.
References beta, Lorene::Metric::cov(), ener, ener_euler, ent, Lorene::Map::flat_met_spher(), gam_euler, gamma, logn, mp, nbar, nn, press, s_euler, set_der_0x0(), Lorene::Tensor::set_etat_zero(), Lorene::Scalar::std_spectral_base(), stress_euler, and u_euler.
Lorene::Star::Star | ( | const Star & | et | ) |
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 226 of file star.C.
References beta, ent, eos, Lorene::Eos::eos_from_file(), Lorene::fread_be(), Lorene::Map::get_mg(), mp, nn, nzet, set_der_0x0(), Lorene::Tensor::set_etat_zero(), stress_euler, and u_euler.
|
virtual |
|
protectedvirtual |
Deletes all the derived quantities.
Reimplemented in Lorene::Star_bin_xcts, Lorene::Star_bin, Lorene::Star_rot, Lorene::Star_bhns, Lorene::Star_rot_Dirac, and Lorene::Star_rot_CFC.
Definition at line 301 of file star.C.
References p_l_surf, p_mass_b, p_mass_g, p_ray_eq, p_ray_eq_3pis2, p_ray_eq_pi, p_ray_eq_pis2, p_ray_pole, p_xi_surf, and set_der_0x0().
|
protectedvirtual |
Sets to ETATNONDEF
(undefined state) the hydrodynamical quantities relative to the Eulerian observer.
Reimplemented in Lorene::Star_bin_xcts, Lorene::Star_bin, Lorene::Star_rot, Lorene::Star_rot_Dirac, and Lorene::Star_rot_CFC.
Definition at line 333 of file star.C.
References del_deriv(), ener_euler, gam_euler, s_euler, Lorene::Scalar::set_etat_nondef(), Lorene::Tensor::set_etat_nondef(), stress_euler, and u_euler.
void Lorene::Star::equation_of_state | ( | ) |
Computes the proper baryon and energy density, as well as pressure from the enthalpy.
Definition at line 465 of file star.C.
References Lorene::Param::add_int(), Lorene::Scalar::allocate_all(), del_deriv(), ener, Lorene::Eos::ener_ent(), ent, eos, Lorene::Mg3d::get_grille3d(), Lorene::Map::get_mg(), Lorene::Mg3d::get_np(), Lorene::Mg3d::get_nr(), Lorene::Mg3d::get_nt(), Lorene::Mg3d::get_nzone(), mp, nbar, Lorene::Eos::nbar_ent(), nzet, press, Lorene::Eos::press_ent(), Lorene::Mtbl::set(), Lorene::Scalar::set_domain(), Lorene::Mtbl::set_etat_qcq(), Lorene::Tbl::set_etat_qcq(), Lorene::Scalar::set_etat_qcq(), Lorene::Scalar::std_spectral_base(), Lorene::Mtbl::t, and Lorene::Grille3d::x.
|
virtual |
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 101 of file star_equil_spher.C.
References Lorene::Star_rot::a_car, Lorene::Map_et::adapt(), Lorene::Param::add_double(), Lorene::Param::add_int(), Lorene::Param::add_int_mod(), Lorene::Param::add_tbl(), Lorene::Scalar::annule(), Lorene::Star_rot::b_car, Lorene::Star_rot::bbb, Lorene::diffrel(), Lorene::Scalar::dsdr(), Lorene::Map_af::dsdr(), Lorene::Star_rot::dzeta, ener, ener_euler, ent, equation_of_state(), Lorene::exp(), gam_euler, gamma, Lorene::Map_af::get_alpha(), Lorene::Map_et::get_alpha(), Lorene::Map_af::get_beta(), Lorene::Map_et::get_beta(), Lorene::Map::get_bvect_cart(), Lorene::Map::get_mg(), Lorene::Mg3d::get_nr(), Lorene::Mg3d::get_nt(), Lorene::Mg3d::get_nzone(), Lorene::Map_af::homothetie(), Lorene::Scalar::integrale(), logn, mass_b(), mass_g(), mp, nn, Lorene::norme(), nzet, Lorene::Map_af::poisson(), press, s_euler, Lorene::Vector::set(), Lorene::Map_af::set_alpha(), Lorene::Map_af::set_beta(), Lorene::Scalar::set_dzpuis(), Lorene::Cmp::set_etat_qcq(), Lorene::Tensor::set_etat_zero(), Lorene::sqrt(), Lorene::Scalar::std_spectral_base(), u_euler, Lorene::Scalar::val_grid_point(), and Lorene::Map::val_r().
|
inline |
|
inline |
|
inline |
Returns the total energy density with respect to the Eulerian observer.
Definition at line 376 of file star.h.
References ener_euler.
|
inline |
|
inline |
|
inline |
|
inline |
|
inline |
|
inline |
|
inline |
|
inline |
|
inline |
|
inline |
|
inline |
|
inline |
Returns the spatial part of the stress-energy tensor with respect to the Eulerian observer.
Definition at line 390 of file star.h.
References stress_euler.
|
inline |
|
virtual |
Computes the hydrodynamical quantities relative to the Eulerian observer from those in the fluid frame (nbar
, ener
and press
).
Reimplemented in Lorene::Star_bin_xcts, Lorene::Star_bin, Lorene::Star_rot, Lorene::Star_rot_Dirac, Lorene::Star_rot_CFC, and Lorene::Star_rot_Dirac_diff.
|
virtual |
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 Lorene::Star_rot.
Definition at line 66 of file star_global.C.
References ent, Lorene::Map::get_mg(), Lorene::Mg3d::get_np(), Lorene::Mg3d::get_nt(), Lorene::Scalar::get_spectral_va(), mp, nzet, p_l_surf, and p_xi_surf.
|
pure virtual |
Baryon mass.
Implemented in Lorene::Star_bin_xcts, Lorene::Star_bin, Lorene::Star_rot, Lorene::Star_bhns, Lorene::Star_rot_Dirac, and Lorene::Star_rot_CFC.
Definition at line 310 of file star_global.C.
References p_mass_b.
|
pure virtual |
Gravitational mass.
Implemented in Lorene::Star_bin_xcts, Lorene::Star_bin, Lorene::Star_rot, Lorene::Star_bhns, Lorene::Star_rot_Dirac, and Lorene::Star_rot_CFC.
Definition at line 330 of file star_global.C.
References p_mass_g.
void Lorene::Star::operator= | ( | const Star & | et | ) |
|
protectedvirtual |
Operator >> (virtual function called by the operator <<).
Reimplemented in Lorene::Star_bin_xcts, Lorene::Star_bin, Lorene::Star_rot, Lorene::Star_bhns, Lorene::Star_rot_Dirac, Lorene::Star_rot_CFC, Lorene::Gravastar, and Lorene::Star_rot_Dirac_diff.
Definition at line 420 of file star.C.
References ener, ent, eos, mass_b(), mass_g(), nbar, nn, nzet, press, ray_eq(), ray_eq_pi(), ray_eq_pis2(), ray_pole(), and Lorene::Scalar::val_grid_point().
double Lorene::Star::ray_eq | ( | ) | const |
Coordinate radius at , [r_unit].
Definition at line 111 of file star_global.C.
References Lorene::Map::get_mg(), Lorene::Mg3d::get_nt(), Lorene::Mg3d::get_type_p(), Lorene::Mg3d::get_type_t(), mp, and p_ray_eq.
double Lorene::Star::ray_eq_3pis2 | ( | ) | const |
Coordinate radius at , [r_unit].
Definition at line 236 of file star_global.C.
References Lorene::Map::get_mg(), Lorene::Mg3d::get_np(), Lorene::Mg3d::get_nt(), Lorene::Mg3d::get_type_p(), Lorene::Mg3d::get_type_t(), mp, and p_ray_eq_3pis2.
double Lorene::Star::ray_eq_pi | ( | ) | const |
Coordinate radius at , [r_unit].
Definition at line 189 of file star_global.C.
References Lorene::Map::get_mg(), Lorene::Mg3d::get_np(), Lorene::Mg3d::get_nt(), Lorene::Mg3d::get_type_p(), Lorene::Mg3d::get_type_t(), mp, and p_ray_eq_pi.
double Lorene::Star::ray_eq_pis2 | ( | ) | const |
Coordinate radius at , [r_unit].
Definition at line 141 of file star_global.C.
References Lorene::Map::get_mg(), Lorene::Mg3d::get_np(), Lorene::Mg3d::get_nt(), Lorene::Mg3d::get_type_p(), Lorene::Mg3d::get_type_t(), mp, and p_ray_eq_pis2.
double Lorene::Star::ray_pole | ( | ) | const |
Coordinate radius at [r_unit].
Definition at line 281 of file star_global.C.
References Lorene::Map::get_mg(), Lorene::Mg3d::get_type_t(), mp, and p_ray_pole.
|
virtual |
Save in a file.
Reimplemented in Lorene::Star_bin_xcts, Lorene::Star_bin, Lorene::Star_rot, Lorene::Star_bhns, Lorene::Star_rot_Dirac, Lorene::Star_rot_CFC, and Lorene::Star_rot_Dirac_diff.
Definition at line 400 of file star.C.
References ent, eos, Lorene::fwrite_be(), logn, nzet, Lorene::Eos::sauve(), and Lorene::Scalar::sauve().
|
protectedvirtual |
Sets to 0x0
all the pointers on derived quantities.
Reimplemented in Lorene::Star_bin_xcts, Lorene::Star_bin, Lorene::Star_rot, Lorene::Star_bhns, Lorene::Star_rot_Dirac, and Lorene::Star_rot_CFC.
Definition at line 319 of file star.C.
References p_l_surf, p_mass_b, p_mass_g, p_ray_eq, p_ray_eq_3pis2, p_ray_eq_pi, p_ray_eq_pis2, p_ray_pole, and p_xi_surf.
void Lorene::Star::set_enthalpy | ( | const Scalar & | ent_i | ) |
Assignment of the enthalpy field.
Definition at line 382 of file star.C.
References del_deriv(), ent, and equation_of_state().
|
inline |
const Tbl & Lorene::Star::xi_surf | ( | ) | const |
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 92 of file star_global.C.
|
friend |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
mutableprotected |
|
mutableprotected |
|
mutableprotected |
|
mutableprotected |
|
mutableprotected |
|
mutableprotected |
|
mutableprotected |
|
mutableprotected |
|
protected |
|
protected |
|
protected |