Etoile (Map &mp_i, int nzet_i, bool relat, const Eos &eos_i)
Standard constructor.
Etoile (const Etoile &)
Copy constructor.
Etoile (Map &mp_i, const Eos &eos_i, FILE *fich)
Constructor from a file (see sauve(FILE*) ).
virtual ~Etoile ()
Destructor.
void operator= (const Etoile &)
Assignment to another Etoile .
Map & set_mp ()
Read/write of the mapping.
void set_enthalpy (const Cmp &)
Assignment of the enthalpy field.
virtual void equation_of_state ()
Computes the proper baryon and energy density, as well as pressure from the enthalpy.
virtual void hydro_euler ()
Computes the hydrodynamical quantities relative to the Eulerian observer from those in the fluid frame (nbar , ener and press ).
virtual void equilibrium_spher (double ent_c, double precis=1.e-14, const Tbl *ent_limit=0x0)
Computes a spherical static configuration.
void equil_spher_regular (double ent_c, double precis=1.e-14)
Computes a spherical static configuration.
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.
const Map & get_mp () const
Returns the mapping.
int get_nzet () const
Returns the number of domains occupied by the star.
bool is_relativistic () const
Returns true for a relativistic star, false for a Newtonian one.
const Eos & get_eos () const
Returns the equation of state.
const Tenseur & get_ent () const
Returns the enthalpy field.
const Tenseur & get_nbar () const
Returns the proper baryon density.
const Tenseur & get_ener () const
Returns the proper total energy density.
const Tenseur & get_press () const
Returns the fluid pressure.
const Tenseur & get_ener_euler () const
Returns the total energy density with respect to the Eulerian observer.
const Tenseur & get_s_euler () const
Returns the trace of the stress tensor in the Eulerian frame.
const Tenseur & get_gam_euler () const
Returns the Lorentz factor between the fluid and Eulerian observers.
const Tenseur & get_u_euler () const
Returns the fluid 3-velocity with respect to the Eulerian observer.
const Tenseur & get_logn_auto () const
Returns the logarithm of the part of the lapse N generated principaly by the star.
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.
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.
const Tenseur & get_d_logn_auto_div () const
Returns the gradient of logn_auto_div.
const Tenseur & get_beta_auto () const
Returns the logarithm of the part of the product AN generated principaly by the star.
const Tenseur & get_nnn () const
Returns the total lapse function N .
const Tenseur & get_shift () const
Returns the total shift vector .
const Tenseur & get_a_car () const
Returns the total conformal factor .
virtual void sauve (FILE *) const
Save in a file.
double ray_eq () const
Coordinate radius at , [r_unit].
double ray_eq_pis2 () const
Coordinate radius at , [r_unit].
double ray_eq_pi () const
Coordinate radius at , [r_unit].
double ray_eq_3pis2 () const
Coordinate radius at , [r_unit].
double ray_pole () const
Coordinate radius at [r_unit].
double ray_eq (int kk) const
Coordinate radius at , [r_unit].
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 on the surface at the collocation points in .
virtual double mass_b () const
Baryon mass.
virtual double mass_g () const
Gravitational mass.
Map & mp
Mapping associated with the star.
int nzet
Number of domains of *mp occupied by the star.
bool relativistic
Indicator of relativity: true for a relativistic star, false for a Newtonian one.
double unsurc2
: unsurc2=1 for a relativistic star, 0 for a Newtonian one.
int k_div
Index of regularity of the gravitational potential logn_auto .
const Eos & eos
Equation of state of the stellar matter.
Tenseur ent
Log-enthalpy (relativistic case) or specific enthalpy (Newtonian case).
Tenseur nbar
Baryon density in the fluid frame.
Tenseur ener
Total energy density in the fluid frame.
Tenseur press
Fluid pressure.
Tenseur ener_euler
Total energy density in the Eulerian frame.
Tenseur s_euler
Trace of the stress tensor in the Eulerian frame.
Tenseur gam_euler
Lorentz factor between the fluid and Eulerian observers.
Tenseur u_euler
Fluid 3-velocity with respect to the Eulerian observer.
Tenseur logn_auto
Total of the logarithm of the part of the lapse N
generated principaly by the star.
Tenseur logn_auto_regu
Regular part of the logarithm of the part of the lapse N
generated principaly by the star.
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.
Tenseur d_logn_auto_div
Gradient of logn_auto_div (if k_div!=0 ).
Tenseur beta_auto
Logarithm of the part of the product AN generated principaly by by the star.
Tenseur nnn
Total lapse function.
Tenseur shift
Total shift vector.
Tenseur a_car
Total conformal factor .
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 on the surface at the collocation points in .
double * p_mass_b
Baryon mass.
double * p_mass_g
Gravitational mass.
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.
Definition at line 424 of file etoile.h .
void Lorene::Etoile::equil_spher_falloff
(
double ent_c ,
double precis = 1.e-14 )
virtual
Computes a spherical static configuration with the outer boundary condition at a finite radius.
Parameters
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 57 of file etoile_eqsph_falloff.C .
References a_car , beta_auto , Lorene::diffrel() , Lorene::Cmp::dsdr() , Lorene::Map_af::dsdr() , ener , ener_euler , ent , equation_of_state() , Lorene::exp() , gam_euler , 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 , and unsurc2 .
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.
Parameters
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 115 of file et_equil_spher_regu.C .
References a_car , beta_auto , d_logn_auto_div , Lorene::diffrel() , Lorene::Map_af::dsdr() , ener , ener_euler , ent , eos , equation_of_state() , Lorene::exp() , gam_euler , Lorene::Mg3d::get_nr() , Lorene::Mg3d::get_nt() , Lorene::Mg3d::get_nzone() , 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() , shift , Lorene::sqrt() , Lorene::Cmp::std_base_scal() , u_euler , and unsurc2 .
void Lorene::Etoile::equilibrium_spher
(
double ent_c ,
double precis = 1.e-14 ,
const Tbl * ent_limit = 0x0 )
virtual
Computes a spherical static configuration.
Parameters
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 87 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() , 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::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 , and unsurc2 .