February 24, 2013

Abstract

This thorn provides a unified EOS (Equation Of State) interface and implements multiple analytic EOS, and also provides table reader and interpolation routines for finite-temperature microphysical EOS available from http://www.stellarcollapse.org. In addition cold and barotropic tabulated EOS are provided for use in initial data thorns. Currently, the implemented analytic EOS are the polytropic EOS, the gamma-law EOS, and a hybrid EOS consisting of a n-piece piecewise-polytrope with a thermal, gamma-law component.

Equations of State (EOS) are crucial for hydrodynamics and hydro codes (as well as other codes needing/providing microphysics) are closely coupled to EOS and call EOS routines many times for each grid point during the calculation of a time update.

EOS_Omni is presently coded for cold and hot EOS, including those based on microphysical models. It does currently assume nuclear statistical equilibrium (NSE) with rest-mass density $\rho $, specific internal energy $\mathit{\epsilon}$ (or temperature $T$), and electron fraction ${Y}_{e}$ being the independent variables. EOS_Omni can be called on arrays or on single grid points.

This thorn uses solar units where $c=G={M}_{\odot}=1$. Temperatures are measured in MeV.

EOS_Omni works via the aliased-function interface, and EOS functions to be used must be declared in interface.ccl. Here is an example interface.ccl entry:

void FUNCTION EOS_Omni_press(CCTK_INT IN eoskey, \

CCTK_INT IN havetemp, \

CCTK_REAL IN rf_precision, \

CCTK_INT IN npoints, \

CCTK_REAL IN ARRAY rho, \

CCTK_REAL INOUT ARRAY eps, \

CCTK_REAL INOUT ARRAY temp, \

CCTK_REAL IN ARRAY ye, \

CCTK_REAL OUT ARRAY press \

CCTK_INT OUT ARRAY keyerr, \

CCTK_INT OUT anyerr)

CCTK_INT IN havetemp, \

CCTK_REAL IN rf_precision, \

CCTK_INT IN npoints, \

CCTK_REAL IN ARRAY rho, \

CCTK_REAL INOUT ARRAY eps, \

CCTK_REAL INOUT ARRAY temp, \

CCTK_REAL IN ARRAY ye, \

CCTK_REAL OUT ARRAY press \

CCTK_INT OUT ARRAY keyerr, \

CCTK_INT OUT anyerr)

Here,

- eoskey is the type of EOS to be used in this call.
- havetemp determines whether the EOS is to be called as a function of $\left(\rho ,\mathit{\epsilon},{Y}_{e}\right)$ (havetemp = 0), or as a function of $\left(\rho ,T,{Y}_{e}\right)$ (havetemp = 1). havetemp = 0 is the method of choice for analytic EOS during evolution, but at the initial data stage one may need to set havetemp = 1 (with $T=0$) to obtain initial values for $\mathit{\epsilon}$. In the case of a finite-temperature EOS (that usually is a function of $\left(\rho ,T,{Y}_{e}\right)$), a call with havetemp = 0 will first lead to the solution of $T(\rho ,\mathit{\epsilon},{Y}_{e}$) via a Newton-Raphson-type iteration (using the supplied value of $T$ as the starting point) and will then calculate the requested dependent variable as a function of $X=X\left(\rho ,T,{Y}_{e}\right)$. Both $X$ and the updated $T$ are returned.
- npoints tells the EOS how many data points are passed in,
- rho,eps,temp,ye,press are obvious,
- rf_precision is a real number telling the root finding algorithm (for finding $T\left(\rho ,\mathit{\epsilon},{Y}_{e}\right)$) at what relative error to terminate the iteration. $1{0}^{-10}$ is a good value.
- keyerr is an array (with $n$ entries for $n$ data points) containing error codes (relevant only for tabulated microphysical EOS),
- anyerr is an integer $>0$ in case any error occured.

The eoskey is obtained by calling EOS_Omni_GetHandle with the name of desired EOS. The currently available EOS are:

- 2D_Polytrope
- Polytropic EOS,
- Ideal_Fluid
- Gamma-Law EOS,
- Hybrid
- Hybrid EOS (N Polys, 1 Gamma-Law), used for stellar core collapse simulations,
- nuc_eos
- Finite-temperature microphysical EOSm
- cold_tabulated
- Cold tabulated EOS with Gamma-Lawm
- barotropic_tabulated
- Tabulated barotropic EOS.

Many hydro codes require a fallback EOS in case something goes wrong. This is also true for the Einstein Toolkit GR hydro code EinsteinEvolve/GRHydro. If you want to use EOS_Omni with GRHydro, you must set parameters for the EOS of your choice and, in addition, the following parameters must be set to sensible values:

eos_omni::poly_gamma

eos_omni::poly_gamma_initial

eos_omni::poly_k

eos_omni::poly_gamma_initial

eos_omni::poly_k

The only non-obvious parameter here is poly_gamma_initial. In most simulations it should not be set. In simulations that are run with a different adiabatic index than what was used to set up the initial data, poly_gamma should be the evolution value, and poly_gamma_initial should be the initial data value. EOS_Omni then rescales poly_k such that the cgs value of poly_k is the same for initial data and evolution. Since the units of poly_k ($\left[K\right]=\frac{N}{{m}^{2}}{\left(\frac{{m}^{3}}{kg}\right)}^{\Gamma}$) depend on poly_gamma this is not fully trivial.

Check param.ccl for parameters for the other EOS.

The poly EOS is a polytropic equation of state, which does not allow for changes in entropy:

$$\begin{array}{rcll}p& =& K{\rho}^{\gamma}& \text{(1)}\text{}\text{}\end{array}$$

where $p$ is the pressure, $\rho $ the density, $K$ the polytropic constant set via poly_k, and $\gamma $ is the adiabatic index set via poly_gamma.

If the internal energy $\mathit{\epsilon}$ is to be “calculated from the temperature” (havetemp = 1), then this is done using the relation

$$\begin{array}{rcll}\mathit{\epsilon}& =& \frac{K}{\gamma -1}{\rho}^{\gamma -1}& \text{(2)}\text{}\text{}\end{array}$$

(which actually ignores the temperature).

Note: This polytropic EOS is also used as fall-back when other EOS fail.

The gl EOS is a gamma-law equation of state, corresponding to an ideal gas:

$$\begin{array}{rcll}p& =& \left(\gamma -1\right)\rho \mathit{\epsilon}& \text{(3)}\text{}\text{}\end{array}$$

where $p$ is the pressure, $\rho $ the density, $\mathit{\epsilon}$ the internal energy, and $\gamma $ is the adiabatic index set via gl_gamma.

At the initial data stage, it may be necessary to set up initial values for $\mathit{\epsilon}$. For this, the gl EOS implements equation (2) just like the poly EOS and the parameters poly_gamma_ini and gl_k must be set for this.

The hybrid EOS was introduced by [1] for use in simplified simulations of stellar collapse to mimic (1) the stiffening of the nuclear EOS at nuclear density and (2) to include thermal pressure in the postbounce phase. It consists of n_pieces polytropes (cold) characterized by (${K}_{0},{K}_{1},\dots ,{K}_{n\text{\_}pieces-1}$, ${\gamma}_{0},{\gamma}_{1},\dots ,{\gamma}_{n\text{\_}pieces-1}$) and a thermal $\gamma -$law component described by ${\gamma}_{th}$. Polytrope $i$ is used below the dividing density hybrid_rho[i] and is smoothly matched to polytrope $i+1$ at the dividing density. Because of this ${K}_{i+1}$ for $i>0$ is completely determined by requiring continuity of ${P}_{cold}\left(\rho \right)$ across hybrid_rho[i], ${K}_{0}$, and ${\gamma}_{0},{\gamma}_{1},\dots ,{\gamma}_{i+1}$.

$$\begin{array}{rcll}{P}_{cold}& =& {K}_{i}{\rho}^{{\gamma}_{i}}\phantom{\rule{2em}{0ex}}\text{for}\rho \text{}\mathtt{\text{hybrid\_rho[i]}}\text{}& \text{(4)}\text{}\text{}\\ {\mathit{\epsilon}}_{cold}& =& {K}_{i}\u2215\left({\gamma}_{i}-1\right){\rho}^{{\gamma}_{i}-1}+{\mathit{\epsilon}}_{i},& \text{(5)}\text{}\text{}\end{array}$$

where ${\mathit{\epsilon}}_{i}$ is determined by requiring continuity of ${\mathit{\epsilon}}_{cold}$. The full functional form of the EOS $P=P\left(\rho ,\mathit{\epsilon}\right)$ includes a thermal component, which takes into account shock heating:

$$\begin{array}{rcll}{P}_{th}& =& \left({\gamma}_{th}-1\right)\phantom{\rule{0.3em}{0ex}}\rho \phantom{\rule{0.3em}{0ex}}\left(\mathit{\epsilon}-{\mathit{\epsilon}}_{cold}\right),& \text{(6)}\text{}\text{}\\ P& =& {P}_{cold}+{P}_{th}& \text{(7)}\text{}\text{}\end{array}$$

The EOS_Omni parameters for the hybrid EOS are the following:

n_pieces | number of piecewise polytropic pieces, |

hybrid_gamma | ${\gamma}_{i}$, $0\le i<\text{}\mathtt{\text{n\_pieces}}$, |

hybrid_gamma_th | ${\gamma}_{th}$, |

hybrid_k0 | ${K}_{0}$, |

hybrid_rho | dividing densities, $0\le i<\text{}\mathtt{\text{n\_pieces}}-1$. |

For example for a simple 2-piece piecewise polytropic EOS mimicking a nuclear equation of state:

hybrid_gamma[0] | ${\gamma}_{0}=1.325$ is an appropriate choice. |

hybrid_gamma[1] | ${\gamma}_{1}=2.5$ is an appropriate choice. |

hybrid_gamma_th | ${\gamma}_{th}$, perhaps $1.5$. |

hybrid_k0 | ${K}_{0}$, $0.4640517$ in solar units for relativistic degenerate e${}^{-}$. |

hybrid_rho[0] | nuclear density, standard is $3.238607\times 1{0}^{-4}$ in solar units. |

Complex microphysical finite-temperature equations of state come usually in tabulated form. EOS_Omni comes
with routines provided as part of the nuc_eos package described in [2] and available at

http://www.stellarcollapse.org. A variety of EOS tables for application in high-density astrophysical
situations (i.e. in stellar collapse or in compact star mergers) are also available from there in HDF5
format.

The parameters controlling the finite-temperature nuclear EOS are the following:

nuceos_read_table | BOOLEAN | Set to yes to read table. |

do_energy_shift | BOOLEAN | Set to yes to subtract the energy shift |

stored in the table to get correctly normalized $\mathit{\epsilon}$. | ||

nuceos_table_name | STRING | Path/Name of the table file. |

Many equations of state for neutron stars are generated under the assumption of zero temperature. This is perfectly appropriate for cold old neutron stars. In simulations of binary mergers, however, shocks will drive the temperature up, adding a thermal pressure component, which can be accounted for approximately with a Gamma-Law: ${P}_{th}=\left({\Gamma}_{th}-1\right)\rho {\mathit{\epsilon}}_{th}$ .

EOS_Omni implements such an equation of state. It reads in an ASCII EOS table (see subdirector tables for an example table for the SLy EOS [3, 4], which was generated according to the prescription in [5, 6]). All EOS parameters are read from the ASCII file, which has the following format:

EoSType = Tabulated

Nrho = 600 NYe = 1 NT = 1

RhoMin = 1e-09 RhoMax = 0.01

HeatCapacityE = 1

GammaTh = 2

Kappa = 1

RhoSpacing = Log

1.57940636422747e-03 1.38773826035349e+00 2.62139412738900e-02

[...]

2.81804006881059e+00 6.89967555695907e-01 8.30537612378975e-01

Nrho = 600 NYe = 1 NT = 1

RhoMin = 1e-09 RhoMax = 0.01

HeatCapacityE = 1

GammaTh = 2

Kappa = 1

RhoSpacing = Log

1.57940636422747e-03 1.38773826035349e+00 2.62139412738900e-02

[...]

2.81804006881059e+00 6.89967555695907e-01 8.30537612378975e-01

The header completely determines the range in baryon rest mass density (in $c=G={M}_{\odot}=1$ units), gives the number of zones (currently only NYe = 1, NT = 1, HeatCapacityE = 1, and RhoSpacing = Log are supported). GammaTh is the ${\Gamma}_{th}$ of the thermal gamma law. The tabulated columns are $\mathit{\epsilon}$, $\Gamma $, ${c}_{s}$ (the speed of sound of the cold component of the EOS). The pressure is obtained via $P=\kappa {\rho}^{\Gamma}$, where $\kappa $ is the Kappa scaling parameter.

Generally, $P=P\left(\rho ,\mathit{\epsilon}\right)$ in this EOS, but note that ${\mathit{\epsilon}}_{th}=\mathit{\epsilon}-{\mathit{\epsilon}}_{cold}$. EOS_Omni uses linear interpolation to first find $P\left(\rho \right)$ and ${\mathit{\epsilon}}_{cold}\left(\rho \right)$ and then computes the thermal component analytically.

coldeos_read_table | BOOLEAN | Set to yes to read table. |

coldeos_use_thermal_gamma_law | BOOLEAN | Set to yes to use the thermal gamma law (default). |

coldeos_table_name | STRING | Path/Name of the table file. |

If you have a parameter file that uses the previous EOS interface in Cactus, you will have to convert it so that it runs with EOS_Omni. The following describes a set of simple rules for this conversion.

- Add EOS_Omni to the thorn list. You can then remove all other EOS_* thorns from the thorn list (or you can leave them in; they are unused).
- Activate EOS_Omni in the parameter file: Add EOS_Omni to your active thorns, and do not activate any other EOS_* thorns.
- Translate all EOS parameters according to table 1.
- All thorns using this EOS interface will have a parameter that determines which EOS to use, typically via a string/keyword parameter specifying an EOS name. Convert these names using table 2.

Old Parameter | Old Value | New Parameter | New Value |

EOS_Polytrope::eos_gamma | EOS_Omni::poly_gamma | ||

EOS_Polytrope::eos_k | EOS_Omni::poly_k | ||

EOS_Polytrope::use_cgs | yes | — | |

EOS_Polytrope::use_cgs | no | unsupported | |

EOS_Polytrope::gamma_ini | EOS_Omni::poly_gamma_ini | ||

EOS_Ideal_Fluid::eos_ideal_fluid_gamma | EOS_Omni::gl_gamma | ||

TODO: complete this table |

EOS | Description | Old Name | New Name |

poly | polytropic | ??? | 2D_Polytrope |

gl | gamma-law | ??? | Ideal_Fluid |

hybrid | hybrid | ??? | Hybrid |

nuc | finite-temperature nuclear | ??? | nuc_eos |

[1] Janka, H.-T., Zwerger, T., & Moenchmeyer, R. 1993, Astron. Astrophys., 268, 360

[2] O’Connor, E., & Ott, C. D. 2010, Class. Quantum Grav., 27, 114103

[3] Douchin, F., & Haensel, P. 2001, Astron. Astrophys., 380, 151

[4] Haensel, P., & Potekhin, A. Y. 2004, Astron. Astrophys., 428, 191

[5] Shibata, M., Taniguchi, K., & Uryū, K. 2005, Phys. Rev. D, 71, 084021

[6] Corvino, G., Rezzolla, L., Bernuzzi, S., De Pietri, R., & Giacomazzo, B. 2010, Class. Quantum Grav., 27, 114104

barotropiceos_gammath | Scope: restricted | REAL |

Description: thermal gamma for barotropic EOS
| ||

Range | Default: 2.0 | |

1.0:* | something
| |

barotropiceos_read_table | Scope: restricted | BOOLEAN |

Description: Read in barotropic EOS table?
| ||

Default: No | ||

barotropiceos_table_name | Scope: restricted | STRING |

Description: table name for barotropic EOS (ASCII)
| ||

Range | Default: blah.asc | |

.+ | Can be anything
| |

barotropiceos_use_thermal_gamma_law | Scope: restricted | BOOLEAN |

Description: use an additional thermal gamma?
| ||

Default: Yes | ||

coldeos_read_table | Scope: restricted | BOOLEAN |

Description: Read in cold EOS table?
| ||

Default: No | ||

coldeos_table_name | Scope: restricted | STRING |

Description: table name for cold EOS (ASCII)
| ||

Range | Default: blah.asc | |

.+ | Can be anything
| |

coldeos_use_thermal_gamma_law | Scope: restricted | BOOLEAN |

Description: use an additional thermal gamma?
| ||

Default: Yes | ||

do_energy_shift | Scope: restricted | BOOLEAN |

Description: shift energies around?
| ||

Default: yes | ||

dump_nuceos_table | Scope: restricted | BOOLEAN |

Description: Dump table in ASCII at beginning
| ||

Default: No | ||

dump_nuceos_table_name | Scope: restricted | STRING |

Description: nuceos dump table name (ASCII)
| ||

Range | Default: blah.asc | |

.+ | Can be anything
| |

gl_gamma | Scope: restricted | REAL |

Description: Adiabatic Index for gamma-law EOS
| ||

Range | Default: 2.0 | |

(0:* | any positive number
| |

gl_k | Scope: restricted | REAL |

Description: Polytropic constant in c=G=Msun=1 for gamma-law EOS
| ||

Range | Default: 100.0 | |

(0:* | any positive number
| |

hybrid_gamma | Scope: restricted | REAL |

Description: Adiabatic index of each piece
| ||

Range | Default: 0.0 | |

unset
| ||

(0:* | Any positive
| |

hybrid_gamma1 | Scope: restricted | REAL |

Description: subnuclear adiabatic Index for hybrid EOS
| ||

Range | Default: 1.325 | |

(0:* | Any positive
| |

hybrid_gamma2 | Scope: restricted | REAL |

Description: subnuclear adiabatic Index for hybrid EOS
| ||

Range | Default: 2.5 | |

(0:* | Any positive
| |

hybrid_gamma_th | Scope: restricted | REAL |

Description: Thermal gamma for hybrid EOS
| ||

Range | Default: 1.5 | |

(0:* | any positive number
| |

hybrid_k0 | Scope: restricted | REAL |

Description: K coefficent of the first piece
| ||

Range | Default: 0.0 | |

unset
| ||

(0:* | Any positive
| |

hybrid_k1 | Scope: restricted | REAL |

Description: Polytropic constant in c=G=Msun=1 for hybrid EOS
| ||

Range | Default: 0.4640517 | |

(0:* | Any positive
| |

hybrid_rho | Scope: restricted | REAL |

Description: Density values at separation points
| ||

Range | Default: 0.0 | |

unset
| ||

(0:* | Any positive
| |

hybrid_rho_nuc | Scope: restricted | REAL |

Description: Density at which to switch between gammas; c=G=Msun=1
| ||

Range | Default: 3.238607e-4 | |

(0:* | Any positive
| |

n_pieces | Scope: restricted | INT |

Description: Number of polytropic pieces
| ||

Range | Default: 2 | |

1:10 | Max 10 pieces
| |

nuceos_read_table | Scope: restricted | BOOLEAN |

Description: Read in EOS table?
| ||

Default: No | ||

nuceos_table_name | Scope: restricted | STRING |

Description: nuceos table name (hdf5)
| ||

Range | Default: blah.h5 | |

.+ | Can be anything
| |

poly_gamma | Scope: restricted | REAL |

Description: Adiabatic Index for poly EOS
| ||

Range | Default: 2.0 | |

(0:* | any positive number
| |

poly_gamma_initial | Scope: restricted | REAL |

Description: Initial Adiabatic Index for poly and hybrid EOS
| ||

Range | Default: -1 | |

-1 | ”use poly_gamma/hybrid_ga mma, ie no change in gamma
during ID”
| |

(0:* | assume that ID used this adiabiatic index, change K
accordingly
| |

poly_k | Scope: restricted | REAL |

Description: Polytropic constant in c=G=Msun=1
| ||

Range | Default: 100.0 | |

: | ||

read_table_on_single_process | Scope: restricted | BOOLEAN |

Description: read table on one process only and bcast data
| ||

Default: no | ||

reader_process | Scope: restricted | INT |

Description: read table on this process and bcast data
| ||

Range | Default: (none) | |

0:* | read on process N
| |

out_dir | Scope: shared from IO | STRING |

Implements:

eos_omni

Provides:

EOS_Omni_GetHandle to

EOS_Omni_press to

EOS_Omni_press_cs2 to

EOS_Omni_pressOMP to

EOS_Omni_DPressByDEps to

EOS_Omni_DPressByDRho to

EOS_Omni_dpderho_dpdrhoe to

EOS_Omni_cs2 to

EOS_Omni_EpsFromPress to

EOS_Omni_RhoFromPressEpsTempEnt to

EOS_Omni_PressEpsTempYe_from_Rho to

EOS_Omni_short to

EOS_Omni_full to

EOS_Omni_DEpsByDRho_DEpsByDPress to

EOS_Omni_press_f_hrho_v2_rhoW to

EOS_Omni_dpdhrho_f_hrho_v2_rhoW to

EOS_Omni_dpdv2_f_hrho_v2_rhoW to

This section lists all the variables which are assigned storage by thorn EinsteinEOS/EOS_Omni. Storage can either last for the duration of the run (Always means that if this thorn is activated storage will be assigned, Conditional means that if this thorn is activated storage will be assigned for the duration of the run if some condition is met), or can be turned on for the duration of a schedule function.

NONE

CCTK_STARTUP

eos_omni_startup

set up conversion factors and other fun stuff

Language: | fortran | |

Options: | global | |

Type: | function | |

CCTK_BASEGRID (conditional)

nuc_eos_readtable_cactus_wrapper

read eos hdf5 table

Language: | c | |

Options: | global | |

Type: | function | |

CCTK_BASEGRID (conditional)

eos_omni_get_energy_shift

setup energy_shift in eos_omni_module

After: | nuc_eos_readtable_cactus_wrapper | |

Language: | fortran | |

Options: | global | |

Type: | function | |

CCTK_BASEGRID (conditional)

eos_omni_dumptable

dump eos hdf5 table in ascii

After: | nuc_eos_readtable_cactus_wrapper | |

Language: | c | |

Options: | global | |

Type: | function | |

CCTK_BASEGRID (conditional)

eos_omni_readcoldtable

read cold eos ascii table

Language: | fortran | |

Options: | global | |

Type: | function | |

CCTK_BASEGRID (conditional)

eos_omni_barotropicreadtable

read barotropic eos ascii table

Language: | fortran | |

Options: | global | |

Type: | function | |