%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\documentclass[12pt]{article}
\renewcommand{\baselinestretch}{1.6}
\textheight 24.cm
\textwidth 17cm
\oddsidemargin -18pt
\evensidemargin -18pt
\topmargin -35pt
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\begin{document}
\vspace*{0.5cm}
\begin{center}
{\large {\bf User Notes for Relativistic EOS Table} }
\vspace*{1.0cm}
H. Shen$^{a,b,}$\footnote{e-mail: shen@rcnp.osaka-u.ac.jp},
H. Toki$^{a,c,}$\footnote{e-mail: toki@rcnp.osaka-u.ac.jp},
K. Oyamatsu$^{c,d,}$\footnote{e-mail: oyak@postman.riken.go.jp}
and K. Sumiyoshi$^{c,}$\footnote{e-mail: sumi@postman.riken.go.jp} \\
\end{center}
$^{a}$Research Center for Nuclear Physics (RCNP), Osaka University,
Ibaraki, Osaka 567-0047, Japan \\
$^{b}$Department of Physics, Nankai University, Tianjin 300071, China \\
$^{c}$The Institute of Physical and Chemical Research (RIKEN),
Wako, Saitama 351-0198, Japan \\
$^{d}$Department of Energy Engineering and Science, Nagoya University,
Nagoya 464-8603, Japan \\
\vspace*{0.5cm}
\begin{abstract}
This is a detailed note for users of the Relativistic Equation of
State of Nuclear Matter for the use of Nuclear Astrophysics.
\end{abstract}
\vspace*{0.5cm}
{\large {\bf Contents} }
\begin{itemize}
\item 1 Introduction
\item 2 Relativistic mean field theory
\item 3 Ideal gas approximation
\item 4 Thomas-Fermi approximation
\item 5 Calculation procedure in details
\item 6 Resulting EOS table
\item 7 Suggestions for using the EOS table
\end{itemize}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\newpage
\section{Introduction}
We write this document as a guide for the use of the relativistic
equation of state (EOS) table.
This EOS table is constructed in the following range
\begin{itemize}
\item temperature $T (MeV)$:
$ -1.0 \leq \log_{10}(T) \leq 2.0; \,\, $
mesh of $\log_{10}(T) \sim 0.1 $
\item proton fraction $Y_p$:
$ -2.00 \leq \log_{10}(Y_p) \leq -0.25; \,\, $
mesh of $\log_{10}(Y_p) \sim 0.025 $
\item baryon mass density $\rho_B (g/cm^3)$:
$ 5.1 \leq \log_{10}(\rho_B) \leq 15.4; \,\, $
mesh of $\log_{10}(\rho_B) \sim 0.1 $
\end{itemize}
We also add the results for zero temperature case ($T=0$)
and pure neutron matter case ($Y_p=0$).
We have worked out consistent calculations for uniform matter and non-uniform
matter in the relativistic mean field (RMF) framework \cite{STOS981,STOS982}.
We use the Thomas-Fermi approximation to describe inhomogeneous
nuclear matter, which can be modeled as a mixture of
free neutrons, free protons, alpha-particles, and a single species
of heavy nuclei.
For extremely low density and finite temperature,
we approximate the nuclear matter as classical ideal gas
of free neutrons, free protons, and alpha-particles.
Antiparticles have some contribution when temperature is very high.
The thermodynamically favorable state is the one that minimizes
the free energy density in this model, and
we determine the most favorable state of nuclear matter
at each temperature, proton fraction, and baryon mass density.
The leptons can be considered as uniform non-interacting relativistic
particles, which is relatively easy to deal with.
Hence, we provide the baryon EOS without the lepton contribution
in this table, users are supposed to add the lepton
contribution to the baryon EOS.
This document is arranged as follows. In Sect.2, we describe the RMF theory.
In Sect.3, we introduce the ideal gas approximation to be used in
low density range.
In Sect.4, we describe the Thomas-Fermi approximation, which is
used for calculating non-uniform matter.
In Sect.5, we make a detailed description of how to work out the EOS table
in this model. In Sect.6, we list the definitions of the physical quantities in
the EOS table. Sect.7 is devoted to suggestions and discussions
for using the EOS table. For users who are not interested in the framework
and calculating processes, they can only read the Sect.6 and Sect.7
for the purpose to use the EOS table correctly.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\section{ Relativistic mean field theory}
We adopt the RMF theory with non-linear $\sigma$ and $\omega$ terms
to describe homogeneous nuclear matter \cite{ST94}.
We start with the lagrangian given by
\begin{eqnarray}
{\cal L}_{RMF} & = & \bar{\psi}\left[i\gamma_{\mu}\partial^{\mu} -M
-g_{\sigma}\sigma-g_{\omega}\gamma_{\mu}\omega^{\mu}
-g_{\rho}\gamma_{\mu}\tau_a\rho^{a\mu}
\right]\psi \\ \nonumber
&& +\frac{1}{2}\partial_{\mu}\sigma\partial^{\mu}\sigma
-\frac{1}{2}m^2_{\sigma}\sigma^2-\frac{1}{3}g_{2}\sigma^{3}
-\frac{1}{4}g_{3}\sigma^{4} \\ \nonumber
&& -\frac{1}{4}W_{\mu\nu}W^{\mu\nu}
+\frac{1}{2}m^2_{\omega}\omega_{\mu}\omega^{\mu}
+\frac{1}{4}c_{3}\left(\omega_{\mu}\omega^{\mu}\right)^2 \\ \nonumber
&& -\frac{1}{4}R^a_{\mu\nu}R^{a\mu\nu}
+\frac{1}{2}m^2_{\rho}\rho^a_{\mu}\rho^{a\mu} .
\end{eqnarray}
Here, $\psi$ denotes a SU(2) baryon field of mass $M$ (proton and neutron);
$\sigma, \omega^{\mu},$ and $\rho^{a\mu}$ are
$\sigma, \omega,$ and $\rho$ meson fields with masses
$m_{\sigma}, m_{\omega},$ and $m_{\rho}$, respectively.
$W^{\mu\nu}$ and $R^{a\mu\nu}$ are the antisymmetric field tensors
for $\omega^{\mu}$ and $\rho^{a\mu}$, which can be written as
\begin{eqnarray}
W^{\mu\nu} & = & \partial^{\mu}\omega^{\nu}
- \partial^{\nu}\omega^{\mu} , \\
R^{a\mu\nu} & = & \partial^{\mu}\rho^{a\nu}- \partial^{\nu}\rho^{a\mu}
+g_{\rho}\epsilon^{abc}\rho^{b\mu}\rho^{c\nu} .
\end{eqnarray}
In the lagrangian, the constants $g_{\sigma}, g_{\omega}, g_{\rho}$
are the coupling constants
for the interactions between mesons and nucleons, the coefficients $g_2$
and $g_3$ are the self-coupling constants for $\sigma$ meson field,
while $c_3$ is the self-coupling constant for $\omega$ meson field.
It is known that the inclusion of the non-linear $\sigma$ terms is
essential to reproduce the properties of nuclei quantitatively and
a reasonable value for the incompressibility, while the
non-linear $\omega$ term is added in order to reproduce the density
dependence of the vector part of the self-energy of the nucleon obtained in
the RBHF theory.
The lagrangian contains the meson masses, the coupling constants, and the
self-coupling constants as parameters.
We adopt the parameter set TM1 \cite{ST94},
which is listed in Table 1.
With the TM1 parameter set, the symmetry energy is 36.9 MeV and
the incompressibility is 281 MeV.
The RMF theory with the TM1 provides
excellent properties of heavy nuclei including unstable nuclei \cite{ST94},
and it is also shown to reproduce satisfactory agreement with experimental data
in the studies of the nuclei with deformed configuration \cite{HS96} and the
giant resonances within the RPA formalism \cite{MT97}.
The properties of the EOS with TM1 at high density were discussed in details
in Ref.\cite{SKT95}, and the neutron star properties were studied using this
EOS in Ref.\cite{STOS982, SKT95}.
\begin{table}
\caption{The parameter set TM1 for the lagrangian in Eq.(1)}
\vspace{0.5cm}
\begin{center}
\begin{tabular}{|l l l|r r r|} \hline
& Parameter & & & TM1 & \\ \hline
& $M (MeV)$ & & & 938.0 & \\
& $m_{\sigma} (MeV)$ & & & 511.19777 & \\
& $m_{\omega} (MeV)$ & & & 783.0 & \\
& $m_{\rho} (MeV)$ & & & 770.0 & \\
& $g_{\sigma}$ & & & 10.02892 & \\
& $g_{\omega}$ & & & 12.61394 & \\
& $g_{\rho} $ & & & 4.63219 & \\
& $g_{2} (fm^{-1})$ & & & -7.23247 & \\
& $g_{3} $ & & & 0.61833 & \\
& $c_{3} $ & & & 71.30747 & \\ \hline
\end{tabular}
\end{center}
\end{table}
Starting with the lagrangian, we derive a set of the Euler-Lagrange
equations. The Dirac equation for the nucleon field is given by
\begin{equation}
\left(i\gamma_{\mu}\partial^{\mu} -M
-g_{\sigma}\sigma-g_{\omega}\gamma_{\mu}\omega^{\mu}
-g_{\rho}\gamma_{\mu}\tau_a\rho^{a\mu}
\right)\psi=0 ,
\end{equation}
and the Klein-Gordon equations for the meson fields are given by
\begin{eqnarray}
\partial_{\nu}\partial^{\nu}\sigma + m_{\sigma}^2\sigma & =
&-g_{\sigma}\bar{\psi}\psi
-g_{2}\sigma^{2}-g_{3}\sigma^{3}, \\
\partial_{\nu}W^{\nu\mu} +m_{\omega}^2\omega^{\mu} & =
& g_{\omega}\bar{\psi}\gamma^{\mu}\psi
-c_{3}\left(\omega_{\nu}\omega^{\nu}\right)\omega^{\mu}, \\
\partial_{\nu}R^{a\nu\mu} +m_{\rho}^2\rho^{a\mu} & =
& g_{\rho}\bar{\psi}\tau_a\gamma^{\mu}\psi
+g_{\rho}\epsilon_{abc}\rho^{b}_{\nu}R^{c\nu\mu} .
\end{eqnarray}
These equations are coupled nonlinear quantum field equations, which
are very difficult to solve exactly. We adopt the relativistic
mean field approximation as described in Ref.\cite{SW86},
in which the meson fields are treated as classical fields,
and the field operators $\sigma, \omega^{\mu},$ and $\rho^{a\mu}$
are replaced by their expectation values
$\langle\sigma\rangle, \langle\omega^{\mu}\rangle,$
and $\langle\rho^{a\mu}\rangle$.
In the present study, we consider first static infinite matter so that
we obtain simplified equations, where the derivative terms
in the Klein-Gordon equations vanish automatically
due to the translational invariance of infinite matter.
The spatial components of the vector meson fields vanish under the rotational
symmetry. For the isovector-vector meson field $\rho^{a\mu}$,
only the third isospin
component has a nonvanishing value because of the charge conservation.
Hence, the equations for meson fields are reduced to
\begin{eqnarray}
\sigma_0 & \equiv & \langle\sigma\rangle =
-\frac{g_{\sigma}}{m_{\sigma}^2} \langle\bar{\psi}\psi\rangle
-\frac{1}{m_{\sigma}^2}\left(g_{2}\sigma_0^{2}+g_{3}\sigma_0^{3}\right)
, \\
\omega_0 & \equiv & \langle\omega^0\rangle =
\frac{g_{\omega}}{m_{\omega}^2} \langle\bar{\psi}\gamma^0\psi\rangle
-\frac{1}{m_{\omega}^2}c_3\omega_0^3
, \\
\rho_0 & \equiv & \langle\rho^{30}\rangle =
\frac{g_{\rho}}{m_{\rho}^2} \langle\bar{\psi}\tau_3\gamma^0\psi\rangle .
\end{eqnarray}
The Dirac equation then becomes
\begin{equation}
\left(-i\alpha_k\nabla^k + \beta M^{*}
+g_{\omega}\omega_0 +g_{\rho}\tau_3\rho_0 \right)\psi_{is}=
\varepsilon_{is} \psi_{is} ,
\end{equation}
where the index $i$ denotes the isospin degree of freedom (proton and neutron)
and $s$ denotes the index of eigenstates of nucleon,
$M^{*} \equiv M+g_{\sigma}\sigma_0$ is the effective mass,
$\varepsilon_{is}$ is the single-particle energy.
Nucleons occupy single-particle orbits with the occupation probability
$f_{is}$. At zero temperature, $f_{is}=1$ under Fermi surface,
while $f_{is}=0$ above Fermi surface.
For finite temperature, the occupation probability is given by Fermi-Dirac
distribution,
\begin{eqnarray}
f_{is}=\frac{1}{1+\exp\left[\left(\varepsilon_{is}-\mu_{i}\right)/T\right]}
=\frac{1}{1+\exp\left[\left(\sqrt{k^2+{M^{*}}^2}-\nu_{i}\right)
/T\right]}
, \\
f_{\bar{i}s}=\frac{1}{1+\exp\left[\left(-\varepsilon_{\bar{i}s}+\mu_{i}
\right)/T\right]}
=\frac{1}{1+\exp\left[\left(\sqrt{k^2+{M^{*}}^2}+\nu_{i}\right)/T\right]} .
\end{eqnarray}
Here the index $i$ and $\bar{i}$ denote nucleons and antinucleons,
$\varepsilon_{is}$ and $\varepsilon_{\bar{i}s}$ are the single-particle
energy of nucleons and antinucleons, respectively.
The relation between the chemical potential $\mu_i$ and the kinetic
part of the chemical potential $\nu_i$ is give by
\begin{equation}
\mu_{i}=\nu_i+g_{\omega}\omega_0 +g_{\rho}\tau_3\rho_0.
\end{equation}
The chemical potential $\mu_i$ is related to the number density of nucleon
$n_i$ as
\begin{equation}
n_{i}=\frac{\gamma}{2\pi^2}
\int_0^{\infty} dk\,k^2\,\left(f_{ik}-f_{\bar{i}k}\right),
\end{equation}
where $\gamma$ is the degeneracy factor for the spin degree of freedom
($\gamma=2$),
the quantum number $s$ is replaced by the momentum $k$ when we do the
integration in the momentum space instead of the summation of eigenstates.
The equations of meson fields can be written as
\begin{eqnarray}
\sigma_0 &=& -\frac{g_{\sigma}}{m_{\sigma}^2} \displaystyle{\sum_i
\frac{\gamma}{2\pi^2} \int_0^{\infty} dk\,k^2\,
\frac{M^{*}}{\sqrt{k^2+{M^*}^2}}
\left(f_{ik}+f_{\bar{i}k}\right) }
-\frac{1}{m_{\sigma}^2}\left(g_{2}\sigma_0^{2}+g_{3}\sigma_0^{3}\right)
, \\
\omega_0 &=&
\frac{g_{\omega}}{m_{\omega}^2}
\left(n_n+n_p\right)
-\frac{c_3}{m_{\omega}^2}\omega_0^3
, \\
\rho_0 &=&
\frac{g_{\rho}}{m_{\rho}^2}
\left(n_n-n_p\right) .
\end{eqnarray}
Here, $n_p$ and $n_n$ are the proton number density and neutron number density
as defined in Eq.(15), and we denote $n_B=n_p+n_n$ as the baryon number
density. We solve these equations self-consistently.
The thermodynamical quantities are described
in Ref.\cite{SW86}, and we just write the expressions here.
The energy density of nuclear matter is given by
\begin{eqnarray}
\epsilon &=& \displaystyle{\sum_i \frac{\gamma}{2\pi^2}
\int_0^{\infty} dk\,k^2\,
\sqrt{k^2+{M^*}^2} \left(f_{ik}+f_{\bar{i}k}\right) }
+\frac{1}{2}m_{\sigma}^2\sigma_0^2+\frac{1}{3}g_{2}\sigma_0^{3}
+\frac{1}{4}g_{3}\sigma_0^{4} \\ \nonumber
& & +g_{\omega}\omega_0 \left(n_n+n_p\right)
-\frac{1}{2}m_{\omega}^2\omega_0^2-\frac{1}{4}c_{3}\omega_0^{4}
\\ \nonumber
& & +g_{\rho}\rho_0 \left(n_n-n_p\right)
-\frac{1}{2}m_{\rho}^2\rho_0^2,
\end{eqnarray}
the pressure of nuclear matter is given by
\begin{eqnarray}
p &=& \displaystyle{\sum_i \frac{\gamma}{6\pi^2} \int_0^{\infty} dk\,k^2\,
\frac{k^2}{\sqrt{k^2+{M^*}^2}} \left(f_{ik}+f_{\bar{i}k}\right) }
-\frac{1}{2}m_{\sigma}^2\sigma_0^2-\frac{1}{3}g_{2}\sigma_0^{3}
-\frac{1}{4}g_{3}\sigma_0^{4} \\ \nonumber
& & +\frac{1}{2}m_{\omega}^2\omega_0^2+\frac{1}{4}c_{3}\omega_0^{4}
+\frac{1}{2}m_{\rho}^2\rho_0^2,
\end{eqnarray}
and the entropy density is calculated by
\begin{equation}
s=\displaystyle{\sum_i \frac{\gamma}{2\pi^2} \int_0^{\infty} dk\,k^2\,
\left[-f_{ik}\ln f_{ik}-\left(1-f_{ik}\right)\ln \left(1-f_{ik}\right)
-f_{\bar{i}k}\ln f_{\bar{i}k}-\left(1-f_{\bar{i}k}\right)
\ln \left(1-f_{\bar{i}k}\right) \right] }.
\end{equation}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\section{Ideal gas approximation}
The interaction between nucleons is negligible at low density,
so that we can treat nucleons as non-interacting Boltzmann particles
at low density. We use the ideal gas approximation
to describe the mixed gas of neutrons, protons, and alpha-particles
at low density and finite temperature. We note that the connection
between the RMF results and the ideal gas approximation is smooth.
For the Boltzmann particles with spin $\frac{1}{2}$, mass $M$, and number
density $n_i$ ($i$ is $n$ or $p$),
we start with the partition function given by
\begin{equation}
Z_i=\left[ (n_Q)^{n_i/2}/(n_i/2)! \right]^2 \, ,
\end{equation}
where we have used the abbreviation
$n_Q=\left(\frac{M\,T}{2\pi\hbar^2}\right)^{3/2}$.
The factor $(n_i/2)!$ comes from avoiding the double counting for
$n_i/2$ indistinguishable particles of spin up or spin down.
We can calculate the free energy density from the partition
function as
\begin{equation}
f_i=-T\ln Z_i=-T\,n_i\left[\ln (2n_Q/n_i)+1\right].
\end{equation}
The entropy density is given by
\begin{equation}
s_i=-\left(\frac{\partial f_i}{\partial T}\right)_{n_i}
=n_i\left[\ln (2n_Q/n_i)+\frac{5}{2}\right],
\end{equation}
the internal energy density can be written as
\begin{equation}
\epsilon_i=f_i+T s_i=\frac{3}{2} T\,n_i,
\end{equation}
the chemical potential is obtained through the relation
\begin{equation}
\mu_i=\left(\frac{\partial f_i}{\partial n_i}\right)_{T}
=-T\ln (2n_Q/n_i),
\end{equation}
and the pressure is given by
\begin{equation}
p_i=\left[n_i^2\frac{\partial \left(f_i/n_i\right)}{\partial n_i}\right]_{T}
=T\,n_i.
\end{equation}
\\
The partition function of alpha-particles can be written as
\begin{equation}
Z_{\alpha}=\left[8 n_Q \exp \left(B_{\alpha}/T\right)\right]^{n_{\alpha}}/
n_{\alpha}! \, ,
\end{equation}
where $n_{\alpha}$ is the alpha-particle number density,
$B_{\alpha}=28.3 MeV$ is the binding energy of an alpha-particle
relative to free neutrons taken from Ref.\cite{LS91},
then the thermodynamical quantities of alpha-particles are given by
\begin{eqnarray}
f_{\alpha} &=&-T\,n_{\alpha}\left[\ln (8n_Q/n_{\alpha})+1\right]
-n_{\alpha} B_{\alpha}, \\
s_{\alpha} &=& n_{\alpha}\left[\ln (8n_Q/n_{\alpha})+\frac{5}{2}\right], \\
\epsilon_{\alpha} &=& \frac{3}{2} T\,n_{\alpha}
-n_{\alpha} B_{\alpha}, \\
\mu_{\alpha} &=& -T\ln (8n_Q/n_{\alpha})-B_{\alpha}, \\
p_{\alpha} &=& T\,n_{\alpha}.
\end{eqnarray}
For the mixed gas with neutron number density $n_n$,
proton number density $n_p$, and alpha-particle number density $n_{\alpha}$,
the free energy density is given by
\begin{equation}
f=f_n+f_p+f_{\alpha}.
\end{equation}
We have to take into account the volume of alpha-particle, otherwise the
alpha-particle fraction will become some big number at high
densities, where the alpha-particles should actually disappear.
After we take into account the volume excluded by the alpha-particles,
the free energy densities in Eq.(34) are given by
\begin{eqnarray}
f_n &=& -(1-u)\,T\,\tilde{n}_n\left[\ln (2n_Q/\tilde{n}_n)+1\right], \\
f_p &=& -(1-u)\,T\,\tilde{n}_p\left[\ln (2n_Q/\tilde{n}_p)+1\right], \\
f_{\alpha} &=& -(1-u)\,\left\{ T\,\tilde{n}_{\alpha}
\left[\ln (8n_Q/\tilde{n}_{\alpha} )
+1\right] -\tilde{n}_{\alpha} B_{\alpha} \right\} ,
\end{eqnarray}
where $u=n_{\alpha}\,v_{\alpha}$ is the fraction of space occupied by
alpha-particles, the effective volume of an alpha-particle,
$v_{\alpha}=24 fm^{-3}$, is taken from Ref.\cite{LS91}.
We denote the effective number densities of neutrons, protons,
and alpha-particles as $\tilde{n}_n=n_n/(1-u)$, $\tilde{n}_p=n_p/(1-u)$,
and $\tilde{n}_{\alpha}=n_{\alpha}/(1-u)$.
The inclusion of the volume excluded by the alpha-particles
has negligible effect in the low density region,
while it is necessary for the calculation at high density.
Actually, we just use Eqs.(35) and (36) for the case of
$\tilde{n}_n+\tilde{n}_p < 10^{-5} fm^{-3}$, we have to use the RMF theory
to describe nucleons when
$\tilde{n}_n+\tilde{n}_p > 10^{-5} fm^{-3}$, where the nucleon free
energy density is obtained by the RMF input table.
We use Eq.(37) to calculate the free energy contributed by
alpha-particles for whole density range,
it is considered as a reasonable approximation because the number
density of alpha-particles is not so large even though the average
baryon density is high, and the alpha-particle fraction tends to
vanish in the high density limit.
For very high temperature, we need to consider the contribution
from antiparticles. The alpha-particle fraction is very small
at high temperature, so we can treat the matter as uniform matter
of nucleons and antinucleons. We have included the freedom
of antiparticles in the RMF theory, but the RMF code meets some difficulty
at extremely low density, so we need to take some approximation
in order to describe the range of low density and high temperature.
We take the non-relativistic approximation of the RMF theory
at low density, where the effective mass $M^{*} \approx M$,
$ \approx 0$,
and $ \sqrt{k^2+{M^{*}}^2} \approx M+\frac{k^2}{2M}$. Then
the occupation probability of the Fermi-Dirac distribution become
\begin{eqnarray}
f_{is}=A C_i \exp\left(-k^2/2MT\right)
, \\
f_{\bar{i}s}=\frac{A}{C_i} \exp\left(-k^2/2MT\right),
\end{eqnarray}
where $A=\exp\left(-M/T\right)$, and $C_i=\exp\left(\mu_i/T\right)$
with $\mu_i$ being the chemical potential of protons $(i=p)$
or neutrons $(i=n)$,
which is related to the number density $n_i$ as
\begin{equation}
n_i=2A\left(C_i-\frac{1}{C_i}\right)\,\left(\frac{MT}{2\pi}\right)^{3/2}.
\end{equation}
The thermodynamical quantities can be approximated as
\begin{eqnarray}
\epsilon_i &=& 2A\left(C_i+\frac{1}{C_i}\right)\,
\left(\frac{MT}{2\pi}\right)^{3/2}\,
\left(\frac{3}{2}T\right) , \\ \nonumber
p_i &=& 2A\left(C_i+\frac{1}{C_i}\right)\,
\left(\frac{MT}{2\pi}\right)^{3/2}\,
T , \\ \nonumber
s_i &=& 2A\left[\frac{5}{2} \left(C_i+\frac{1}{C_i}\right)
-C_i\ln \left(AC_i\right)-\frac{1}{C_i}
\ln \left(\frac{A}{C_i}\right) \right]\,
\left(\frac{MT}{2\pi}\right)^{3/2}.
\end{eqnarray}
The contribution from antiparticles is negligible when $C_i\gg \frac{1}{C_i}$,
where the above expressions agree with the ideal gas approximation.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\section{Thomas-Fermi approximation}
In the range, $T < 15 MeV$ and $\rho_B < 10^{14.2}g/cm^3$, where heavy nuclei
may be formed in order to lower the free energy,
we perform the Thomas-Fermi calculation based on the work done
by Oyamatsu \cite{O93}. In this case, the non-uniform matter can be modeled
as a mixture of free neutrons, free protons, alpha-particles, and a single
species of heavy nuclei, while the leptons can be treated as uniform
non-interacting particles separately.
For the system with fixed proton fraction, the leptons play no role in the
minimization of the free energy. Hence we mainly
pay attention to baryon contribution in this calculation.
Hereafter, we will not mention the leptons frequently, while we should
keep in mind that there exists an uniform lepton gas existing everywhere.
We assume that each heavy spherical nucleus is
located in the center of a charge-neutral cell consisting of a
vapor of neutrons, protons and alpha-particles. The nuclei form the
body-centered-cubic (BCC) lattice to minimize the Coulomb lattice energy.
It is useful to introduce the Wigner-Seitz cell to simplify the energy of
the unit cell. The Wigner-Seitz cell is a sphere whose volume is
the same as the unit cell in the BCC lattice.
We define the lattice constant $a$ as the cubic
root of the cell volume. Then, we have
\begin{equation}
V_{cell}=a^3=N_B / n_B,
\end{equation}
where $N_B$ and $n_{B}$ are the
baryon number per cell and the average baryon number density, respectively.
We define the baryon mass density as $\rho_B=m_{u} n_B$ with
$m_{u}$ being the atomic mass unit.
We calculate the Coulomb energy using this Wigner-Seitz approximation and
add a correction energy for the BCC lattice \cite{O93}. This correction
energy is negligible unless the nuclear size is comparable to the cell size.
We assume the nucleon distribution functions $n_i(r)$
($i=n$ for neutron, $i=p$ for proton)
in the Wigner-Seitz cell as
\begin{equation}
n_i\left(r\right)=\left\{
\begin{array}{ll}
\left(n_i^{in}-n_i^{out}\right) \left[1-\left(\frac{r}{R_i}\right)^{t_i}
\right]^3 +n_i^{out}, & 0 \leq r \leq R_i \\
n_i^{out}, & R_i \leq r \leq R_{cell} \\
\end{array} \right. ,
\end{equation}
where $r$ represents the distance from the center of the nucleus and
$R_{cell}$ is the radius of the Wigner-Seitz cell defined by the relation,
\begin{equation}
V_{cell} \equiv \frac{4 \pi}{3} R_{cell}^3.
\end{equation}
The density parameters
$n_i^{in}$ and $n_i^{out}$ are
the densities at $r=0$ and $ r \geq R_i $, respectively. The parameters
$R_i$ and $t_i$ determine the boundary and the relative surface thickness
of the heavy nucleus.
For the distribution function of alpha-particle $n_{\alpha}(r)$,
which should decrease as $r$ approaches the center of the heavy nucleus,
we assume
\begin{equation}
n_{\alpha}\left(r\right)=\left\{
\begin{array}{ll}
-n_{\alpha}^{out} \left[1-\left(\frac{r}{R_p}\right)^{t_p}
\right]^3 +n_{\alpha}^{out}, & 0 \leq r \leq R_p \\
n_{\alpha}^{out}, & R_{p} \leq r \leq R_{cell} \\
\end{array} \right. ,
\end{equation}
which will give $n_{\alpha} (r=0) =0$
and $n_{\alpha} (r>R_p) =n_{\alpha}^{out}$.
Here we use the same parameters $R_p$ and $t_p$ for both proton distribution
function and alpha-particle distribution function in order to avoid
too many parameters in the minimization procedure.
The parameters $R_n$ and $t_n$ may be a little different from
$R_p$ and $t_p$ due to the additional neutrons forming a neutron skin
in the surface region.
For the system with fixed temperature $T$, proton fraction $Y_p$,
and baryon mass density $\rho_B$, there are eight independent parameters
in the ten variables; $a, n_n^{in}, n_n^{out}, R_n, t_n,
n_p^{in}, n_p^{out}, R_p, t_p, n_{\alpha}^{out}$.
The thermodynamically favorable state is the one that minimizes
the free energy density with respect to these eight independent parameters.
In this model the free energy density contributed by baryons is given as
\begin{equation}
f\,=\,F_{cell}\,/\,a^3
=\,\left(\,E_{cell}\,-\,T\,S_{cell}\,\right)\,/\,a^3 ,
\end{equation}
where the free energy per cell $F_{cell}$ can be written as
\begin{equation}
F_{cell}=(E_{bulk}+E_s+E_C)- T S_{cell}=F_{bulk}+E_s+E_C.
\end{equation}
The bulk energy $E_{bulk}$, entropy $S_{cell}$, and bulk free energy
$F_{bulk}$ are calculated by
\begin{eqnarray}
E_{bulk} &=&\int_{cell} \epsilon \left( \, n_n\left(r\right),
\, n_p\left(r\right), \, n_{\alpha}\left(r\right) \, \right) d^3r, \\
S_{cell} &=&\int_{cell} s \left(\, n_n\left(r\right),
\, n_p\left(r\right), \, n_{\alpha}\left(r\right) \, \right) d^3r, \\
F_{bulk} &=&\int_{cell} f \left( \, n_n\left(r\right),
\, n_p\left(r\right), \, n_{\alpha}\left(r\right) \, \right) d^3r.
\end{eqnarray}
Here $\epsilon \left( \, n_n\left(r\right),
\, n_p\left(r\right), \, n_{\alpha}\left(r\right) \, \right)$,
$s \left( \, n_n\left(r\right),
\, n_p\left(r\right), \, n_{\alpha}\left(r\right) \, \right)$,
and $f \left( \, n_n\left(r\right),
\, n_p\left(r\right), \, n_{\alpha}\left(r\right) \, \right)$
are the local energy density, entropy density, and free energy density
at the radius $r$, where it can be considered as a mixed uniform matter
of neutrons, protons, and alpha-particles.
Note that these local densities at each radius
are calculated by treating neutrons and protons in RMF theory
as described in Sect.2
for the case of $n_n(r)+n_p(r) > 10^{-5} fm^{-3}$,
while we can use the method described in the Sect.3 for low density case of
$n_n(r)+n_p(r) < 10^{-5} fm^{-3}$.
As for the surface energy term $E_s$ due to the inhomogeneity of nucleon
distribution, we take a simple form as,
\begin{equation}
E_s=\int_{cell} F_0 \mid \nabla \left( \, n_n\left(r\right)+
n_p\left(r\right) \, \right) \mid^2 d^3r.
\end{equation}
The parameter $F_0=70 \, MeV\cdot fm^5$ is determined by doing the Thomas-Fermi
calculations of finite nuclei so as to reproduce the gross properties of
nuclear mass, charge radii and the beta stability line as described in
the appendix in Ref.\cite{O93}.
The Coulomb energy per cell $E_C$ is calculated using the Wigner-Seitz
approximation and add a correction term for the BCC lattice \cite{O93}
\begin{equation}
E_C=\frac{1}{2}\int_{cell} e \left[n_p\left(r\right)
+2n_{\alpha}\left(r\right)-n_e\right]\,\phi(r) d^3r
\,+\,\triangle E_C,
\end{equation}
where $\phi(r)$ stands for the electrostatic potential calculated
in the Wigner-Seitz approximation,
$n_e$ is the electron number density of uniform electron gas ($n_e=Y_p\,n_B$),
$\triangle E_C$ is the correction term for the BCC lattice,
which can be approximated as
\begin{equation}
\triangle E_C=C_{BCC}\frac{(Z_{non}e)^2} {a},
\end{equation}
here $a$ is the lattice constant as defined in Eq.(42),
the coefficient $C_{BCC}=0.006562$ is taken from Ref.\cite{O93},
$Z_{non}$ is the non-uniform part of the charge number per cell
given by
\begin{equation}
Z_{non}=\int_{0}^{R_p}
(n_p^{in}-n_p^{out}-2 n_{\alpha}^{out})
\left[1-\left(\frac{r}{R_p}\right)^{t_p}\right]^3
4\pi r^2 dr.
\end{equation}
Because of the long-range nature of the Coulomb interaction,
the Coulomb energy will depend on the lattice type.
This dependence was discussed in more details in Ref.\cite{O93}.
The system prefers the BCC lattice because the BCC lattice
gives the lowest Coulomb energy.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\section{Calculation procedure in details}
We calculate the EOS table in the range mentioned in Sect.1.
For each $T$, $Y_p$, and $\rho_B$, we perform the minimization of the free
energy for non-uniform matter using the Thomas-Fermi method described
in Sect.4.
We also do the minimization for the uniform matter with respect to
converting two protons and two neutrons into an alpha-particle,
there is only one independent parameter, alpha-particle fraction $X_{\alpha}$,
in this minimization procedure.
Here the phase of heavy nuclei formed together with free nucleons
and alpha-particles is referred to as non-uniform matter
while the phase of nucleons mixed with alpha-particles without
heavy nuclei is referred to as uniform matter.
By comparing the free energies of the uniform matter and non-uniform matter,
we determine the most favorable state of nuclear matter
at each $T$, $Y_p$, and $\rho_B$.
The density of the phase transition between uniform matter and non-uniform
matter depends on both $T$ and $Y_p$.
The non-uniform matter phase exists only in low temperature range ($T<15 MeV$).
We construct the EOS table in three dimensions ($T$, $Y_p$, and $\rho_B$),
first fix $T$, second fix $Y_p$, third fix $\rho_B$.
We should finish the following steps in order to work out one EOS table
with fixed $T$ where the non-uniform matter phase exists.
\begin{itemize}
\item 1 calculate input table in RMF theory for the use of minimization
\item 2 do minimization for uniform matter at low density
\item 3 do minimization for non-uniform matter at medium density
\item 4 calculate RMF result for uniform matter without alpha-particles
at high density
\item 5 determine the phase transition and connect the free energy
\item 6 calculate pressure and chemical potentials from free energy
\item 7 construct the finial EOS table at $T$
\end{itemize}
Now, we explain each step in more details.
$\bullet$ 1
The calculation of the minimization procedure needs the free energy
per baryon as a function of $n_n$ and $n_p$ as input table,
and this input table should be big enough in order to get good numerical
accuracy in doing linear approximation.
We calculate this input table in the RMF theory which describes
the homogeneous matter of protons and neutrons.
This input table is covering the following range
\begin{eqnarray}
\nonumber
\begin{array}{lllll}
Y_p=\frac{n_p}{n_n+n_p} & from\,\, 0.0 & to\,\, 0.5
& mesh=0.0005 & point=1001 \\
n_B=n_n+n_p (fm^{-3})
& from\,\, 0.000010 & to\,\, 0.000100 & mesh=0.000001 & point=91 \\
& from\,\, 0.000110 & to\,\, 0.001000 & mesh=0.000010 & point=90 \\
& from\,\, 0.001100 & to\,\, 0.060000 & mesh=0.000100 & point=590 \\
& from\,\, 0.061000 & to\,\, 0.160000 & mesh=0.001000 & point=100 \\
\end{array}
\end{eqnarray}
The density meshes are determined by the behavior of the free energy
per baryon, which is necessary to get good accuracy for the
minimization results.
In this input table, we list the quantities of the energy per baryon $E$
and the entropy per baryon $S$, then free energy per baryon is
$F=E-TS$.
The data at $n_n+n_p>0.16$ is not necessary for the minimization codes,
while we use ideal gas approximation to describe neutrons and protons
when $n_n+n_p<0.00001 fm^{-3}$.
$\bullet$ 2
For low density and finite temperature, the thermodynamically favorable
state might be the uniform matter, which is a mixed gas of neutrons,
protons, and alpha-particles. For each $Y_p$ and $\rho_B$, we have to do
minimization in order to find out the alpha-particle fraction having
the minimum free energy. The free energy of the mixed gas is calculated
by Eq.(34) as a function of $n_n$, $n_p$, and $n_{\alpha}$, where
the RMF input table is used instead of Eqs.(35) and (36) when
$n_n+n_p > 10^{-5} fm^{-3}$.
Because $n_B=n_n+n_p+4n_{\alpha}$ and $Y_p=(n_p+2n_{\alpha})/n_B$
are fixed in this case, there is only one independent parameter,
the alpha-particle fraction $X_{\alpha}=4n_{\alpha}/n_B$ as
the independent parameter in this minimization procedure.
$\bullet$ 3
The non-uniform matter phase exists in some middle density range
where the heavy nuclei can be formed in order to lower the free
energy. The starting density of the non-uniform matter phase
depends on the temperature strongly, while the ending density
is nearly independent of the temperature.
Both of them have some weak dependence on $Y_p$.
So we should first try to find out the starting density at $T$,
then calculate the non-uniform matter results in this range by
using the minimization codes with respect to some independent
parameters as described in Sect.4.
Actually, we have to use a few codes to do this calculation
because of the numerical difficulty in minimization procedure.
As the density $\rho_B$ increases with fixed $T$ and $Y_p$,
the heavy nucleus fraction
increases, while the nucleon fractions and alpha-particle
fraction decrease. When one of them decreases to some small number
like $10^{-5}$, it will bring some difficulty for the minimization code,
and lose good accuracy in the results.
Then we have to use another code which take this composition out
by deleting one independent parameter.
It may run into difficulty again as density increases, then we need to
delete one more parameter. We list all codes used in this procedure in table 2.
We determine which code we should use by comparing the free energies,
we choose the one with lowest free energy.
Generally speaking, the alpha-particles will disappear first,
then free neutrons ($Y_p>0.45$) or free protons ($Y_p<0.45$) disappear,
at the end only heavy nuclei are formed without free particles outside
(this pure nucleus phase doesn't occur at $Y_p < 0.3$).
When $T\geq 5 MeV$, the calculation becomes relatively easy because
the favorable state is always the one including all compositions.
As temperature increases, the heavy nucleus fraction decreases,
and it will disappear completely when $T \geq 15 MeV$.
\begin{table}
\caption{The minimization codes for non-uniform matter phase}
\vspace{0.5cm}
\begin{center}
\begin{tabular}{|c|llll|c|} \hline
code & compositions & of non-uniform &matter & & parameter \\
name & & & & & number \\ \hline
npah &free neutrons, &free protons, &alpha particles, & heavy nuclei & 8 \\
pah & &free protons, &alpha particles, & heavy nuclei & 7 \\
nah &free neutrons, & &alpha particles, & heavy nuclei & 7 \\
nph &free neutrons, &free protons, & & heavy nuclei & 7 \\
nh &free neutrons, & & & heavy nuclei & 6 \\
ph & &free protons, & & heavy nuclei & 6 \\
ah & & &alpha particles, & heavy nuclei & 6 \\
h & & & & heavy nuclei & 5 \\
\hline
\end{tabular}
\end{center}
\end{table}
$\bullet$ 4
As density $\rho_B$ increases beyond $10^{14.2} g/cm^3$,
the non-uniform matter phase disappears, it becomes homogeneous matter
of neutrons and protons with less alpha-particles.
We switch off alpha-particle freedom when $X_{\alpha}< 10^{-4}$,
then we can use RMF code to do the calculation in this high density range.
$\bullet$ 5
We judge the phase transitions in whole density range at each $Y_p$ with
fixed $T$ by comparing the free energies calculated with the uniform matter
codes and non-uniform matter codes in table 2,
all the coming results in step 2, 3, and 4 will be the input files for this
judging code. We get two output files from it, eos.fe and eos.tab,
eos.fe includes only free energy per particle $F$ which will be used
for calculating pressure and chemical potentials, eos.tab includes all
quantities which are needed for the finial EOS table.
$\bullet$ 6
We calculate pressure and chemical potentials through the following
thermodynamical relations
\begin{eqnarray}
p &=&\left[\,n_B^2 (\partial F/\partial n_B)\, \right]_{T,Y_p}, \\
\mu_n&=&\left[\,\partial (n_B F) /\partial n_n \, \right]_{T,n_p}, \\
\mu_p&=&\left[\,\partial (n_B F) /\partial n_p \, \right]_{T,n_n},
\end{eqnarray}
where F is the free energy per baryon from the file eos.fe obtained in
step 5, $n_B$ is the baryon number density related to the baryon mass
density $\rho_B$ as $\rho_B=m_{u} n_B$ with $m_{u}=931.49432 MeV$
being the atomic mass unit, $n_n=(1-Y_p)n_B$ and $n_p=Y_p n_B$ are the
neutron number density and proton number density.
The numerical derivatives are calculated by five point differentiation
method, where the five points $F$ around $\rho_B$
and $Y_p$ are obtained by spline method.
For pressure, we use one dimension spline method because the differentiation
in Eq.(55) is under the condition of fixed $Y_p$. We have to use two
dimension spline method to get the five point $F$ for calculating
chemical potential which is under the condition of fixed $n_p$ or fixed $n_n$.
In high density range where we switch off alpha-particles freedom,
we take the exact results calculated in the RMF theory
by Eq.(20) and Eq.(14) instead of the numerical differentiations.
The output files of this step are pre.tab for pressure and
che.tab for chemical potentials.
$\bullet$ 7
We combine the files eos.tab, pre.tab, and che.tab in order to
get the finial EOS table at $T$.
Because of the use of many codes in table 2 for non-uniform matter
phase, some fractions might be equal to zero at some densities in the file
eos.tab, this is due to some compositions taken out for keeping good accuracy.
However, we know the fraction should be some finite number
at finite temperature if the chemical potential is finite.
In this case, we calculate the fraction $X_i$ by chemical potential $\mu_i$
and write it in the finial EOS table instead of $X_i=0$ in eos.tab
\begin{equation}
X_i=(n_i^{out} V^{out})/(n_B V_{cell}),
\end{equation}
where $V_{cell}$ is the cell volume as defined in Eq.(42),
$V^{out}=V_{cell}-\frac{4\pi}{3}R_{A}^3$ is the volume outside
the heavy nucleus with $R_{A}$ being the maximum between $R_p$ and $R_n$
considered as the boundary of the heavy nucleus.
$n_i^{out}$ is the free neutron number density ($i=n$) or
free proton number density ($i=p$). For alpha-particle ($i=\alpha$),
$n_i^{out}$ should be 4 times of the alpha-particle number density outside,
because there are 2 protons and 2 neutrons inside one alpha-particle.
The number densities outside the heavy nucleus can be obtained through
the relations,
\begin{eqnarray}
n_n &=& 2 \left(\frac{M\,T}{2\pi\hbar^2}\right)^{3/2}
\exp (\frac{\mu_n}{T}), \\
n_p &=& 2 \left(\frac{M\,T}{2\pi\hbar^2}\right)^{3/2}
\exp (\frac{\mu_p}{T}), \\
n_{\alpha} &=& 8 \left(\frac{M\,T}{2\pi\hbar^2}\right)^{3/2}
\exp (\frac{\mu_{\alpha}+B_{\alpha}}{T}),
\end{eqnarray}
here the alpha-particle chemical potential is given by
the equilibrium condition, $\mu_{\alpha}=2\mu_n+2\mu_p$.
The finial EOS table at this $T$ is the output file of this step.
\\
We have to go through these steps for each $T$ where non-uniform matter
phase exists. We also need to add the $T=0$ case, which is calculated
by similar procedure as the above steps. There are two difference between
$T=0$ and $T \neq 0$ cases. (1) There are no free protons and alpha-particles
outside the heavy nucleus at $T=0$, so we needn't calculate so many codes
in table 2. The starting density of the non-uniform matter phase is below
$10^{5} g/cm^3$ at $T \leq 0.4 MeV$ and $T=0$, so we can skip the step 2
in this case (2) The ideal gas approximation is not available for $T=0$,
we use a function $C k_f^2$ to express the non-relativistic
kinetic energy at low density instead of $\frac{3}{2}T$ in ideal gas
approximation, note that the free energy at $T=0$ is equal to the
internal energy, the coefficient $C$ is fixed by the first point
in the RMF table in order to get smooth connection.
When the non-uniform matter phase disappears at high temperature,
the calculation becomes relatively easy, the step 3 is skipped which is
the most extensive calculation.
As temperature increases, the alpha-particle contribution gets less and less,
we actually do this uniform matter calculation dependent on $T$.
For $T=15, 20 MeV$, we still use the large RMF input table as described in
step 1 to calculate the mixed gas phase including alpha-particles.
When we write the finial EOS table, we switch off the alpha-particle
contribution in the range of $X_{\alpha}<10^{-4}$ due to accuracy requirement.
For $T=25, 32, 40 MeV$, we use a small RMF input table for the calculation
of the fractions ($X_n$, $X_p$, and $X_{\alpha}$) without alpha-particle
contribution for other quantities, actually the alpha-particle contribution
is really small in this high temperature range.
As for $T=50, 63, 80, 100 MeV$, we only use a small RMF input table to
calculate $X_{\alpha}$, while all other quantities are coming from
the ideal gas approximation (including antiparticles) or the RMF theory.
The connection between ideal gas approximation and RMF theory is also
dependent on $T$, and we try to find out the density
where they can be connected with each other smoothly.
At the end, we should add the pure neutron matter results ($Y_p=0$),
there is only one composition in this case,
we use the RMF theory to do this calculation in high density range,
while adopt ideal gas approximation for low density range,
then combine them to get finial EOS table for $Y_p=0$.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\section{Resulting EOS table}
In this section, we will explain the resulting EOS table and
list the definitions of the quantities tabulated.
We give the resulting EOS by three tables, which are named
(1) eos.tab (main EOS table, size: 54.5MB)
\begin{itemize}
\item temperature $T (MeV)$:
$ -1.0 \leq \log_{10}(T) \leq 2.0; \,\, $
mesh of $\log_{10}(T) \sim 0.1 $
\item proton fraction $Y_p$:
$ -2.00 \leq \log_{10}(Y_p) \leq -0.25; \,\, $
mesh of $\log_{10}(Y_p) \sim 0.025 $
\item baryon mass density $\rho_B (g/cm^3)$:
$ 5.1 \leq \log_{10}(\rho_B) \leq 15.4; \,\, $
mesh of $\log_{10}(\rho_B) \sim 0.1 $
\end{itemize}
(2) eos.t00 (EOS at $T=0$, same range of $Y_p$ and $\rho_B$, size: 1.76MB)
(3) eos.yp0 (EOS at $Y_p=0$, same range of $T$ and $\rho_B$, size: 0.79MB) \\
We putted them in the website after compressed by 'gzip' \\
http://www-server.rcnp.osaka-u.ac.jp/$^{\sim}$shen/table/ \\
In the model used in this calculation, we consider four compositions
which are free neutrons, free protons, alpha-particles, and a single species
of heavy nuclei.
Generally speaking, there are two phases existing in the range
covered by this table,
we call the phase where the heavy nuclei are formed as non-uniform matter
phase, and call the phase without heavy nuclei as uniform matter phase.
The phase transition can be found by the variation
of the heavy nucleus fraction between $X_A=0$ and $X_A \neq 0$.
\\
We write the table in the following order, first fix $T$ which is noted
at the top of each block, second fix $Y_p$, third fix $\rho_B$.
The blocks with different $T$ are divided by line 'ccccccc'.
For each $T$, $Y_p$, and $\rho_B$, we tabulate all the quantities
in one line in the order
\begin{itemize}
\item
(1) logarithm of baryon mass density: $\log_{10}(\rho_B)$ [$g/cm^3$] \\
\item
(2) baryon number density: $n_B$ [$fm^{-3}$] \\
The baryon number density is related to the baryon mass density as
$\rho_B=m_{u} n_B$ with $m_{u}=931.49432 MeV$
being the atomic mass unit taken from Ref.\cite{cons}.
\item
(3) logarithm of proton fraction: $\log_{10}(Y_p)$ \\
The proton fraction $Y_p$ of uniform matter is defined by
\begin{equation}
Y_p=\frac{n_p+2n_{\alpha}}{n_B}
=\frac{n_p+2n_{\alpha}}{n_n+n_p+4n_{\alpha}}
\end{equation}
where $n_p$ is the proton number density, $n_n$ is the neutron number density,
$n_{\alpha}$ is the alpha-particle number density,
and $n_B$ is the baryon number density.
For non-uniform matter case, $Y_p$ is the average proton fraction defined by
\begin{equation}
Y_p=\frac{N_p}{N_B}
\end{equation}
where $N_p$ is the proton number per cell, $N_B$ is the baryon number per cell;
\begin{eqnarray}
N_p &=& \int_{cell} \left[\, n_p\left(r\right) + 2n_{\alpha}\left(r\right)
\,\right] d^3r, \\
N_B &=& \int_{cell} \left[\, n_n\left(r\right) +
n_p\left(r\right) + 4n_{\alpha}\left(r\right)
\,\right] d^3r,
\end{eqnarray}
$n_p(r)$ and $n_n(r)$ are the proton and neutron density distribution function
given in Eq.(43),
$n_{\alpha}(r)$ is the alpha-particle density
distribution function given in Eq.(45).
\item
(4) proton fraction: $Y_p$ \\
\item
(5) free energy per baryon: $F$ [$MeV$] \\
The free energy per baryon is defined as relative to
the free nucleon mass $M=938.0 MeV$ in TM1 parameter set as
\begin{equation}
F=\frac{f}{n_B}-M.
\end{equation}
\item
(6) internal energy per baryon: $E_{int}$ [$MeV$] \\
The internal energy per baryon is defined as relative to
the atomic mass unit $m_{u}$ as
\begin{equation}
E_{int}=\frac{\epsilon}{n_B}-m_{u}.
\end{equation}
\item
(7) entropy per baryon: $S$ [$k_B$] \\
The entropy per baryon is related to the entropy density via
\begin{equation}
S=\frac{s}{n_B}.
\end{equation}
\item
(8) mass number of heavy nucleus: $A$ \\
The mass number of heavy nucleus is defined by
\begin{equation}
A=\int_{0}^{R_A} \left[\, n_n\left(r\right) + n_p\left(r\right)
\,\right] 4\pi r^2 dr,
\end{equation}
where $n_n(r)$ and $n_p(r)$ are the density
distribution functions in the Thomas-Fermi approximation as given in Eq.(43),
$R_{A}$ is the maximum between $R_p$ and $R_n$, which is considered as
the boundary of the heavy nucleus.
\item
(9) charge number of heavy nucleus: $Z$ \\
The charge number of heavy nucleus is defined by
\begin{equation}
Z=\int_{0}^{R_A} \, n_p\left(r\right) \, 4\pi r^2 dr,
\end{equation}
\item
(10) effective mass: $M^{*}$ [$MeV$] \\
The effective mass is obtained in the RMF theory for uniform matter.
In non-uniform matter case, the effective mass is a function of space
due to inhomogeneity of nucleon distribution, so it is meaningless
to list the effective mass for non-uniform matter.
We replace the effective mass $M^{*}$ by the free nucleon mass $M$
in the non-uniform matter phase.
\item
(11) free neutron fraction: $X_n$ \\
The free neutron fraction is given by
\begin{equation}
X_n=(n_n^{out} V^{out})/(n_B V_{cell})
\end{equation}
where $V_{cell}$ is the cell volume as defined in Eq.(42),
$V^{out}=V_{cell}-\frac{4\pi}{3}R_{A}^3$ is the volume outside
heavy nucleus, $n_n^{out}$ is the free neutron number density
outside heavy nucleus, while $n_B$ is the average baryon number density.
\item
(12) free proton fraction: $X_p$ \\
The free proton fraction is given by
\begin{equation}
X_p=(n_p^{out} V^{out})/(n_B V_{cell})
\end{equation}
with $n_p^{out}$ is the free proton number density
outside heavy nucleus.
\item
(13) alpha-particle fraction: $X_{\alpha}$ \\
The alpha-particle fraction is defined by
\begin{equation}
X_{\alpha}=4N_{\alpha}/(n_B V_{cell})
\end{equation}
where $N_{\alpha}$ is the alpha-particle number per cell obtained by
\begin{equation}
N_{\alpha}=\int_{cell} n_{\alpha}\left(r\right) d^3r,
\end{equation}
$n_{\alpha}(r)$ is the alpha-particle density
distribution function given in Eq.(45).
\item
(14) heavy nucleus fraction: $X_A$ \\
The heavy nucleus fraction is defined by
\begin{equation}
X_A=A/(n_B V_{cell})
\end{equation}
where $A$ is the mass number of heavy nucleus as defined in Eq.(69).
\item
(15) pressure: $p$ [$MeV/fm^3$] \\
The pressure is calculated through the following
thermodynamical relation
\begin{equation}
p =\left[\,n_B^2 (\partial F/\partial n_B) \,\right]_{T,Y_p}
\end{equation}
\item
(16) chemical potential of neutron: $\mu_n$ [$MeV$] \\
The chemical potential of neutron relative to
the free nucleon mass $M$ is calculated through
the thermodynamical relation
\begin{equation}
\mu_n=\left[\,\partial (n_B F) /\partial n_n \,\right]_{T,n_p}
\end{equation}
here $n_n=(1-Y_p)\,n_B$.
\item
(17) chemical potential of proton: $\mu_p$ [$MeV$] \\
The chemical potential of proton relative to
the free nucleon mass $M$ is calculated through
the thermodynamical relation
\begin{equation}
\mu_p=\left[\,\partial (n_B F) /\partial n_p \,\right]_{T,n_n}
\end{equation}
here $n_p=Y_p\,n_B$.
\end{itemize}
We have done the following check for the EOS table,\\
(1) the consistency of the thermodynamical quantities,
\begin{equation}
F + \frac{p}{n_B}=\mu_n (1-Y_p)+\mu_p Y_p,
\end{equation}
(2) the consistency of the fractions,
\begin{equation}
X_n+X_p+X_{\alpha}+X_A=1,
\end{equation}
(3) the consistency of the relation between $F$, $E_{int}$, and $S$,
\begin{equation}
F=E_{int}-TS +m_{u} -M.
\end{equation}
Generally these consistency relations can be satisfied
within the accuracy of 0.001.
Physical constants to convert units are taken from Ref.\cite{cons}.
\\
We also put the minimization resulting parameters in one table
named par.tab (size: 33.5MB) in the website \\
http://www-server.rcnp.osaka-u.ac.jp/$^{\sim}$shen/parameter/ \\
We write the parameters in the following order,
(1)$dnn=n_n^{in}-n_n^{out}$, (2)$t_n$, (3)$R_n$,
(4)$dnp=n_p^{in}-n_p^{out}$, (5)$t_p$, (6)$R_p$,
(7)$a$, (8)$n_{\alpha}^{out}$, (9)$n_n^{out}$, (10)$n_p^{out}$.
This file was written in the same range of $Y_p$ and $\rho_B$,
but $T$ is below $15 MeV$ where the non-uniform matter phase exists.
The parameters,
(1)$dnn=n_n^{in}-n_n^{out}$, (2)$t_n$, (3)$R_n$,
(4)$dnp=n_p^{in}-n_p^{out}$, (5)$t_p$, (6)$R_p$,
(7)$a$, are written as zero in uniform matter phase.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\section{Suggestions for using the EOS table}
This relativistic EOS table is designed for the use
of supernova simulations, so we do the calculation at each
$T$, $Y_p$, and $\rho_B$. We note here that the meshes of
$\log_{10}(\rho_B)$ and $\log_{10}(T)$ are not exactly equal meshes
in the EOS table because of the difficulty in numerical calculation.
In order to perform the minimization of the free
energy for non-uniform matter, we have to parametrize
the nucleon distributions, so that some quantities
in the EOS table like $A$, $Z$, $X_n$, $X_p$, $X_{\alpha}$, and $X_A$ are
dependent on this parametrization method.
One should keep in mind their definitions when these quantities
are used. We suggest users pay more attention to the
thermodynamical quantities like $F$, $E_{int}$, $S$, $P$, $\mu_n$, and $\mu_p$,
which are supposed not to be sensitive to the parametrization method.
The thermodynamical quantities are found more smooth in
the resulting EOS, while the use of many codes in table 2 may bring some
fluctuations for the fractions, especially in the range of
$2MeV