LORENE
Lorene::Eos_multi_poly Class Reference

Base class for a multiple polytropic equation of state. More...

#include <eos_multi_poly.h>

Inheritance diagram for Lorene::Eos_multi_poly:
Lorene::Eos

Public Member Functions

 Eos_multi_poly (int npoly, double *gamma_i, double kappa0_i, double logP1_i, double *logRho_i, double *decInc_i)
 Standard constructor (sets m0 to 1).
 
 Eos_multi_poly (const Eos_multi_poly &)
 Copy constructor.
 
virtual ~Eos_multi_poly ()
 Destructor.
 
void operator= (const Eos_multi_poly &)
 Assignment to another Eos_multi_poly.
 
virtual bool operator== (const Eos &) const
 Read/write kappa.
 
virtual bool operator!= (const Eos &) const
 Comparison operator (difference)
 
virtual int identify () const
 Returns a number to identify the sub-classe of Eos the object belongs to.
 
const intget_npeos () const
 Returns the number of polytropes npeos.
 
const doubleget_gamma (int n) const
 Returns the adiabatic index $\gamma$.
 
const doubleget_kappa0 () const
 Returns the pressure coefficient for the crust.
 
const doubleget_logP1 () const
 Returns the exponent of the pressure at the fiducial density.
 
const doubleget_logRho (int n) const
 Returns the exponent of fiducial densities.
 
const doubleget_kappa (int n) const
 Returns the pressure coefficient $\kappa$ [unit: $\rho_{\rm nuc} c^2 / n_{\rm nuc}^\gamma$], where $\rho_{\rm nuc} := 1.66\ 10^{17} \ {\rm kg/m}^3$ and $n_{\rm nuc} := 0.1 \ {\rm fm}^{-3}$.
 
const doubleget_nbCrit (int n) const
 Returns the critical number density.
 
const doubleget_entCrit (int n) const
 Returns the critical enthalpy.
 
virtual void sauve (FILE *) const
 Save in a file.
 
virtual double nbar_ent_p (double ent, const Param *par=0x0) const
 Computes the baryon density from the log-enthalpy.
 
virtual double ener_ent_p (double ent, const Param *par=0x0) const
 Computes the total energy density from the log-enthalpy.
 
virtual double press_ent_p (double ent, const Param *par=0x0) const
 Computes the pressure from the log-enthalpy.
 
virtual double der_nbar_ent_p (double ent, const Param *par=0x0) const
 Computes the logarithmic derivative $d\ln n/d\ln H$ from the log-enthalpy.
 
virtual double der_ener_ent_p (double ent, const Param *par=0x0) const
 Computes the logarithmic derivative $d\ln e/d\ln H$ from the log-enthalpy.
 
virtual double der_press_ent_p (double ent, const Param *par=0x0) const
 Computes the logarithmic derivative $d\ln p/d\ln H$ from the log-enthalpy.
 
virtual double der_press_nbar_p (double ent, const Param *par=0x0) const
 Computes the logarithmic derivative $d\ln p/d\ln n$ from the log-enthalpy.
 
const charget_name () const
 Returns the EOS name.
 
void set_name (const char *name_i)
 Sets the EOS name.
 
Cmp nbar_ent (const Cmp &ent, int nzet, int l_min=0, const Param *par=0x0) const
 Computes the baryon density field from the log-enthalpy field and extra parameters.
 
Scalar nbar_ent (const Scalar &ent, int nzet, int l_min=0, const Param *par=0x0) const
 Computes the baryon density field from the log-enthalpy field and extra parameters.
 
Cmp ener_ent (const Cmp &ent, int nzet, int l_min=0, const Param *par=0x0) const
 Computes the total energy density from the log-enthalpy and extra parameters.
 
Scalar ener_ent (const Scalar &ent, int nzet, int l_min=0, const Param *par=0x0) const
 Computes the total energy density from the log-enthalpy and extra parameters.
 
Cmp press_ent (const Cmp &ent, int nzet, int l_min=0, const Param *par=0x0) const
 Computes the pressure from the log-enthalpy and extra parameters.
 
Scalar press_ent (const Scalar &ent, int nzet, int l_min=0, const Param *par=0x0) const
 Computes the pressure from the log-enthalpy and extra parameters.
 
Cmp der_nbar_ent (const Cmp &ent, int nzet, int l_min=0, const Param *par=0x0) const
 Computes the logarithmic derivative $d\ln n/d\ln H$ from the log-enthalpy and extra parameters.
 
Scalar der_nbar_ent (const Scalar &ent, int nzet, int l_min=0, const Param *par=0x0) const
 Computes the logarithmic derivative $d\ln n/d\ln H$ from the log-enthalpy and extra parameters.
 
Cmp der_ener_ent (const Cmp &ent, int nzet, int l_min=0, const Param *par=0x0) const
 Computes the logarithmic derivative $d\ln e/d\ln H$ from the log-enthalpy and extra parameters.
 
Scalar der_ener_ent (const Scalar &ent, int nzet, int l_min=0, const Param *par=0x0) const
 Computes the logarithmic derivative $d\ln e/d\ln H$ from the log-enthalpy and extra parameters.
 
Cmp der_press_ent (const Cmp &ent, int nzet, int l_min=0, const Param *par=0x0) const
 Computes the logarithmic derivative $d\ln p/d\ln H$ from the log-enthalpy and extra parameters.
 
Scalar der_press_ent (const Scalar &ent, int nzet, int l_min=0, const Param *par=0x0) const
 Computes the logarithmic derivative $d\ln p/d\ln H$ from the log-enthalpy and extra parameters.
 

Static Public Member Functions

static Eoseos_from_file (FILE *)
 Construction of an EOS from a binary file.
 
static Eoseos_from_file (ifstream &)
 Construction of an EOS from a formatted file.
 

Protected Member Functions

 Eos_multi_poly (FILE *)
 Constructor from a binary file (created by the function sauve(FILE*) ).
 
 Eos_multi_poly (ifstream &)
 Constructor from a formatted file.
 
void set_auxiliary ()
 Computes the auxiliary quantities.
 
virtual ostreamoperator>> (ostream &) const
 Operator >>
 
void calcule (const Cmp &thermo, int nzet, int l_min, double(Eos::*fait)(double, const Param *) const, const Param *par, Cmp &resu) const
 General computational method for Cmp 's.
 
void calcule (const Scalar &thermo, int nzet, int l_min, double(Eos::*fait)(double, const Param *) const, const Param *par, Scalar &resu) const
 General computational method for Scalar 's.
 

Protected Attributes

int npeos
 Number of polytropic equations of state.
 
doublegamma
 Array (size: npeos) of adiabatic index $\gamma$.
 
double kappa0
 Pressure coefficient for the crust [unit: $({\rm g/cm^3})^{1-\gamma_0}$].
 
double logP1
 Exponent of the pressure at the fiducial density $\rho_1$.
 
doublelogRho
 Array (size: npeos - 1) of the exponent of fiducial densities.
 
doublekappa
 Array (size: npeos) of pressure coefficient $\kappa$ [unit: $\rho_{\rm nuc} c^2 / n_{\rm nuc}^\gamma$], where $\rho_{\rm nuc} := 1.66\ 10^{17} \ {\rm kg/m}^3$ and $n_{\rm nuc} := 0.1 \ {\rm fm}^{-3}$.
 
doublenbCrit
 Array (size npeos - 1) of the number density at which the polytropic EOS changes its index and constant.
 
doubleentCrit
 Array (size npeos - 1) of the critical enthalpy at which the polytropic EOS changes its index and constant.
 
doubledecInc
 Array (size npeos - 1) of the percentage which detemines the terminating enthalpy for lower density and the starting enthalpy for higher density.
 
double m0
 Individual particule mass $m0$ [unit: $m_B = 1.66\ 10^{-27} \ {\rm kg}$].
 
doublemu0
 Array (size: npeos) of the relativistic chemical potential at zero pressure [unit: $m_B c^2$, with $m_B = 1.66\ 10^{-27} \ {\rm kg}$].
 
char name [100]
 EOS name.
 

Friends

EosEos::eos_from_file (FILE *)
 The construction functions from a file.
 
EosEos::eos_from_file (ifstream &)
 

Detailed Description

Base class for a multiple polytropic equation of state.

This equation of state mimics some realistic, tabulated EOSs. ()

Definition at line 84 of file eos_multi_poly.h.

Constructor & Destructor Documentation

◆ Eos_multi_poly() [1/4]

Lorene::Eos_multi_poly::Eos_multi_poly ( int  npoly,
double gamma_i,
double  kappa0_i,
double  logP1_i,
double logRho_i,
double decInc_i 
)

Standard constructor (sets m0 to 1).

The individual particle mass $m0$ is set to the mean baryon mass $m_B = 1.66\ 10^{-27} \ {\rm kg}$.

Parameters
npolynumber of polytropes
gamma_iarray of adiabatic index $\gamma$
kappa0_ipressure coefficient for the crust
logP1_iexponent of the pressure at the fiducial density
logRho_iarray of the exponent of fiducial densities
decInc_iarray of percentage

Definition at line 94 of file eos_multi_poly.C.

References decInc, gamma, logRho, npeos, and set_auxiliary().

◆ Eos_multi_poly() [2/4]

Lorene::Eos_multi_poly::Eos_multi_poly ( const Eos_multi_poly eosmp)

Copy constructor.

Definition at line 126 of file eos_multi_poly.C.

References decInc, entCrit, gamma, kappa, logRho, mu0, nbCrit, and npeos.

◆ Eos_multi_poly() [3/4]

Lorene::Eos_multi_poly::Eos_multi_poly ( FILE fich)
protected

Constructor from a binary file (created by the function sauve(FILE*) ).

This constructor is protected because any EOS construction from a binary file must be done via the function Eos::eos_from_file(FILE*) .

Definition at line 176 of file eos_multi_poly.C.

References decInc, Lorene::fread_be(), gamma, kappa0, logP1, logRho, m0, npeos, and set_auxiliary().

◆ Eos_multi_poly() [4/4]

Lorene::Eos_multi_poly::Eos_multi_poly ( ifstream fich)
protected

Constructor from a formatted file.

This constructor is protected because any EOS construction from a formatted file must be done via the function Eos::eos_from_file(ifstream&) .

Definition at line 208 of file eos_multi_poly.C.

References decInc, gamma, kappa0, logP1, logRho, m0, npeos, and set_auxiliary().

◆ ~Eos_multi_poly()

Lorene::Eos_multi_poly::~Eos_multi_poly ( )
virtual

Destructor.

Definition at line 242 of file eos_multi_poly.C.

References decInc, entCrit, gamma, kappa, logRho, mu0, and nbCrit.

Member Function Documentation

◆ calcule() [1/2]

void Lorene::Eos::calcule ( const Cmp thermo,
int  nzet,
int  l_min,
double(Eos::*)(double, const Param *) const  fait,
const Param par,
Cmp resu 
) const
protectedinherited

General computational method for Cmp 's.

Parameters
thermo[input] thermodynamical quantity (for instance the enthalpy field)from which the thermodynamical quantity resu is to be computed.
nzet[input] number of domains where resu is to be computed.
l_min[input] index of the innermost domain is which resu is to be computed [default value: 0]; resu is computed only in domains whose indices are in [l_min,l_min+nzet-1] . In the other domains, it is set to zero.
fait[input] pointer on the member function of class Eos which performs the pointwise calculation.
parpossible extra parameters of the EOS
resu[output] result of the computation.

Definition at line 203 of file eos.C.

◆ calcule() [2/2]

void Lorene::Eos::calcule ( const Scalar thermo,
int  nzet,
int  l_min,
double(Eos::*)(double, const Param *) const  fait,
const Param par,
Scalar resu 
) const
protectedinherited

General computational method for Scalar 's.

Parameters
thermo[input] thermodynamical quantity (for instance the enthalpy field)from which the thermodynamical quantity resu is to be computed.
nzet[input] number of domains where resu is to be computed.
l_min[input] index of the innermost domain is which resu is to be computed [default value: 0]; resu is computed only in domains whose indices are in [l_min,l_min+nzet-1] . In the other domains, it is set to zero.
fait[input] pointer on the member function of class Eos which performs the pointwise calculation.
parpossible extra parameters of the EOS
resu[output] result of the computation.

Definition at line 268 of file eos.C.

◆ der_ener_ent() [1/2]

Cmp Lorene::Eos::der_ener_ent ( const Cmp ent,
int  nzet,
int  l_min = 0,
const Param par = 0x0 
) const
inherited

Computes the logarithmic derivative $d\ln e/d\ln H$ from the log-enthalpy and extra parameters.

Parameters
ent[input, unit: $c^2$] log-enthalpy H defined by $H = c^2 \ln\left( {e+p \over m_B c^2 n} \right) $, where e is the (total) energy density, p the pressure, n the baryon density, and $m_B$ the baryon mass
nzetnumber of domains where the derivative dln(e)/dln(H) is to be computed.
l_minindex of the innermost domain is which the coefficient dln(n)/dln(H) is to be computed [default value: 0]; the derivative dln(e)/dln(H) is computed only in domains whose indices are in [l_min,l_min+nzet-1] . In the other domains, it is set to zero.
parpossible extra parameters of the EOS
Returns
dln(e)/dln(H)

Definition at line 430 of file eos.C.

References Lorene::Cmp::get_mp().

◆ der_ener_ent() [2/2]

Scalar Lorene::Eos::der_ener_ent ( const Scalar ent,
int  nzet,
int  l_min = 0,
const Param par = 0x0 
) const
inherited

Computes the logarithmic derivative $d\ln e/d\ln H$ from the log-enthalpy and extra parameters.

Parameters
ent[input, unit: $c^2$] log-enthalpy H defined by $H = c^2 \ln\left( {e+p \over m_B c^2 n} \right) $, where e is the (total) energy density, p the pressure, n the baryon density, and $m_B$ the baryon mass
nzetnumber of domains where the derivative dln(e)/dln(H) is to be computed.
l_minindex of the innermost domain is which the coefficient dln(n)/dln(H) is to be computed [default value: 0]; the derivative dln(e)/dln(H) is computed only in domains whose indices are in [l_min,l_min+nzet-1] . In the other domains, it is set to zero.
parpossible extra parameters of the EOS
Returns
dln(e)/dln(H)

Definition at line 440 of file eos.C.

References Lorene::Tensor::get_mp().

◆ der_ener_ent_p()

double Lorene::Eos_multi_poly::der_ener_ent_p ( double  ent,
const Param par = 0x0 
) const
virtual

Computes the logarithmic derivative $d\ln e/d\ln H$ from the log-enthalpy.

Parameters
ent[input, unit: $c^2$] log-enthalpy H
parpossible extra parameters of the EOS
Returns
dln(e)/dln(H)

Implements Lorene::Eos.

Definition at line 857 of file eos_multi_poly.C.

References decInc, der_press_ent_p(), der_press_nbar_p(), entCrit, Lorene::exp(), gamma, kappa, Lorene::log10(), m0, mu0, npeos, and Lorene::pow().

◆ der_nbar_ent() [1/2]

Cmp Lorene::Eos::der_nbar_ent ( const Cmp ent,
int  nzet,
int  l_min = 0,
const Param par = 0x0 
) const
inherited

Computes the logarithmic derivative $d\ln n/d\ln H$ from the log-enthalpy and extra parameters.

Parameters
ent[input, unit: $c^2$] log-enthalpy H defined by $H = c^2 \ln\left( {e+p \over m_B c^2 n} \right) $, where e is the (total) energy density, p the pressure, n the baryon density, and $m_B$ the baryon mass
nzetnumber of domains where the derivative dln(n)/dln(H) is to be computed.
l_minindex of the innermost domain is which the coefficient dln(n)/dln(H) is to be computed [default value: 0]; the derivative dln(n)/dln(H) is computed only in domains whose indices are in [l_min,l_min+nzet-1] . In the other domains, it is set to zero.
parpossible extra parameters of the EOS
Returns
dln(n)/dln(H)

Definition at line 407 of file eos.C.

References Lorene::Cmp::get_mp().

◆ der_nbar_ent() [2/2]

Scalar Lorene::Eos::der_nbar_ent ( const Scalar ent,
int  nzet,
int  l_min = 0,
const Param par = 0x0 
) const
inherited

Computes the logarithmic derivative $d\ln n/d\ln H$ from the log-enthalpy and extra parameters.

Parameters
ent[input, unit: $c^2$] log-enthalpy H defined by $H = c^2 \ln\left( {e+p \over m_B c^2 n} \right) $, where e is the (total) energy density, p the pressure, n the baryon density, and $m_B$ the baryon mass
nzetnumber of domains where the derivative dln(n)/dln(H) is to be computed.
l_minindex of the innermost domain is which the coefficient dln(n)/dln(H) is to be computed [default value: 0]; the derivative dln(n)/dln(H) is computed only in domains whose indices are in [l_min,l_min+nzet-1] . In the other domains, it is set to zero.
parpossible extra parameters of the EOS
Returns
dln(n)/dln(H)

Definition at line 417 of file eos.C.

References Lorene::Tensor::get_mp().

◆ der_nbar_ent_p()

double Lorene::Eos_multi_poly::der_nbar_ent_p ( double  ent,
const Param par = 0x0 
) const
virtual

Computes the logarithmic derivative $d\ln n/d\ln H$ from the log-enthalpy.

Parameters
ent[input, unit: $c^2$] log-enthalpy H
parpossible extra parameters of the EOS
Returns
dln(n)/dln(H)

Implements Lorene::Eos.

Definition at line 798 of file eos_multi_poly.C.

References decInc, der_press_ent_p(), der_press_nbar_p(), entCrit, Lorene::exp(), gamma, m0, mu0, and npeos.

◆ der_press_ent() [1/2]

Cmp Lorene::Eos::der_press_ent ( const Cmp ent,
int  nzet,
int  l_min = 0,
const Param par = 0x0 
) const
inherited

Computes the logarithmic derivative $d\ln p/d\ln H$ from the log-enthalpy and extra parameters.

Parameters
ent[input, unit: $c^2$] log-enthalpy H defined by $H = c^2 \ln\left( {e+p \over m_B c^2 n} \right) $, where e is the (total) energy density, p the pressure, n the baryon density, and $m_B$ the baryon mass
nzetnumber of domains where the derivative dln(p)/dln(H) is to be computed.
parpossible extra parameters of the EOS
l_minindex of the innermost domain is which the coefficient dln(n)/dln(H) is to be computed [default value: 0]; the derivative dln(p)/dln(H) is computed only in domains whose indices are in [l_min,l_min+nzet-1] . In the other domains, it is set to zero.
Returns
dln(p)/dln(H)

Definition at line 452 of file eos.C.

References Lorene::Cmp::get_mp().

◆ der_press_ent() [2/2]

Scalar Lorene::Eos::der_press_ent ( const Scalar ent,
int  nzet,
int  l_min = 0,
const Param par = 0x0 
) const
inherited

Computes the logarithmic derivative $d\ln p/d\ln H$ from the log-enthalpy and extra parameters.

Parameters
ent[input, unit: $c^2$] log-enthalpy H defined by $H = c^2 \ln\left( {e+p \over m_B c^2 n} \right) $, where e is the (total) energy density, p the pressure, n the baryon density, and $m_B$ the baryon mass
nzetnumber of domains where the derivative dln(p)/dln(H) is to be computed.
parpossible extra parameters of the EOS
l_minindex of the innermost domain is which the coefficient dln(n)/dln(H) is to be computed [default value: 0]; the derivative dln(p)/dln(H) is computed only in domains whose indices are in [l_min,l_min+nzet-1] . In the other domains, it is set to zero.
Returns
dln(p)/dln(H)

Definition at line 462 of file eos.C.

References Lorene::Tensor::get_mp().

◆ der_press_ent_p()

double Lorene::Eos_multi_poly::der_press_ent_p ( double  ent,
const Param par = 0x0 
) const
virtual

Computes the logarithmic derivative $d\ln p/d\ln H$ from the log-enthalpy.

Parameters
ent[input, unit: $c^2$] log-enthalpy H
parpossible extra parameters of the EOS
Returns
dln(p)/dln(H)

Implements Lorene::Eos.

Definition at line 990 of file eos_multi_poly.C.

References decInc, entCrit, Lorene::exp(), gamma, kappa, Lorene::log10(), m0, mu0, npeos, and Lorene::pow().

◆ der_press_nbar_p()

double Lorene::Eos_multi_poly::der_press_nbar_p ( double  ent,
const Param par = 0x0 
) const
virtual

Computes the logarithmic derivative $d\ln p/d\ln n$ from the log-enthalpy.

Parameters
ent[input, unit: $c^2$] log-enthalpy H
parpossible extra parameters of the EOS
Returns
dln(p)/dln(n)

Definition at line 1088 of file eos_multi_poly.C.

References decInc, entCrit, gamma, Lorene::log10(), and npeos.

◆ ener_ent() [1/2]

Cmp Lorene::Eos::ener_ent ( const Cmp ent,
int  nzet,
int  l_min = 0,
const Param par = 0x0 
) const
inherited

Computes the total energy density from the log-enthalpy and extra parameters.

Parameters
ent[input, unit: $c^2$] log-enthalpy H defined by $H = c^2 \ln\left( {e+p \over m_B c^2 n} \right) $, where e is the (total) energy density, p the pressure, n the baryon density, and $m_B$ the baryon mass
nzetnumber of domains where the energy density is to be computed.
l_minindex of the innermost domain is which the energy density is to be computed [default value: 0]; the energy density is computed only in domains whose indices are in [l_min,l_min+nzet-1] . In the other domains, it is set to zero.
parpossible extra parameters of the EOS
Returns
energy density [unit: $\rho_{\rm nuc} c^2$], where $\rho_{\rm nuc} := 1.66\ 10^{17} \ {\rm kg/m}^3$

Definition at line 363 of file eos.C.

References Lorene::Cmp::get_mp().

◆ ener_ent() [2/2]

Scalar Lorene::Eos::ener_ent ( const Scalar ent,
int  nzet,
int  l_min = 0,
const Param par = 0x0 
) const
inherited

Computes the total energy density from the log-enthalpy and extra parameters.

Parameters
ent[input, unit: $c^2$] log-enthalpy H defined by $H = c^2 \ln\left( {e+p \over m_B c^2 n} \right) $, where e is the (total) energy density, p the pressure, n the baryon density, and $m_B$ the baryon mass
nzetnumber of domains where the energy density is to be computed.
l_minindex of the innermost domain is which the energy density is to be computed [default value: 0]; the energy density is computed only in domains whose indices are in [l_min,l_min+nzet-1] . In the other domains, it is set to zero.
parpossible extra parameters of the EOS
Returns
energy density [unit: $\rho_{\rm nuc} c^2$], where $\rho_{\rm nuc} := 1.66\ 10^{17} \ {\rm kg/m}^3$

Definition at line 373 of file eos.C.

References Lorene::Tensor::get_mp().

◆ ener_ent_p()

double Lorene::Eos_multi_poly::ener_ent_p ( double  ent,
const Param par = 0x0 
) const
virtual

Computes the total energy density from the log-enthalpy.

Parameters
ent[input, unit: $c^2$] log-enthalpy H
parpossible extra parameters of the EOS
Returns
energy density e [unit: $\rho_{\rm nuc} c^2$], where $\rho_{\rm nuc} := 1.66\ 10^{17} \ {\rm kg/m}^3$

Implements Lorene::Eos.

Definition at line 591 of file eos_multi_poly.C.

References decInc, entCrit, Lorene::exp(), gamma, kappa, Lorene::log10(), m0, mu0, npeos, and Lorene::pow().

◆ eos_from_file() [1/2]

Eos * Lorene::Eos::eos_from_file ( FILE fich)
staticinherited

Construction of an EOS from a binary file.

The file must have been created by the function sauve(FILE*) .

Definition at line 177 of file eos_from_file.C.

References Lorene::fread_be().

◆ eos_from_file() [2/2]

Eos * Lorene::Eos::eos_from_file ( ifstream fich)
staticinherited

Construction of an EOS from a formatted file.

The fist line of the file must start by the EOS number, according to the following conventions:

  • 1 = relativistic polytropic EOS (class Eos_poly ).
  • 2 = Newtonian polytropic EOS (class Eos_poly_newt ).
  • 3 = Relativistic incompressible EOS (class Eos_incomp ).
  • 4 = Newtonian incompressible EOS (class Eos_incomp_newt ).
  • 5 = Strange matter (MIT Bag model)
  • 6 = Strange matter (MIT Bag model) with crust
  • 10 = SLy4 (Douchin & Haensel 2001)
  • 11 = FPS (Friedman-Pandharipande + Skyrme)
  • 12 = BPAL12 (Bombaci et al. 1995)
  • 13 = AkmalPR (Akmal, Pandharipande & Ravenhall 1998)
  • 14 = BBB2 (Baldo, Bombaci & Burgio 1997)
  • 15 = BalbN1H1 (Balberg 2000)
  • 16 = GlendNH3 (Glendenning 1985, case 3)
  • 17 = Compstar (Tabulated EOS for 2010 CompStar school)
  • 18 = magnetized (tabulated) equation of state
  • 19 = relativistic ideal Fermi gas at zero temperature (class Eos_Fermi)
  • 100 = Multi-domain EOS (class MEos )
  • 110 = Multi-polytropic EOS (class Eos_multi_poly )
  • 120 = Fitted SLy4 (Shibata 2004)
  • 121 = Fitted FPS (Shibata 2004)
  • 122 = Fitted AkmalPR (Taniguchi 2005)

The second line in the file should contain a name given by the user to the EOS. The following lines should contain the EOS parameters (one parameter per line), in the same order than in the class declaration.

Definition at line 314 of file eos_from_file.C.

◆ get_entCrit()

const double & Lorene::Eos_multi_poly::get_entCrit ( int  n) const
inline

Returns the critical enthalpy.

Definition at line 252 of file eos_multi_poly.h.

References entCrit, and npeos.

◆ get_gamma()

const double & Lorene::Eos_multi_poly::get_gamma ( int  n) const
inline

Returns the adiabatic index $\gamma$.

Definition at line 217 of file eos_multi_poly.h.

References gamma, and npeos.

◆ get_kappa()

const double & Lorene::Eos_multi_poly::get_kappa ( int  n) const
inline

Returns the pressure coefficient $\kappa$ [unit: $\rho_{\rm nuc} c^2 / n_{\rm nuc}^\gamma$], where $\rho_{\rm nuc} := 1.66\ 10^{17} \ {\rm kg/m}^3$ and $n_{\rm nuc} := 0.1 \ {\rm fm}^{-3}$.

Definition at line 240 of file eos_multi_poly.h.

References kappa, and npeos.

◆ get_kappa0()

const double & Lorene::Eos_multi_poly::get_kappa0 ( ) const
inline

Returns the pressure coefficient for the crust.

Definition at line 223 of file eos_multi_poly.h.

References kappa0.

◆ get_logP1()

const double & Lorene::Eos_multi_poly::get_logP1 ( ) const
inline

Returns the exponent of the pressure at the fiducial density.

Definition at line 226 of file eos_multi_poly.h.

References logP1.

◆ get_logRho()

const double & Lorene::Eos_multi_poly::get_logRho ( int  n) const
inline

Returns the exponent of fiducial densities.

Definition at line 229 of file eos_multi_poly.h.

References logRho, and npeos.

◆ get_name()

const char * Lorene::Eos::get_name ( ) const
inherited

Returns the EOS name.

Definition at line 169 of file eos.C.

References Lorene::Eos::name.

◆ get_nbCrit()

const double & Lorene::Eos_multi_poly::get_nbCrit ( int  n) const
inline

Returns the critical number density.

Definition at line 246 of file eos_multi_poly.h.

References nbCrit, and npeos.

◆ get_npeos()

const int & Lorene::Eos_multi_poly::get_npeos ( ) const
inline

Returns the number of polytropes npeos.

Definition at line 214 of file eos_multi_poly.h.

References npeos.

◆ identify()

int Lorene::Eos_multi_poly::identify ( ) const
virtual

Returns a number to identify the sub-classe of Eos the object belongs to.

Implements Lorene::Eos.

Definition at line 165 of file eos_from_file.C.

◆ nbar_ent() [1/2]

Cmp Lorene::Eos::nbar_ent ( const Cmp ent,
int  nzet,
int  l_min = 0,
const Param par = 0x0 
) const
inherited

Computes the baryon density field from the log-enthalpy field and extra parameters.

Parameters
ent[input, unit: $c^2$] log-enthalpy H defined by $H = c^2 \ln\left( {e+p \over m_B c^2 n} \right) $, where e is the (total) energy density, p the pressure, n the baryon density, and $m_B$ the baryon mass
nzetnumber of domains where the baryon density is to be computed.
l_minindex of the innermost domain is which the baryon density is to be computed [default value: 0]; the baryon density is computed only in domains whose indices are in [l_min,l_min+nzet-1] . In the other domains, it is set to zero.
parpossible extra parameters of the EOS
Returns
baryon density [unit: $n_{\rm nuc} := 0.1 \ {\rm fm}^{-3}$]

Definition at line 338 of file eos.C.

References Lorene::Cmp::get_mp().

◆ nbar_ent() [2/2]

Scalar Lorene::Eos::nbar_ent ( const Scalar ent,
int  nzet,
int  l_min = 0,
const Param par = 0x0 
) const
inherited

Computes the baryon density field from the log-enthalpy field and extra parameters.

Parameters
ent[input, unit: $c^2$] log-enthalpy H defined by $H = c^2 \ln\left( {e+p \over m_B c^2 n} \right) $, where e is the (total) energy density, p the pressure, n the baryon density, and $m_B$ the baryon mass
nzetnumber of domains where the baryon density is to be computed.
l_minindex of the innermost domain is which the baryon density is to be computed [default value: 0]; the baryon density is computed only in domains whose indices are in [l_min,l_min+nzet-1] . In the other domains, it is set to zero.
parpossible extra parameters of the EOS
Returns
baryon density [unit: $n_{\rm nuc} := 0.1 \ {\rm fm}^{-3}$]

Definition at line 348 of file eos.C.

References Lorene::Tensor::get_mp().

◆ nbar_ent_p()

double Lorene::Eos_multi_poly::nbar_ent_p ( double  ent,
const Param par = 0x0 
) const
virtual

Computes the baryon density from the log-enthalpy.

Parameters
ent[input, unit: $c^2$] log-enthalpy H
parpossible extra parameters of the EOS
Returns
baryon density n [unit: $n_{\rm nuc} := 0.1 \ {\rm fm}^{-3}$]

Implements Lorene::Eos.

Definition at line 493 of file eos_multi_poly.C.

References decInc, entCrit, Lorene::exp(), gamma, kappa, Lorene::log10(), m0, mu0, npeos, and Lorene::pow().

◆ operator!=()

bool Lorene::Eos_multi_poly::operator!= ( const Eos eos_i) const
virtual

Comparison operator (difference)

Implements Lorene::Eos.

Definition at line 385 of file eos_multi_poly.C.

References operator==().

◆ operator=()

void Lorene::Eos_multi_poly::operator= ( const Eos_multi_poly )

Assignment to another Eos_multi_poly.

Definition at line 258 of file eos_multi_poly.C.

◆ operator==()

bool Lorene::Eos_multi_poly::operator== ( const Eos eos_i) const
virtual

Read/write kappa.

Comparison operator (egality)

Implements Lorene::Eos.

Definition at line 341 of file eos_multi_poly.C.

References gamma, get_gamma(), get_kappa(), get_npeos(), identify(), Lorene::Eos::identify(), kappa, and npeos.

◆ operator>>()

ostream & Lorene::Eos_multi_poly::operator>> ( ostream ost) const
protectedvirtual

Operator >>

Implements Lorene::Eos.

Definition at line 419 of file eos_multi_poly.C.

References entCrit, gamma, kappa, logP1, logRho, mu0, nbCrit, npeos, and Lorene::pow().

◆ press_ent() [1/2]

Cmp Lorene::Eos::press_ent ( const Cmp ent,
int  nzet,
int  l_min = 0,
const Param par = 0x0 
) const
inherited

Computes the pressure from the log-enthalpy and extra parameters.

Parameters
ent[input, unit: $c^2$] log-enthalpy H defined by $H = c^2 \ln\left( {e+p \over m_B c^2 n} \right) $, where e is the (total) energy density, p the pressure, n the baryon density, and $m_B$ the baryon mass
nzetnumber of domains where the pressure is to be computed.
l_minindex of the innermost domain is which the pressure is to be computed [default value: 0]; the pressure is computed only in domains whose indices are in [l_min,l_min+nzet-1] . In the other domains, it is set to zero.
parpossible extra parameters of the EOS
Returns
pressure [unit: $\rho_{\rm nuc} c^2$], where $\rho_{\rm nuc} := 1.66\ 10^{17} \ {\rm kg/m}^3$

Definition at line 385 of file eos.C.

References Lorene::Cmp::get_mp().

◆ press_ent() [2/2]

Scalar Lorene::Eos::press_ent ( const Scalar ent,
int  nzet,
int  l_min = 0,
const Param par = 0x0 
) const
inherited

Computes the pressure from the log-enthalpy and extra parameters.

Parameters
ent[input, unit: $c^2$] log-enthalpy H defined by $H = c^2 \ln\left( {e+p \over m_B c^2 n} \right) $, where e is the (total) energy density, p the pressure, n the baryon density, and $m_B$ the baryon mass
nzetnumber of domains where the pressure is to be computed.
l_minindex of the innermost domain is which the pressure is to be computed [default value: 0]; the pressure is computed only in domains whose indices are in [l_min,l_min+nzet-1] . In the other domains, it is set to zero.
parpossible extra parameters of the EOS
Returns
pressure [unit: $\rho_{\rm nuc} c^2$], where $\rho_{\rm nuc} := 1.66\ 10^{17} \ {\rm kg/m}^3$

Definition at line 395 of file eos.C.

References Lorene::Tensor::get_mp().

◆ press_ent_p()

double Lorene::Eos_multi_poly::press_ent_p ( double  ent,
const Param par = 0x0 
) const
virtual

Computes the pressure from the log-enthalpy.

Parameters
ent[input, unit: $c^2$] log-enthalpy H
parpossible extra parameters of the EOS
Returns
pressure p [unit: $\rho_{\rm nuc} c^2$], where $\rho_{\rm nuc} := 1.66\ 10^{17} \ {\rm kg/m}^3$

Implements Lorene::Eos.

Definition at line 700 of file eos_multi_poly.C.

References decInc, entCrit, Lorene::exp(), gamma, kappa, Lorene::log10(), m0, mu0, npeos, and Lorene::pow().

◆ sauve()

void Lorene::Eos_multi_poly::sauve ( FILE fich) const
virtual

Save in a file.

Reimplemented from Lorene::Eos.

Definition at line 395 of file eos_multi_poly.C.

References decInc, Lorene::fwrite_be(), gamma, kappa0, logP1, logRho, npeos, and Lorene::Eos::sauve().

◆ set_auxiliary()

void Lorene::Eos_multi_poly::set_auxiliary ( )
protected

Computes the auxiliary quantities.

Definition at line 270 of file eos_multi_poly.C.

References entCrit, gamma, kappa, kappa0, Lorene::log(), logP1, logRho, m0, mu0, nbCrit, npeos, and Lorene::pow().

◆ set_name()

void Lorene::Eos::set_name ( const char name_i)
inherited

Sets the EOS name.

Definition at line 163 of file eos.C.

References Lorene::Eos::name.

Friends And Related Symbol Documentation

◆ Eos::eos_from_file

Eos * Eos::eos_from_file ( FILE )
friend

The construction functions from a file.

Member Data Documentation

◆ decInc

double* Lorene::Eos_multi_poly::decInc
protected

Array (size npeos - 1) of the percentage which detemines the terminating enthalpy for lower density and the starting enthalpy for higher density.

Definition at line 129 of file eos_multi_poly.h.

◆ entCrit

double* Lorene::Eos_multi_poly::entCrit
protected

Array (size npeos - 1) of the critical enthalpy at which the polytropic EOS changes its index and constant.

Definition at line 123 of file eos_multi_poly.h.

◆ gamma

double* Lorene::Eos_multi_poly::gamma
protected

Array (size: npeos) of adiabatic index $\gamma$.

Definition at line 94 of file eos_multi_poly.h.

◆ kappa

double* Lorene::Eos_multi_poly::kappa
protected

Array (size: npeos) of pressure coefficient $\kappa$ [unit: $\rho_{\rm nuc} c^2 / n_{\rm nuc}^\gamma$], where $\rho_{\rm nuc} := 1.66\ 10^{17} \ {\rm kg/m}^3$ and $n_{\rm nuc} := 0.1 \ {\rm fm}^{-3}$.

Definition at line 113 of file eos_multi_poly.h.

◆ kappa0

double Lorene::Eos_multi_poly::kappa0
protected

Pressure coefficient for the crust [unit: $({\rm g/cm^3})^{1-\gamma_0}$].

Definition at line 99 of file eos_multi_poly.h.

◆ logP1

double Lorene::Eos_multi_poly::logP1
protected

Exponent of the pressure at the fiducial density $\rho_1$.

Definition at line 102 of file eos_multi_poly.h.

◆ logRho

double* Lorene::Eos_multi_poly::logRho
protected

Array (size: npeos - 1) of the exponent of fiducial densities.

Definition at line 105 of file eos_multi_poly.h.

◆ m0

double Lorene::Eos_multi_poly::m0
protected

Individual particule mass $m0$ [unit: $m_B = 1.66\ 10^{-27} \ {\rm kg}$].

Definition at line 134 of file eos_multi_poly.h.

◆ mu0

double* Lorene::Eos_multi_poly::mu0
protected

Array (size: npeos) of the relativistic chemical potential at zero pressure [unit: $m_B c^2$, with $m_B = 1.66\ 10^{-27} \ {\rm kg}$].

(The value for the EOS which covers the lowest density: 1)

Definition at line 142 of file eos_multi_poly.h.

◆ name

char Lorene::Eos::name[100]
protectedinherited

EOS name.

Definition at line 196 of file eos.h.

◆ nbCrit

double* Lorene::Eos_multi_poly::nbCrit
protected

Array (size npeos - 1) of the number density at which the polytropic EOS changes its index and constant.

Definition at line 118 of file eos_multi_poly.h.

◆ npeos

int Lorene::Eos_multi_poly::npeos
protected

Number of polytropic equations of state.

Definition at line 91 of file eos_multi_poly.h.


The documentation for this class was generated from the following files: