LORENE
|
Base class for stars *** DEPRECATED : use class Star
instead ***.
More...
#include <etoile.h>
Public Member Functions | |
Etoile (Map &mp_i, int nzet_i, bool relat, const Eos &eos_i) | |
Standard constructor. More... | |
Etoile (const Etoile &) | |
Copy constructor. More... | |
Etoile (Map &mp_i, const Eos &eos_i, FILE *fich) | |
Constructor from a file (see sauve(FILE*) ). More... | |
virtual | ~Etoile () |
Destructor. More... | |
void | operator= (const Etoile &) |
Assignment to another Etoile . More... | |
Map & | set_mp () |
Read/write of the mapping. More... | |
void | set_enthalpy (const Cmp &) |
Assignment of the enthalpy field. More... | |
virtual 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 *ent_limit=0x0) |
Computes a spherical static configuration. More... | |
void | equil_spher_regular (double ent_c, double precis=1.e-14) |
Computes a spherical static configuration. More... | |
virtual void | equil_spher_falloff (double ent_c, double precis=1.e-14) |
Computes a spherical static configuration with the outer boundary condition at a finite radius. More... | |
const Map & | get_mp () const |
Returns the mapping. More... | |
int | get_nzet () const |
Returns the number of domains occupied by the star. More... | |
bool | is_relativistic () const |
Returns true for a relativistic star, false for a Newtonian one. More... | |
const Eos & | get_eos () const |
Returns the equation of state. More... | |
const Tenseur & | get_ent () const |
Returns the enthalpy field. More... | |
const Tenseur & | get_nbar () const |
Returns the proper baryon density. More... | |
const Tenseur & | get_ener () const |
Returns the proper total energy density. More... | |
const Tenseur & | get_press () const |
Returns the fluid pressure. More... | |
const Tenseur & | get_ener_euler () const |
Returns the total energy density with respect to the Eulerian observer. More... | |
const Tenseur & | get_s_euler () const |
Returns the trace of the stress tensor in the Eulerian frame. More... | |
const Tenseur & | get_gam_euler () const |
Returns the Lorentz factor between the fluid and Eulerian observers. More... | |
const Tenseur & | get_u_euler () const |
Returns the fluid 3-velocity with respect to the Eulerian observer. More... | |
const Tenseur & | get_logn_auto () const |
Returns the logarithm of the part of the lapse N generated principaly by the star. More... | |
const Tenseur & | get_logn_auto_regu () const |
Returns the regular part of the logarithm of the part of the lapse N generated principaly by the star. More... | |
const Tenseur & | get_logn_auto_div () const |
Returns the divergent part of the logarithm of the part of the lapse N generated principaly by the star. More... | |
const Tenseur & | get_d_logn_auto_div () const |
Returns the gradient of logn_auto_div . More... | |
const Tenseur & | get_beta_auto () const |
Returns the logarithm of the part of the product AN generated principaly by the star. More... | |
const Tenseur & | get_nnn () const |
Returns the total lapse function N. More... | |
const Tenseur & | get_shift () const |
Returns the total shift vector . More... | |
const Tenseur & | get_a_car () const |
Returns the total conformal factor . 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... | |
double | ray_eq (int kk) 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 |
Baryon mass. More... | |
virtual double | mass_g () const |
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... | |
bool | relativistic |
Indicator of relativity: true for a relativistic star, false for a Newtonian one. More... | |
double | unsurc2 |
: unsurc2=1 for a relativistic star, 0 for a Newtonian one. More... | |
int | k_div |
Index of regularity of the gravitational potential logn_auto . More... | |
const Eos & | eos |
Equation of state of the stellar matter. More... | |
Tenseur | ent |
Log-enthalpy (relativistic case) or specific enthalpy (Newtonian case) More... | |
Tenseur | nbar |
Baryon density in the fluid frame. More... | |
Tenseur | ener |
Total energy density in the fluid frame. More... | |
Tenseur | press |
Fluid pressure. More... | |
Tenseur | ener_euler |
Total energy density in the Eulerian frame. More... | |
Tenseur | s_euler |
Trace of the stress tensor in the Eulerian frame. More... | |
Tenseur | gam_euler |
Lorentz factor between the fluid and Eulerian observers. More... | |
Tenseur | u_euler |
Fluid 3-velocity with respect to the Eulerian observer. More... | |
Tenseur | logn_auto |
Total of the logarithm of the part of the lapse N generated principaly by the star. More... | |
Tenseur | logn_auto_regu |
Regular part of the logarithm of the part of the lapse N generated principaly by the star. More... | |
Tenseur | logn_auto_div |
Divergent part (if k_div!=0 ) of the logarithm of the part of the lapse N generated principaly by the star. More... | |
Tenseur | d_logn_auto_div |
Gradient of logn_auto_div (if k_div!=0 ) More... | |
Tenseur | beta_auto |
Logarithm of the part of the product AN generated principaly by by the star. More... | |
Tenseur | nnn |
Total lapse function. More... | |
Tenseur | shift |
Total shift vector. More... | |
Tenseur | a_car |
Total conformal factor . 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 Etoile &) |
Display. More... | |
Base class for stars *** DEPRECATED : use class Star
instead ***.
()
An Etoile
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 Tenseur
) which describle the hydrodynamical quantities as well as the gravitational field (spacetime metric).
According to the 3+1 formalism, the spacetime metric is written
where is a 3-metric, the exact form of which is specified in the derived classes of Etoile
. The base class Etoile
by itself provides only storage for the lapse function N (member nnn
), the shift vector (member shift
) and the conformal factor (member a_car
).
The 3+1 formalism introduces two kinds of priviledged 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
.
A star of class Etoile
can be either relativistic or Newtonian, depending on the boolean indicator relativistic
. For a Newtonian star, the metric coefficients N and A are set to 1, and is set to zero; the only relevant gravitational quantity in this case is logn_auto
which represents the (regular part of) the Newtonian gravitational potential generated by the star.
Standard constructor.
mp_i | Mapping on which the star will be defined |
nzet_i | Number of domains occupied by the star |
relat | should be true for a relativistic star, false for a Newtonian one |
eos_i | Equation of state of the stellar matter |
Definition at line 146 of file etoile.C.
References a_car, beta_auto, d_logn_auto_div, ener, ener_euler, ent, eos, gam_euler, logn_auto, logn_auto_div, logn_auto_regu, nbar, nnn, press, relativistic, s_euler, set_der_0x0(), Lorene::Tenseur::set_std_base(), shift, u_euler, and unsurc2.
Lorene::Etoile::Etoile | ( | const Etoile & | 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 278 of file etoile.C.
References beta_auto, d_logn_auto_div, ent, eos, Lorene::Eos::eos_from_file(), Lorene::fread_be(), Lorene::Map::get_bvect_spher(), k_div, logn_auto, logn_auto_div, logn_auto_regu, mp, nzet, relativistic, set_der_0x0(), shift, and unsurc2.
|
virtual |
|
protectedvirtual |
Deletes all the derived quantities.
Reimplemented in Lorene::Etoile_rot, Lorene::Etoile_bin, Lorene::Et_rot_bifluid, and Lorene::Et_rot_mag.
Definition at line 381 of file etoile.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::Etoile_rot, Lorene::Etoile_bin, Lorene::Et_rot_bifluid, and Lorene::Et_rot_mag.
Definition at line 413 of file etoile.C.
References del_deriv(), ener_euler, gam_euler, s_euler, Lorene::Tenseur::set_etat_nondef(), and u_euler.
|
virtual |
Computes the proper baryon and energy density, as well as pressure from the enthalpy.
Reimplemented in Lorene::Et_magnetisation, and Lorene::Et_rot_bifluid.
Definition at line 569 of file etoile.C.
References Lorene::Param::add_int(), Lorene::Cmp::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::Tenseur::set(), Lorene::Cmp::set(), Lorene::Mtbl::set_etat_qcq(), Lorene::Tbl::set_etat_qcq(), Lorene::Tenseur::set_etat_qcq(), Lorene::Tenseur::set_std_base(), Lorene::Cmp::std_base_scal(), Lorene::Mtbl::t, and Lorene::Grille3d::x.
|
virtual |
Computes a spherical static configuration with the outer boundary condition at a finite radius.
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) |
Definition at line 60 of file etoile_eqsph_falloff.C.
References a_car, Lorene::Tenseur::annule(), beta_auto, Lorene::diffrel(), Lorene::Cmp::dsdr(), Lorene::Map_af::dsdr(), ener, ener_euler, ent, equation_of_state(), Lorene::exp(), gam_euler, Lorene::Map::get_mg(), Lorene::Mg3d::get_nr(), Lorene::Mg3d::get_nt(), Lorene::Mg3d::get_nzone(), Lorene::Map_af::homothetie(), logn_auto, mass_b(), mass_g(), mp, nbar, nnn, Lorene::norme(), nzet, press, relativistic, s_euler, Lorene::Tenseur::set(), Lorene::Tenseur::set_etat_qcq(), Lorene::Tenseur::set_std_base(), shift, Lorene::sqrt(), u_euler, unsurc2, and Lorene::Map::val_r().
void Lorene::Etoile::equil_spher_regular | ( | double | ent_c, |
double | precis = 1.e-14 |
||
) |
Computes a spherical static configuration.
The sources for Poisson equations are regularized by extracting analytical diverging parts.
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) |
Definition at line 118 of file et_equil_spher_regu.C.
References a_car, Lorene::Tenseur::annule(), beta_auto, d_logn_auto_div, Lorene::Eos::der_ener_ent_p(), Lorene::Eos::der_nbar_ent_p(), Lorene::diffrel(), Lorene::Map_af::dsdr(), ener, ener_euler, ent, eos, equation_of_state(), Lorene::exp(), gam_euler, Lorene::Map::get_bvect_spher(), Lorene::Map::get_mg(), Lorene::Mg3d::get_nr(), Lorene::Mg3d::get_nt(), Lorene::Mg3d::get_nzone(), Lorene::Tenseur::gradient_spher(), Lorene::Map_af::homothetie(), k_div, logn_auto, logn_auto_div, logn_auto_regu, mass_b(), mass_g(), mp, nbar, nnn, Lorene::norme(), nzet, Lorene::Map_af::poisson(), Lorene::Map_af::poisson_regular(), press, relativistic, s_euler, Lorene::Tenseur::set(), Lorene::Tenseur::set_etat_qcq(), Lorene::Tenseur::set_std_base(), Lorene::Tenseur::set_triad(), shift, Lorene::sqrt(), Lorene::Cmp::std_base_scal(), u_euler, unsurc2, and Lorene::Map::val_r().
|
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 90 of file etoile_equil_spher.C.
References a_car, Lorene::Map_et::adapt(), Lorene::Param::add_double(), Lorene::Param::add_int(), Lorene::Param::add_int_mod(), Lorene::Param::add_tbl(), Lorene::Tenseur::annule(), beta_auto, Lorene::diffrel(), Lorene::Cmp::dsdr(), Lorene::Map_af::dsdr(), ener, ener_euler, ent, equation_of_state(), Lorene::exp(), gam_euler, Lorene::Map_af::get_alpha(), Lorene::Map_et::get_alpha(), Lorene::Map_af::get_beta(), Lorene::Map_et::get_beta(), get_ent(), Lorene::Map::get_mg(), Lorene::Mg3d::get_nr(), Lorene::Mg3d::get_nt(), Lorene::Mg3d::get_nzone(), get_press(), Lorene::Map_af::homothetie(), logn_auto, mass_b(), mass_g(), mp, nbar, nnn, Lorene::norme(), nzet, Lorene::Map_af::poisson(), press, relativistic, s_euler, Lorene::Tenseur::set(), Lorene::Map_af::set_alpha(), Lorene::Map_af::set_beta(), Lorene::Tenseur::set_etat_qcq(), Lorene::Tenseur::set_std_base(), shift, Lorene::sqrt(), u_euler, unsurc2, and Lorene::Map::val_r().
|
inline |
|
inline |
|
inline |
Returns the gradient of logn_auto_div
.
Definition at line 722 of file etoile.h.
References d_logn_auto_div.
|
inline |
|
inline |
Returns the total energy density with respect to the Eulerian observer.
Definition at line 688 of file etoile.h.
References ener_euler.
|
inline |
|
inline |
|
inline |
|
inline |
|
inline |
Returns the divergent part of the logarithm of the part of the lapse N generated principaly by the star.
In the Newtonian case, this is the diverging part of the Newtonian gravitational potential (in units of ).
Definition at line 718 of file etoile.h.
References logn_auto_div.
|
inline |
Returns the regular part of the logarithm of the part of the lapse N generated principaly by the star.
In the Newtonian case, this is the Newtonian gravitational potential (in units of ).
Definition at line 711 of file etoile.h.
References logn_auto_regu.
|
inline |
|
inline |
|
inline |
|
inline |
|
inline |
|
inline |
|
inline |
|
inline |
|
virtual |
Computes the hydrodynamical quantities relative to the Eulerian observer from those in the fluid frame (nbar
, ener
and press
).
Reimplemented in Lorene::Etoile_rot, Lorene::Etoile_bin, Lorene::Et_rot_bifluid, and Lorene::Et_rot_diff.
|
inline |
Returns true
for a relativistic star, false
for a Newtonian one.
Definition at line 670 of file etoile.h.
References relativistic.
|
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::Etoile_rot, and Lorene::Et_rot_bifluid.
Definition at line 78 of file etoile_global.C.
References ent, Lorene::Map::get_mg(), Lorene::Mg3d::get_np(), Lorene::Mg3d::get_nt(), mp, nzet, p_l_surf, and p_xi_surf.
|
virtual |
Baryon mass.
Reimplemented in Lorene::Etoile_rot, Lorene::Etoile_bin, and Lorene::Et_rot_bifluid.
Definition at line 528 of file etoile_global.C.
References Lorene::Tenseur::get_etat(), nbar, p_mass_b, and relativistic.
|
virtual |
Gravitational mass.
Reimplemented in Lorene::Etoile_rot, Lorene::Etoile_bin, Lorene::Et_magnetisation, Lorene::Et_rot_bifluid, and Lorene::Et_rot_mag.
Definition at line 554 of file etoile_global.C.
References mass_b(), p_mass_g, and relativistic.
void Lorene::Etoile::operator= | ( | const Etoile & | et | ) |
Assignment to another Etoile
.
Definition at line 433 of file etoile.C.
References a_car, beta_auto, d_logn_auto_div, del_deriv(), ener, ener_euler, ent, eos, gam_euler, k_div, logn_auto, logn_auto_div, logn_auto_regu, mp, nbar, nnn, nzet, press, relativistic, s_euler, shift, u_euler, and unsurc2.
|
protectedvirtual |
Operator >> (virtual function called by the operator <<).
Reimplemented in Lorene::Etoile_rot, Lorene::Etoile_bin, Lorene::Et_magnetisation, Lorene::Et_rot_bifluid, Lorene::Et_bin_nsbh, Lorene::Et_rot_mag, and Lorene::Et_rot_diff.
Definition at line 514 of file etoile.C.
References a_car, ener, ent, eos, k_div, mass_b(), mass_g(), nbar, nnn, nzet, press, ray_eq(), ray_eq_pi(), ray_eq_pis2(), ray_pole(), and relativistic.
double Lorene::Etoile::ray_eq | ( | ) | const |
Coordinate radius at , [r_unit].
Definition at line 123 of file etoile_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::Etoile::ray_eq | ( | int | kk | ) | const |
Coordinate radius at , [r_unit].
Definition at line 443 of file etoile_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(), and mp.
double Lorene::Etoile::ray_eq_3pis2 | ( | ) | const |
Coordinate radius at , [r_unit].
Definition at line 338 of file etoile_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::Etoile::ray_eq_pi | ( | ) | const |
Coordinate radius at , [r_unit].
Definition at line 259 of file etoile_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::Etoile::ray_eq_pis2 | ( | ) | const |
Coordinate radius at , [r_unit].
Definition at line 172 of file etoile_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::Etoile::ray_pole | ( | ) | const |
Coordinate radius at [r_unit].
Definition at line 418 of file etoile_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::Etoile_rot, Lorene::Etoile_bin, Lorene::Et_magnetisation, Lorene::Et_rot_bifluid, Lorene::Et_bin_nsbh, Lorene::Et_rot_mag, Lorene::Et_rot_diff, and Lorene::Et_bin_bhns_extr.
Definition at line 486 of file etoile.C.
References beta_auto, d_logn_auto_div, ent, eos, Lorene::fwrite_be(), k_div, logn_auto, logn_auto_div, nzet, relativistic, Lorene::Eos::sauve(), and Lorene::Tenseur::sauve().
|
protectedvirtual |
Sets to 0x0
all the pointers on derived quantities.
Reimplemented in Lorene::Etoile_rot, Lorene::Etoile_bin, Lorene::Et_rot_bifluid, and Lorene::Et_rot_mag.
Definition at line 399 of file etoile.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::Etoile::set_enthalpy | ( | const Cmp & | ent_i | ) |
Assignment of the enthalpy field.
Definition at line 468 of file etoile.C.
References del_deriv(), ent, and equation_of_state().
|
inline |
const Tbl & Lorene::Etoile::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 104 of file etoile_global.C.
|
friend |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
mutableprotected |
|
mutableprotected |
|
mutableprotected |
|
mutableprotected |
|
mutableprotected |
|
mutableprotected |
|
mutableprotected |
|
mutableprotected |
|
mutableprotected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |