eos_had_base.h
Go to the documentation of this file.
1 /*
2  -------------------------------------------------------------------
3 
4  Copyright (C) 2006-2017, Andrew W. Steiner
5 
6  This file is part of O2scl.
7 
8  O2scl is free software; you can redistribute it and/or modify
9  it under the terms of the GNU General Public License as published by
10  the Free Software Foundation; either version 3 of the License, or
11  (at your option) any later version.
12 
13  O2scl is distributed in the hope that it will be useful,
14  but WITHOUT ANY WARRANTY; without even the implied warranty of
15  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16  GNU General Public License for more details.
17 
18  You should have received a copy of the GNU General Public License
19  along with O2scl. If not, see <http://www.gnu.org/licenses/>.
20 
21  -------------------------------------------------------------------
22 */
23 /** \file eos_had_base.h
24  \brief File defining \ref o2scl::eos_had_base
25 */
26 #ifndef O2SCL_HADRONIC_EOS_H
27 #define O2SCL_HADRONIC_EOS_H
28 
29 #include <iostream>
30 #include <string>
31 
32 #include <boost/numeric/ublas/vector.hpp>
33 
34 #include <o2scl/deriv_gsl.h>
35 #include <o2scl/mroot.h>
36 #include <o2scl/mroot_hybrids.h>
37 #include <o2scl/mm_funct.h>
38 #include <o2scl/eos_base.h>
39 #include <o2scl/fermion_eff.h>
40 #include <o2scl/part.h>
41 #include <o2scl/lib_settings.h>
42 
43 #ifndef DOXYGEN_NO_O2NS
44 namespace o2scl {
45 #endif
46 
47  /** \brief Hadronic equation of state [abstract base]
48 
49  \note This class and all of its children expect the neutron and
50  proton properties to be specified in powers of \f$ \mathrm{fm}
51  \f$ so that masses and chemical potentials are in powers of \f$
52  \mathrm{fm}^{-1} \f$ and energy densities and pressures are in
53  powers of \f$ \mathrm{fm}^{-4} \f$ .
54 
55  Denote the number density of neutrons as \f$ n_n \f$, the number
56  density of protons as \f$ n_p \f$, the total baryon density \f$
57  n_B = n_n + n_p \f$, the asymmetry \f$ \delta \equiv
58  (n_n-n_p)/n_B \f$, the nuclear saturation density as \f$ n_0
59  \approx 0.16~\mathrm{fm}^{-3} \f$, and the quantity \f$ \epsilon
60  \equiv (n_B-n_0)/3n_0 \f$. (Note that some authors define
61  \f$ \delta \f$ as \f$ n_n - n_p \f$, which is not the same as
62  the definition above.) Then the energy per baryon of nucleonic
63  matter can be written as an expansion around
64  \f$ \epsilon =\delta = 0 \f$
65  \f[
66  E(n_B,\delta) = -B + \frac{\tilde{K}}{2!} {\epsilon}^2 +
67  \frac{Q}{3!} {\epsilon}^3 + \delta^2 \left(S + L \epsilon +
68  \frac{K_{\mathrm{sym}}}{2!} {\epsilon}^2
69  + \frac{Q_{\mathrm{sym}}}{3!} \epsilon^3 \right) +
70  E_4(n_B,\delta) + {\cal O}(\delta^6)
71  \qquad \left(\mathrm{Eq.}~1\right)
72  \f]
73  where \f$ E_4 \f$ represents the quartic terms
74  \f[
75  E_4(n_B,\delta)
76  = \delta^4 \left(S_4 + L_4 \epsilon + \frac{K_4}{2!} {\epsilon}^2
77  + \frac{Q_4}{3!} \epsilon^3 \right) \qquad
78  \left(\mathrm{Eq.}~2\right)
79  \f]
80  (Adapted slightly from \ref Piekarewicz09). From this, one can
81  compute the energy density of nuclear matter \f$
82  \varepsilon(n_B,\delta) = n_B E(n_B,\delta) \f$, the chemical
83  potentials \f$ \mu_i \equiv (\partial \varepsilon) / (\partial
84  n_i ) \f$ and the pressure \f$ P = -\varepsilon + \mu_n n_n +
85  \mu_p n_p \f$. This expansion motivates the definition of
86  several separate terms. The binding energy \f$ B \f$ of
87  symmetric nuclear matter (\f$ \delta = 0 \f$) is around 16 MeV.
88 
89  The compression modulus is usually defined by \f$ \chi = -1/V
90  (dV/dP) = 1/n (dP/dn)^{-1} \f$ . In nuclear physics it has
91  become common to use the incompressibility (or bulk) modulus
92  with an extra factor of 9 (motivated by the 3 in the denominator
93  in the definition of \f$ \epsilon \f$), \f$ K=9/(n \chi) \f$ and
94  refer to \f$ K \f$ simply as the incompressibility. Here, we
95  define the function
96  \f[
97  K(n_B,\delta) \equiv 9 \left( \frac{\partial P}{\partial n_B}
98  \right) = 9 n_B \left(\frac{\partial^2 \varepsilon}
99  {\partial n_B^2}\right)
100  \f]
101  This quantity is computed by the function \ref fcomp() by
102  computing the first derivative of the pressure, which is more
103  numerically stable than the second derivative of the energy
104  density (since most \o2 EOSs compute the pressure exactly). This
105  function is typically evaluated at the point \f$
106  (n_B=n_0,\delta=0) \f$ and is stored in \ref comp by the
107  function \ref saturation(). This quantity is not always the same
108  as \f$ \tilde{K} \f$, defined here as
109  \f[
110  \tilde{K}(n_B,\delta) = 9 n_B^2 \left(\frac{\partial^2 E}{\partial
111  n_B^2}\right) = K(n_B,\delta) - \frac{1}{n_B} 18 P(n_B,\delta)
112  \f]
113  We denote \f$ K \equiv K(n_B=n_0,\delta=0) \f$ and similarly for
114  \f$ \tilde{K} \f$, the quantity in Eq. 1 above. In nuclear
115  matter at saturation, the pressure is zero and \f$ K = \tilde{K}
116  \f$. See \ref Chabanat97 for further discussion of the distinction
117  between \f$ K \f$ and \f$ \tilde{K} \f$.
118 
119  The symmetry energy, \f$ S(n_B,\delta), \f$ can be defined as
120  \f[
121  S(n_B,\delta) \equiv \frac{1}{2 n_B}\frac{\partial^2 \varepsilon}
122  {\partial \delta^2}
123  \f]
124  and the parameter \f$ S \f$ in Eq. 1 is just \f$ S(n_0,0) \f$.
125  Using
126  \f[
127  \left(\frac{\partial \varepsilon}{\partial \delta}\right)_{n_B} =
128  \frac{\partial \varepsilon}{\partial n_n}
129  \left(\frac{\partial n_n}{\partial \delta}\right)_{n_B} +
130  \frac{\partial \varepsilon}{\partial n_p}
131  \left(\frac{\partial n_p}{\partial \delta}\right)_{n_B}
132  = \frac{n_B}{2} \left(\mu_n - \mu_p \right)
133  \f]
134  this can be rewritten
135  \f[
136  S(n_B,\delta) = \frac{1}{4} \frac{\partial}{\partial \delta}
137  \left(\mu_n - \mu_p\right)
138  \f]
139  where the dependence of the chemical potentials on \f$ n_B \f$
140  and \f$ \delta \f$ is not written explicitly. This quantity is
141  computed by function \ref fesym(). Note that many of the
142  functions in this class are written in terms of the proton
143  fraction \f$ x_p = (1-\delta)/2 \f$ denoted as <tt>'pf'</tt>
144  instead of as functions of \f$ \delta \f$. Frequently, \f$
145  S(n_B,\delta) \f$ is evaluated at \f$ \delta=0 \f$ to give a
146  univariate function of the baryon density. It is sometimes also
147  evaluated at the point \f$ (n_B=n_0, \delta=0) \f$, and this
148  value is denoted by \f$ S \f$ above and is typically stored in
149  \ref esym.
150  Alternatively, one can define the symmetry energy by
151  \f[
152  \tilde{S}(n_B) = E(n_B,\delta=1)-E(n_B,\delta=0)
153  \f]
154  which is computed by function \ref fesym_diff() . The
155  functions \f$ S(n_B,\delta=0) \f$ and \f$ \tilde{S}(n_B) \f$
156  are equal when \f$ {\cal O}(\delta^4) \f$ terms are zero.
157  In this case, \f$ \mu_n - \mu_p \f$ is proportional to
158  \f$ \delta \f$ and so
159  \f[
160  S(n_B) = \tilde{S}(n_B) = \frac{1}{4}
161  \frac{(\mu_n-\mu_p)}{\delta} \, .
162  \f]
163 
164  These functions can also be generalized to finite temperature
165  \f[
166  S(n_B,\delta,T) = \frac{1}{4} \frac{\partial}{\partial \delta}
167  \left[\mu_n(n_B,\delta,T) - \mu_p(n_B,\delta,T)\right] \, ,
168  \f]
169  and
170  \f[
171  \tilde{S}(n_B,T) = F(n_B,\delta=1,T)-F(n_B,\delta=0,T) \, .
172  \f]
173 
174  The symmetry energy slope parameter \f$ L \f$, can be defined
175  by
176  \f[
177  L(n_B,\delta) \equiv 3 n_B \frac{\partial S(n_B,\delta)}
178  {\partial n_B} = 3 n_B \frac{\partial}{\partial n_B} \left[
179  \frac{1}{2 n_B} \frac{\partial^2 \varepsilon}{\partial \delta^2}
180  \right]
181  \f]
182  This can be rewritten as
183  \f[
184  L(n_B,\delta) = \frac{3 n_B}{4} \frac{\partial}{\partial n_B}
185  \frac{\partial}{\partial \delta} \left(\mu_n - \mu_p\right)
186  \f]
187  (where the derivatives can be evaluated in either order)
188  and this is the method used to compute this function
189  in \ref fesym_slope(). Alternatively, using
190  \f[
191  \left(\frac{\partial \varepsilon}{\partial n_B}\right)_{\delta} =
192  \frac{\partial \varepsilon}{\partial n_n}
193  \left(\frac{\partial n_n}{\partial n_B}\right)_{\delta} +
194  \frac{\partial \varepsilon}{\partial n_p}
195  \left(\frac{\partial n_p}{\partial n_B}\right)_{\delta}
196  = \frac{1}{2} \left(\mu_n + \mu_p \right)
197  \f]
198  \f$ L \f$ can be rewritten
199  \f{eqnarray*}
200  L(n_B,\delta) &=& 3 n_B \left[\frac{-1}{2 n_B^2}
201  \frac{\partial^2 \varepsilon}{\partial \delta^2} +
202  \frac{1}{4 n_B} \frac{\partial^2}{\partial \delta^2}
203  \left(\mu_n + \mu_p\right)\right] \\
204  &=& \frac{3}{4}\frac{\partial^2}{\partial \delta^2}
205  \left(\mu_n + \mu_p\right) - 3 S(n_B,\delta) \, .
206  \f}
207 
208  The third derivative with respect to the baryon density is
209  sometimes called the skewness. Here, we define
210  \f[
211  Q(n_B,\delta) = 27 n_B^3 \frac{\partial^3}{\partial n_B^3}
212  \left(\frac{\varepsilon}{n_B}\right) =
213  27 n_B^3 \frac{\partial^2}{\partial n_B^2}
214  \left(\frac{P}{n_B^2}\right)
215  \f]
216  and this function is computed in \ref fkprime() by computing
217  the second derivative of the pressure.
218 
219  The second derivative of the symmetry energy with respect
220  to the baryon density is
221  \f[
222  K_{\mathrm{sym}}(n_B,\delta) = 9 n_B^2
223  \frac{\partial^2}{\partial n_B^2} S(n_B,\delta)
224  \f]
225  and this function is computed in \ref fesym_curve().
226 
227  The third derivative of the symmetry energy with respect
228  to the baryon density is
229  \f[
230  Q_{\mathrm{sym}}(n_B,\delta) = 27 n_B^3
231  \frac{\partial^3}{\partial n_B^3} S(n_B,\delta)
232  \f]
233  and this function is computed in \ref fesym_skew(). Note that
234  the numerical evaluation of higher derivatives can make \ref
235  eos_had_base::fesym_curve() and \ref eos_had_base::fesym_skew()
236  inaccurate.
237 
238  Note that assuming terms of order \f$ \epsilon^3 \f$ and higher
239  are zero and solving for the baryon density for which \f$ P=0 \f$
240  gives, to order \f$ \delta^2 \f$ (\ref Piekarewicz09),
241  \f[
242  n_{B,\mathrm{sat}} = n_0 \left[ 1 - \frac{3 L \delta^2}{K} \right]
243  \f]
244  this implies a new 'incompressibility' around the saturation
245  point, i.e.
246  \f[
247  K(n_B=n_{B,\mathrm{sat}},\delta)=
248  K+\delta^2 \left( K_{\mathrm{sym}}-6 L- \frac{L Q}{K} \right)
249  + {\cal O}\left(\delta^4\right)
250  \f]
251  The quantity in parenthesis is referred to by some authors as
252  \f$ K_{\tau} \f$. Note that, because one is evaluating this at
253  \f$ n_B=n_{B,\mathrm{sat}} \f$, this is distinct from
254  \f[
255  \tilde{K}_{\tau} \equiv \frac{1}{2} \frac{\partial^2 K(n_B,\delta)}
256  {\partial \delta^2}
257  \f]
258  which is equal to \f$ K_{\mathrm{sym}} + 6 L \f$ to lowest order
259  in \f$ \delta \f$ at \f$ n_B = n_0 \f$.
260 
261  The quartic symmetry energy \f$ S_4(n_B,\delta) \f$ can be defined as
262  \f[
263  S_4(n_B,\delta) \equiv \frac{1}{24 n_B}\frac{\partial^4 \varepsilon}
264  {\partial \delta^4}
265  \f]
266  However, fourth derivatives are difficult numerically, and so an
267  alternative quantity is preferable. Instead, one can evaluate
268  the extent to which \f$ {\cal O}(\delta^4) \f$ terms are
269  important from
270  \f[
271  \eta(n_B) \equiv \frac{E(n_B,1)-E(n_B,1/2)}
272  {3 \left[E(n_B,1/2)-E(n_B,0)\right]}
273  \f]
274  as described in \ref Steiner06 . This function can be expressed
275  either in terms of \f$ \tilde{S} \f$ or \f$ S_4 \f$
276  \f[
277  \eta(n_B) = \frac{5 \tilde{S}(n_B) - S(n_B,0)}
278  {\tilde{S}(n_B) + 3 S(n_B,0)} =
279  \frac{5 S_4(n_B,0) + 4 S(n_B,0)}
280  {S_4(n_B,0) + 4 S(n_B,0)}
281  \f]
282  Alternatively, \f$ S_4 \f$ can be written
283  \f[
284  4 S(n_B) \left[ \frac{1-\eta(n_B)}{\eta(n_B)-5} \right] \, .
285  \f]
286 
287  Evaluating this function at the saturation density gives
288  \f[
289  \eta(n_0) = \frac{4 S + 5 S_4}{4 S + S_4}
290  \f]
291  (Note that \f$ S_4 \f$ is referred to as \f$ Q \f$ in
292  \ref Steiner06). Sometimes it is useful to separate out
293  the kinetic and potential parts of the energy density when
294  computing \f$ \eta \f$, and the class \ref eos_had_sym4_base
295  is useful for this purpose.
296 
297  The function \f$ L_4 \f$ can also be rewritten in
298  \f$ \eta^{\prime} \f$ (now suppressing the dependence
299  on \f$ n_B \f$),
300  \f[
301  \eta^{\prime} = \frac{16 \left( L_4 S - L S_4 \right)}
302  {3 n_B \left(4 S +S_4 \right)^2}
303  \f]
304  then using the expression for \f$ S_4 \f$,
305  \f[
306  \eta^{\prime} = \frac{\left(\eta -5\right)}{48 n_B S }
307  \left[ L_4 \left(\eta -5\right) + 4 L
308  \left(\eta -1\right)\right]
309  \f]
310 
311  \future Replace fmsom() with f_effm_scalar(). This has to wait
312  until f_effm_scalar() has a sensible definition when mn is
313  not equal to mp
314  \future Could write a function to compute the "symmetry free energy"
315  or the "symmetry entropy"
316  \future Compute the speed of sound or the number susceptibilities?
317  \future A lot of the numerical derivatives here might possibly
318  request negative number densities for the nucleons, which
319  may cause exceptions, espescially at very low densities.
320  Since the default EOS objects are GSL derivatives, one can
321  get around these issues by setting the GSL derivative object
322  step size, but this is a temporary solution.
323  */
324  class eos_had_base : public eos_base {
325 
326  public:
327 
329 
330  eos_had_base();
331 
332  virtual ~eos_had_base() {};
333 
334  /** \brief Binding energy (without the rest mass) in
335  \f$ \mathrm{fm}^{-1} \f$
336  */
337  double eoa;
338 
339  /// Compression modulus in \f$ \mathrm{fm}^{-1} \f$
340  double comp;
341 
342  /// Symmetry energy in \f$ \mathrm{fm}^{-1} \f$
343  double esym;
344 
345  /// Saturation density in \f$ \mathrm{fm}^{-3} \f$
346  double n0;
347 
348  /// Effective mass (neutron)
349  double msom;
350 
351  /// Skewness in \f$ \mathrm{fm}^{-1} \f$
352  double kprime;
353 
354  /// \name Equation of state
355  //@{
356  /** \brief Equation of state as a function of the chemical potentials
357  */
358  virtual int calc_p(fermion &n, fermion &p, thermo &th)=0;
359 
360  /** \brief Equation of state as a function of density
361  */
362  virtual int calc_e(fermion &n, fermion &p, thermo &th)=0;
363  //@}
364 
365  /// \name EOS properties
366  //@{
367  /** \brief Calculate the incompressibility in \f$ \mathrm{fm}^{-1} \f$
368  using calc_e()
369 
370  This function computes \f$ K (n_B,\delta) = 9 n_B \partial^2
371  \varepsilon /(\partial n_B^2) = 9 \partial P / (\partial n_B)
372  \f$ . The value of \f$ K(n_0,0) \f$, often referred to as the
373  "compressibility", is stored in \ref comp by \ref saturation()
374  and is about 240 MeV at saturation density.
375  */
376  virtual double fcomp(double nb, double delta=0.0);
377 
378  /** \brief Compute the incompressibility and its uncertainty
379 
380  This function works like \ref fcomp(), except that it also
381  returns the uncertainty in \c unc.
382  */
383  virtual double fcomp_err(double nb, double delta, double &unc);
384 
385  /** \brief Calculate the energy per baryon in \f$ \mathrm{fm}^{-1} \f$
386  using calc_e()
387 
388  This function computes the energy per baryon of matter
389  without the nucleon rest masses at the specified baryon
390  density, \c nb, and isospin asymmetry \c delta.
391  */
392  virtual double feoa(double nb, double delta=0.0);
393 
394  /** \brief Calculate symmetry energy of matter in
395  \f$ \mathrm{fm}^{-1} \f$ using \ref calc_dmu_delta() .
396 
397  This function computes the symmetry energy,
398  \f[
399  \left(\frac{1}{2 n_B}\frac{d^2 \varepsilon}{d \delta^2}
400  \right) = \frac{1}{4} \frac{\partial}{\partial \delta}
401  \left(\mu_n - \mu_p \right)
402  \f]
403  at the value of \f$ n_B \f$ given in \c nb and \f$ \delta \f$
404  given in \c delta. The symmetry energy at \f$ \delta=0 \f$ at
405  the saturation density and is stored in \ref esym by
406  \ref saturation().
407  */
408  virtual double fesym(double nb, double delta=0.0);
409 
410  /** \brief Calculate symmetry energy of matter and its
411  uncertainty in \f$ \mathrm{fm}^{-1} \f$
412 
413  This estimates the uncertainty due to the numerical
414  differentiation, assuming that difference betwen the neutron
415  and proton chemical potentials is computed exactly by \ref
416  calc_dmu_delta() .
417  */
418  virtual double fesym_err(double nb, double delta, double &unc);
419 
420  /** \brief The symmetry energy slope parameter
421  in \f$ \mathrm{fm}^{-1} \f$
422 
423  This returns the value of the "slope parameter" of the
424  symmetry energy as a function of baryon density \c nb and
425  isospin asymmetry \c delta. This ranges between about zero and
426  200 MeV for most equations of state.
427  */
428  virtual double fesym_slope(double nb, double delta=0.0);
429 
430  /** \brief The curvature of the symmetry energy
431  in \f$ \mathrm{fm}^{-1} \f$
432  */
433  virtual double fesym_curve(double nb, double delta=0.0);
434 
435  /** \brief The skewness of the symmetry energy
436  in \f$ \mathrm{fm}^{-1} \f$
437  */
438  virtual double fesym_skew(double nb, double delta=0.0);
439 
440  /** \brief Calculate symmetry energy of matter as energy of
441  neutron matter minus the energy of nuclear matter
442  in \f$ \mathrm{fm}^{-1} \f$
443 
444  This function returns the energy per baryon of neutron matter
445  minus the energy per baryon of nuclear matter. This will
446  deviate significantly from the results from fesym() only if
447  the dependence of the symmetry energy on \f$ \delta \f$ is not
448  quadratic.
449  */
450  virtual double fesym_diff(double nb);
451 
452  /** \brief The strength parameter for quartic terms in the
453  symmetry energy
454  */
455  virtual double feta(double nb);
456 
457  /** \brief The derivative of the strength parameter for quartic
458  terms in the symmetry energy
459  */
460  virtual double feta_prime(double nb);
461 
462  /** \brief Calculate skewness of nuclear matter
463  in \f$ \mathrm{fm}^{-1} \f$ using calc_e()
464 
465  The skewness is defined to be
466  \f$ 27 n^3 d^3 (\varepsilon/n)/(d n^3) =
467  27 n^3 d^2 (P/n^2) / (d n^2) \f$
468  and is denoted 'kprime'. This definition seems to be ambiguous
469  for densities other than the saturation density and is not
470  quite analogous to the compression modulus.
471  */
472  virtual double fkprime(double nb, double delta=0.0);
473 
474  /** \brief Calculate reduced neutron effective mass using calc_e()
475 
476  Neutron effective mass (as stored in <tt>part::ms</tt>)
477  divided by vacuum mass (as stored in <tt>part::m</tt>) in
478  nuclear matter at saturation density. Note that this simply
479  uses the value of n.ms from calc_e(), so that this effective
480  mass could be either the Landau or Dirac mass depending on the
481  context. Note that this may not be equal to the reduced proton
482  effective mass.
483  */
484  virtual double fmsom(double nb, double delta=0.0);
485 
486  /** \brief Neutron (reduced) effective mass
487  */
488  virtual double f_effm_neut(double nb, double delta=0.0);
489 
490  /** \brief Proton (reduced) effective mass
491  */
492  virtual double f_effm_prot(double nb, double delta=0.0);
493 
494  /** \brief Scalar effective mass
495 
496  Given the reduced nucleon effective masses, \f$ m_n^{*} \f$
497  and \f$ m_p^{*} \f$, the scalar and vector effective masses
498  are defined by (see e.g. \ref Farine01)
499  \f[
500  \frac{1}{m^{*}_n} = (1+\delta) \frac{1}{m^{*}_s} -
501  \delta \frac{1}{m^{*}_v}
502  \f]
503  \f[
504  \frac{1}{m^{*}_p} = (1-\delta) \frac{1}{m^{*}_s} +
505  \delta \frac{1}{m^{*}_v}
506  \f]
507  this implies
508  \f[
509  m_{\mathrm{scalar}}^{*} =
510  \frac{2 m^{*}_n m^{*}_p}{m^{*}_n+m^{*}_p}
511  \f]
512  and
513  \f[
514  m_{\mathrm{vector}}^{*} =
515  \frac{2 m^{*}_n m^{*}_p \delta}{m^{*}_n - m^{*}_p
516  + \delta(m^{*}_n + m^{*}_p)}
517  \f]
518  */
519  virtual double f_effm_scalar(double nb, double delta=0.0);
520 
521  /** \brief Vector effective mass
522 
523  See documentation for \ref eos_had_base::f_effm_scalar().
524 
525  Note that the vector effective mass diverges when \f$ m^{*}_n
526  = m^{*}_p \f$ and \f$ \delta = 0 \f$, but many models have
527  vector effective masses which are independent of \f$ \delta
528  \f$. For now, we set \f$ \delta =1 \f$ to be the default
529  value, corresponding to neutron matter.
530  */
531  virtual double f_effm_vector(double nb, double delta=1.0);
532 
533  /** \brief Calculate saturation density using calc_e()
534 
535  This function finds the baryon density for which the pressure
536  vanishes.
537  */
538  virtual double fn0(double delta, double &leoa);
539 
540  /** \brief Compute the number susceptibilities as a function of
541  the chemical potentials, \f$ \partial^2 P / \partial \mu_i
542  \mu_j \f$
543  */
544  virtual void f_number_suscept(double mun, double mup, double &dPdnn,
545  double &dPdnp, double &dPdpp);
546 
547  /** \brief Compute the 'inverse' number susceptibilities as a
548  function of the densities, \f$ \partial^2 \varepsilon /
549  \partial n_i n_j \f$
550  */
551  virtual void f_inv_number_suscept(double mun, double mup, double &dednn,
552  double &dednp, double &dedpp);
553 
554  /** \brief Calculates some of the EOS properties at the saturation
555  density
556 
557  This computes the saturation density, and the
558  incompressibility, the symmetry energy, the binding energy,
559  the reduced neutron effective mass at the saturation density,
560  and the skewness in isospin-symmetric matter. The results are
561  stored in \ref n0, \ref comp, \ref esym, \ref eoa, \ref msom,
562  and \ref kprime, respectively.
563 
564  \future It would be great to provide numerical uncertainties
565  in the saturation properties.
566  */
567  virtual void saturation();
568  //@}
569 
570  /// \name Functions for calculating physical properties
571  //@{
572  /** \brief Compute the neutron chemical potential at fixed
573  density
574 
575  This function uses \ref neutron, \ref proton, \ref
576  eos_base::eos_thermo, and \ref calc_e() .
577  */
578  double calc_mun_e(double nn, double np);
579 
580  /** \brief Compute the energy density as a function of the nucleon
581  densities
582  */
583  double calc_ed(double nn, double np);
584 
585  /** \brief Compute the pressure as a function of the nucleon
586  chemical potentials
587  */
588  double calc_pr(double nn, double np);
589 
590  /** \brief Compute the proton chemical potential at fixed
591  density
592 
593  This function uses \ref neutron, \ref proton, \ref
594  eos_base::eos_thermo, and \ref calc_e() .
595  */
596  double calc_mup_e(double nn, double np);
597 
598  /** \brief Compute the neutron density at fixed
599  chemical potential
600 
601  This function uses \ref neutron, \ref proton, \ref
602  eos_base::eos_thermo, and \ref calc_e() .
603  */
604  double calc_nn_p(double mun, double mup);
605 
606  /** \brief Compute the proton density at fixed
607  chemical potential
608 
609  This function uses \ref neutron, \ref proton, \ref
610  eos_base::eos_thermo, and \ref calc_e() .
611  */
612  double calc_np_p(double mun, double mup);
613 
614  /** \brief Compute the difference between neutron and proton chemical
615  potentials as a function of the isospin asymmetry
616 
617  This function uses \ref neutron, \ref proton, \ref
618  eos_base::eos_thermo, and \ref calc_e() .
619  */
620  double calc_dmu_delta(double delta, double nb);
621 
622  /** \brief Compute the sum of the neutron and proton chemical
623  potentials as a function of the isospin asymmetry
624 
625  This uses \ref neutron, \ref proton, \ref eos_base::eos_thermo,
626  and \ref calc_e() .
627  */
628  double calc_musum_delta(double delta, double nb);
629 
630  /** \brief Compute the pressure as a function of baryon density
631  at fixed isospin asymmetry
632 
633  Used by fcomp().
634  */
635  double calc_pressure_nb(double nb, double delta=0.0);
636 
637  /** \brief Compute the energy density as a function of baryon density
638  at fixed isospin asymmetry
639 
640  This uses \ref neutron, \ref proton, \ref eos_base::eos_thermo,
641  and \ref calc_e() .
642  */
643  double calc_edensity_nb(double nb, double delta=0.0);
644 
645  /** \brief Compute derivatives at constant proton fraction
646  */
647  void const_pf_derivs(double nb, double pf,
648  double &dednb_pf, double &dPdnb_pf);
649 
650  /** \brief Calculate pressure / baryon density squared in nuclear
651  matter as a function of baryon density at fixed isospin asymmetry
652 
653  Used by fkprime().
654 
655  This uses \ref neutron, \ref proton, \ref eos_base::eos_thermo,
656  and \ref calc_e() .
657  */
658  double calc_press_over_den2(double nb, double delta=0.0);
659 
660  /** \brief Calculate energy density as a function of the
661  isospin asymmetry at fixed baryon density
662 
663  Used by fesym().
664 
665  This function calls \ref eos_had_base::calc_e() with the
666  internally stored neutron and proton objects.
667  */
668  double calc_edensity_delta(double delta, double nb);
669  //@}
670 
671  /// \name Nuclear matter functions
672  //@{
673  /** \brief Solve for the chemical potentials given the densities
674 
675  The neutron and proton chemical potentials should be stored in
676  <tt>x[0]</tt> and <tt>x[1]</tt> and the neutron and proton
677  densities should be stored in <tt>pa[0]</tt> and
678  <tt>pa[1]</tt>.
679 
680  Because this function is designed to be used in a solver,
681  it returns <tt>exc_efailed</tt> without calling the error
682  handler if the densities are not finite.
683 
684  This function is used by \ref eos_had_pres_base::calc_e().
685  */
686  int nuc_matter_p(size_t nv, const ubvector &x, ubvector &y,
687  double nn0, double np0);
688 
689  /** \brief Solve for the densities given the chemical potentials
690 
691  The neutron and proton densities should be stored in
692  <tt>x[0]</tt> and <tt>x[1]</tt> and the neutron and proton
693  chemical potentials should be stored in <tt>pa[0]</tt> and
694  <tt>pa[1]</tt>.
695 
696  Because this function is designed to be used in a solver,
697  it returns <tt>exc_efailed</tt> without calling the error
698  handler if the chemical potentials are not finite.
699 
700  This function is used by \ref eos_had_eden_base::calc_p().
701  */
702  int nuc_matter_e(size_t nv, const ubvector &x, ubvector &y,
703  double mun0, double mup0);
704  //@}
705 
706  /// \name Set auxilliary objects
707  //@{
708  /** \brief Set class mroot object for use in calculating chemical
709  potentials from densities
710 
711  \note While in principle this allows one to use any \ref mroot
712  object, in practice some of the current EOSs require \ref
713  mroot_hybrids because it automatically avoids regions
714  where the equations are undefined.
715  */
716  virtual void set_mroot(mroot<> &mr);
717 
718  /** \brief Set class mroot object for use calculating saturation density
719 
720  \note While in principle this allows one to use any \ref mroot
721  object, in practice some of the current EOSs require \ref
722  mroot_hybrids because it automatically avoids regions
723  where the equations are undefined.
724  */
725  virtual void set_sat_root(root<> &mr);
726 
727  /// Set \ref deriv_base object to use to find saturation properties
728  virtual void set_sat_deriv(deriv_base<> &de);
729 
730  /** \brief Set the second \ref deriv_base object to use to find
731  saturation properties
732 
733  Computing the slope of the symmetry energy at the saturation
734  density requires two derivative objects, because it has to
735  take an isospin derivative and a density derivative. Thus this
736  second \ref deriv_base object is used in the function
737  fesym_slope().
738  */
739  virtual void set_sat_deriv2(deriv_base<> &de);
740 
741  /// Set neutron and proton
742  virtual void set_n_and_p(fermion &n, fermion &p);
743  //@}
744 
745  /** \brief The defaut neutron
746 
747  By default this has a spin degeneracy of 2 and a mass of \ref
748  o2scl_mks::mass_neutron . Also the value of
749  <tt>part::non_interacting</tt> is set to <tt>false</tt>.
750  */
752 
753  /** \brief The defaut proton
754 
755  By default this has a spin degeneracy of 2 and a mass of \ref
756  o2scl_mks::mass_proton . Also the value of
757  <tt>part::non_interacting</tt> is set to <tt>false</tt>.
758  */
760 
761  /// \name Default solvers and derivative classes
762  //@{
763  /** \brief The default object for derivatives
764 
765  The value of deriv_gsl::h is set to \f$ 10^{-3} \f$ in
766  the eos_had_base constructor.
767  */
769 
770  /** \brief The second default object for derivatives
771 
772  The value of deriv_gsl::h is set to \f$ 10^{-3} \f$ in
773  the eos_had_base constructor.
774  */
776 
777  /** \brief The default solver
778 
779  Used by calc_e() to solve nuc_matter_p() (2 variables) and by
780  calc_p() to solve nuc_matter_e() (2 variables).
781  */
783 
784  /** \brief The default solver for calculating the saturation
785  density
786 
787  Used by fn0() (which is called by saturation()) to solve
788  saturation_matter_e() (1 variable).
789  */
791  //@}
792 
793  /// \name Other functions
794  //@{
795  /** \brief Calculate coefficients for \gradient \part of Hamiltonian
796 
797  \note This is still somewhat experimental.
798 
799  We want the \gradient \part of the Hamiltonian in the form
800  \f[
801  {\cal H}_{\mathrm{grad}} = \frac{1}{2} \sum_{i=n,p}
802  \sum_{j=n,p} Q_{ij}
803  \vec{\nabla} n_i \cdot \vec{\nabla} n_j
804  \f]
805 
806  The expression for the \gradient terms from \ref Pethick95 is
807  \f{eqnarray*}
808  {\cal H}_{\mathrm{grad}}&=&-\frac{1}{4}
809  \left(2 P_1 + P_{1;f}-P_{2;f}\right)
810  \nonumber \\
811  && +\frac{1}{2} \left( Q_1+Q_2 \right)
812  \left(n_n \nabla^2 n_n+n_p \nabla^2 n_p\right) \nonumber \\
813  && + \frac{1}{4}\left( Q_1-Q_2 \right)
814  \left[\left(\nabla n_n\right)^2+\left(\nabla n_p\right)^2
815  \right] \nonumber \\
816  && + \frac{1}{2} \frac{d Q_2}{d n}
817  \left( n_n \nabla n_n + n_p \nabla n_p \right) \cdot \nabla n
818  \f}
819  This can be rewritten
820  \f{eqnarray*}
821  {\cal H}_{\mathrm{grad}}&=&\frac{1}{2}\left(\nabla n\right)^2
822  \left[ \frac{3}{2} P_1+n \frac{d P_1}{d n}\right] \nonumber \\
823  && - \frac{3}{4} \left[ \left( \nabla n_n\right)^2 +
824  \left( \nabla n_p \right)^2 \right] \nonumber \\
825  && -\frac{1}{2} \left[ \right] \cdot \nabla n \frac{d Q_1}{d n}
826  \nonumber \\ && - \frac{1}{4} \left( \nabla n\right)^2 P_2
827  - \frac{1}{4} \left[ \left( \nabla n_n\right)^2 +
828  \left( \nabla n_p \right)^2 \right] Q_2
829  \f}
830  or
831  \f{eqnarray*}
832  {\cal H}_{\mathrm{grad}}&=&\frac{1}{4} \left( \nabla n\right)^2
833  \left[3 P_1 + 2 n \frac{d P_1}{d n}-P_2\right] \nonumber \\
834  && - \frac{1}{4} \left( 3 Q_1+Q_2 \right)
835  \left[ \left( \nabla n_n\right)^2 +
836  \left( \nabla n_p \right)^2 \right] \nonumber \\
837  && - \frac{1}{2} \frac{d Q_1}{d n}
838  \left[ n_n \nabla n_n + n_p \nabla n_p \right]
839  \cdot \nabla n
840  \f}
841  or
842  \f{eqnarray*}
843  {\cal H}_{\mathrm{grad}}&=&\frac{1}{4} \left( \nabla n\right)^2
844  \left[3 P_1 + 2 n \frac{d P_1}{d n}-P_2\right] \nonumber \\
845  && - \frac{1}{4} \left( 3 Q_1+Q_2 \right)
846  \left[ \left( \nabla n_n\right)^2 +
847  \left( \nabla n_p \right)^2 \right] \nonumber \\
848  && - \frac{1}{2} \frac{d Q_1}{d n}
849  \left[ n_n \left( \nabla n_n \right)^2 +
850  n_p \left( \nabla n_p \right)^2 + n \nabla n_n \cdot
851  \nabla n_p \right]
852  \f}
853 
854  Generally, for Skyrme-like interactions
855  \f{eqnarray*}
856  P_i &=& \frac{1}{4} t_i \left(1+\frac{1}{2} x_i \right) \nonumber \\
857  Q_i &=& \frac{1}{4} t_i \left(\frac{1}{2} + x_i \right) \, .
858  \f}
859  for \f$ i=1,2 \f$ .
860 
861  This function uses the assumption \f$ x_1=x_2=0 \f$ to
862  calculate \f$ t_1 \f$ and \f$ t_2 \f$ from the neutron
863  and proton effective masses assuming the Skyrme form. The
864  values of \f$ Q_{ij} \f$ and their derivatives are then computed.
865 
866  The functions set_n_and_p() and set_thermo() will be called by
867  gradient_qij(), to facilitate the use of the \c n, \c p, and
868  \c th parameters.
869 
870  */
871  void gradient_qij(fermion &n, fermion &p, thermo &th,
872  double &qnn, double &qnp, double &qpp,
873  double &dqnndnn, double &dqnndnp,
874  double &dqnpdnn, double &dqnpdnp,
875  double &dqppdnn, double &dqppdnp);
876 
877  /// Return string denoting type ("eos_had_base")
878  virtual const char *type() { return "eos_had_base"; }
879  //@}
880 
881  /// \name Consistency checks
882  //@{
883  /** \brief Check the chemical potentials by computing the
884  derivatives numerically
885  */
886  void check_mu(fermion &n, fermion &p, thermo &th,
887  double &mun_deriv,
888  double &mup_deriv,
889  double &mun_err, double &mup_err);
890 
891  /** \brief Check the densities by computing the
892  derivatives numerically
893  */
894  void check_den(fermion &n, fermion &p, thermo &th,
895  double &nn_deriv, double &np_deriv,
896  double &nn_err, double &np_err);
897  //@}
898 
899 #ifndef DOXYGEN_INTERNAL
900 
901  protected:
902 
903  /// Compute t1 for gradient_qij().
904  double t1_fun(double barn);
905 
906  /// Compute t2 for gradient_qij().
907  double t2_fun(double barn);
908 
909  /// The EOS solver
911 
912  /// The solver to compute saturation properties
914 
915  /// The derivative object for saturation properties
917 
918  /// The second derivative object for saturation properties
920 
921  /// The neutron object
923 
924  /// The proton object
926 
927 #endif
928 
929  };
930 
931  /// A hadronic EOS based on a function of the densities [abstract base]
933  public:
934 
935  /** \brief Equation of state as a function of density
936  */
937  virtual int calc_e(fermion &n, fermion &p, thermo &th)=0;
938 
939  /** \brief Equation of state as a function of the chemical potentials
940  */
941  virtual int calc_p(fermion &n, fermion &p, thermo &th);
942 
943  };
944 
945  /** \brief A hadronic EOS based on a function of the chemical
946  potentials [abstract base]
947  */
949  public:
950 
951  /** \brief Equation of state as a function of the chemical potentials
952  */
953  virtual int calc_p(fermion &n, fermion &p, thermo &th)=0;
954 
955  /** \brief Equation of state as a function of density
956  */
957  virtual int calc_e(fermion &n, fermion &p, thermo &th);
958 
959  };
960 
961  /// A finite temperature hadronic EOS [abstract base]
963 
964 #ifndef DOXYGEN_INTERNAL
965 
966  protected:
967 
968  /// Fermion thermodynamics (default is \ref def_fet)
970 
971  /// Solve for nuclear matter at finite temperature given density
972  int nuc_matter_temp_e(size_t nv, const ubvector &x,
973  ubvector &y, double nn0, double np0, double T);
974 
975  /// Solve for nuclear matter at finite temperature given mu
976  int nuc_matter_temp_p(size_t nv, const ubvector &x,
977  ubvector &y, double mun0, double mup0, double T);
978 
979  /** \brief Solve for the liquid gas phase transition as
980  a function of the densities
981  */
982  int liqgas_dens_solve(size_t nv, const ubvector &x,
983  ubvector &y, fermion &n1, fermion &p1,
984  fermion &n2, fermion &p2, double T,
985  thermo &th1, thermo &th2);
986 
987  /** \brief Solve for the liquid-gas phase transition at fixed
988  baryon density and electron fraction
989  */
990  int liqgas_solve(size_t nv, const ubvector &x,
991  ubvector &y, fermion &n1, fermion &p1,
992  fermion &n2, fermion &p2, double nB0,
993  double Ye0, double T,
994  thermo &th1, thermo &th2);
995 
996  /** \brief Solve for the liquid-gas phase transition in
997  beta-equilibrium
998  */
999  int liqgas_beta_solve(size_t nv, const ubvector &x,
1000  ubvector &y, fermion &n1, fermion &p1,
1001  fermion &n2, fermion &p2,
1002  double nB0, double T,
1003  thermo &th1, thermo &th2, fermion &e);
1004 
1005  /** \brief Compute the entropy
1006  */
1007  double calc_entropy_delta(double delta, double nb, double T);
1008 
1009  /** \brief Compute the difference between the neutron and
1010  proton chemical potentials
1011  */
1012  double calc_dmu_delta_T(double delta, double nb, double T);
1013 
1014 #endif
1015 
1016  public:
1017 
1018  eos_had_temp_base() {
1019  fet=&def_fet;
1020  }
1021 
1022  virtual ~eos_had_temp_base() {}
1023 
1024  /// \name Basic usage
1025  //@{
1026  /** \brief Equation of state as a function of density
1027  */
1028  virtual int calc_e(fermion &n, fermion &p, thermo &th)=0;
1029 
1030  /** \brief Equation of state as a function of densities at
1031  finite temperature
1032  */
1033  virtual int calc_temp_e(fermion &n, fermion &p, double T,
1034  thermo &th)=0;
1035 
1036  /** \brief Equation of state as a function of the chemical potentials
1037  */
1038  virtual int calc_p(fermion &n, fermion &p, thermo &th)=0;
1039 
1040  /** \brief Equation of state as a function of the chemical potentials
1041  at finite temperature
1042  */
1043  virtual int calc_temp_p(fermion &n, fermion &p, double T,
1044  thermo &th)=0;
1045  //@}
1046 
1047  /// Computing finite-temperature integrals
1048  //@{
1049  /** \brief Set the object for computing finite-temperature fermions
1050  (default is \ref def_fet)
1051  */
1053  fet=&f;
1054  }
1055 
1056  /// Default fermion thermodynamics object
1058  //@}
1059 
1060  /// \name Liquid-gas transition functions
1061  //@{
1062  /** \brief Compute liquid-gas phase transition densities using
1063  \ref eos_had_temp_base::calc_temp_e() .
1064 
1065  At fixed baryon number density for \c n1, this determines the
1066  baryon number densities for \c p1, \c n2, and \c p2 which give
1067  chemical and mechanical equilibrium at a fixed temperature
1068  \c T. The thermodynamic quantities assuming bulk matter for
1069  each set is stored in \c th1 and \c th2.
1070  */
1071  virtual int calc_liqgas_dens_temp_e
1072  (fermion &n1, fermion &p1, fermion &n2, fermion &p2,
1073  double T, thermo &th1, thermo &th2);
1074 
1075  /** \brief Compute the liquid-gas phase transition using
1076  \ref eos_had_temp_base::calc_temp_e() .
1077 
1078  At fixed baryon number density, \c nB, fixed electron
1079  fraction, \c Ye, and fixed temperature, \c T, this function
1080  determines the baryon number densities for \c n1, \c p1, \c
1081  n2, and \c p2 which give chemical and mechanical equilibrium.
1082  The thermodynamic quantities assuming bulk matter for each set
1083  is stored in \c th1 and \c th2, and the volume fraction of
1084  phase 1 is stored in \c chi.
1085  */
1086  virtual int calc_liqgas_temp_e
1087  (fermion &n1, fermion &p1, fermion &n2, fermion &p2,
1088  double nB, double Ye, double T, thermo &th1, thermo &th2,
1089  double &chi);
1090 
1091  /** \brief Compute the liquid-gas phase transition in
1092  beta-equilibrium using \ref eos_had_temp_base::calc_temp_e() .
1093 
1094  At fixed baryon number density, \c nB, and fixed temperature,
1095  \c T, this function determines the baryon number densities for
1096  \c n1, \c p1, \c n2, and \c p2 which give chemical and
1097  mechanical equilibrium assuming beta-equilibrium with
1098  electrons. The thermodynamic quantities assuming bulk matter
1099  for each set is stored in \c th1 and \c th2, the electron
1100  fraction is stored in \c Ye, and the volume fraction of phase
1101  1 is stored in \c chi.
1102  */
1103  virtual int calc_liqgas_beta_temp_e
1104  (fermion &n1, fermion &p1, fermion &n2, fermion &p2,
1105  double nB, double T, thermo &th1, thermo &th2,
1106  double &Ye, double &chi);
1107  //@}
1108 
1109  /// \name Functions related to the symmetry energy
1110  //@{
1111  /** \brief Compute the symmetry energy at finite temperature
1112  */
1113  virtual double fesym_T(double nb, double T, double delta=0.0);
1114 
1115  /** \brief Compute the symmetry entropy at finite temperature
1116  */
1117  virtual double fsyment_T(double nb, double T, double delta=0.0);
1118  //@}
1119 
1120  /// \name Helper functions
1121  //@{
1122  /// Neutron chemical potential as a function of the densities
1123  virtual double calc_temp_mun_e(double nn, double np, double T);
1124  /// Proton chemical potential as a function of the densities
1125  virtual double calc_temp_mup_e(double nn, double np, double T);
1126  /// Neutron density as a function of the chemical potentials
1127  virtual double calc_temp_nn_p(double mun, double mup, double T);
1128  /// Proton density as a function of the chemical potentials
1129  virtual double calc_temp_np_p(double mun, double mup, double T);
1130 
1131  /** \brief Compute the free energy as a function of the temperature
1132  and the densities
1133  */
1134  double calc_fr(double nn, double np, double T);
1135  //@}
1136 
1137  /// \name Susceptibilities
1138  //@{
1139  /** \brief Compute the number susceptibilities as a function of
1140  the chemical potentials, \f$ \partial^2 P / \partial \mu_i
1141  \mu_j \f$ at a fixed temperature
1142  */
1143  virtual void f_number_suscept_T
1144  (double mun, double mup, double T, double &dPdnn,
1145  double &dPdnp, double &dPdpp);
1146 
1147  /** \brief Compute the 'inverse' number susceptibilities as a
1148  function of the densities, \f$ \partial^2 \varepsilon /
1149  \partial n_i n_j \f$ at a fixed temperature
1150  */
1151  virtual void f_inv_number_suscept_T
1152  (double mun, double mup, double T, double &dednn,
1153  double &dednp, double &dedpp);
1154  //@}
1155 
1156  /// \name Consistency check
1157  //@{
1158  /** \brief Check the entropy by computing the
1159  derivative numerically
1160  */
1161  void check_en(fermion &n, fermion &p, double T, thermo &th,
1162  double &en_deriv, double &en_err);
1163 
1164  /** \brief Check the chemical potentials at finite temperature
1165  by computing the derivative numerically
1166  */
1167  void check_mu_T(fermion &n, fermion &p, double T, thermo &th,
1168  double &mun_deriv, double &mup_deriv,
1169  double &mun_err, double &mup_err);
1170  //@}
1171 
1172  };
1173 
1174  /** \brief A hadronic EOS at finite temperature
1175  based on a function of the densities [abstract base]
1176 
1177  This class provides automatic computation of \ref calc_e() and
1178  \ref calc_temp_e() assuming that \ref calc_p() and \ref
1179  calc_temp_p() are specified.
1180  */
1182  public:
1183 
1184  /** \brief Equation of state as a function of density
1185  */
1186  virtual int calc_e(fermion &n, fermion &p, thermo &th)=0;
1187 
1188  /** \brief Equation of state as a function of densities at
1189  finite temperature
1190  */
1191  virtual int calc_temp_e(fermion &n, fermion &p, double T,
1192  thermo &th)=0;
1193 
1194  /** \brief Equation of state as a function of the chemical potentials
1195  */
1196  virtual int calc_p(fermion &n, fermion &p, thermo &th);
1197 
1198  /** \brief Equation of state as a function of the chemical potentials
1199  at finite temperature
1200  */
1201  virtual int calc_temp_p(fermion &n, fermion &p, double T,
1202  thermo &th);
1203 
1204  };
1205 
1206  /** \brief A hadronic EOS at finite temperature based on a function
1207  of the chemical potentials [abstract base]
1208 
1209  This class provides automatic computation of \ref calc_p() and
1210  \ref calc_temp_p() assuming that \ref calc_e() and \ref
1211  calc_temp_e() are specified.
1212  */
1214  public:
1215 
1216  /** \brief Equation of state as a function of the chemical potentials
1217  */
1218  virtual int calc_p(fermion &n, fermion &p, thermo &th)=0;
1219 
1220  /** \brief Equation of state as a function of the chemical potentials
1221  at finite temperature
1222  */
1223  virtual int calc_temp_p(fermion &n, fermion &p, double T,
1224  thermo &th)=0;
1225 
1226  /** \brief Equation of state as a function of density
1227  */
1228  virtual int calc_e(fermion &n, fermion &p, thermo &th);
1229 
1230  /** \brief Equation of state as a function of densities at
1231  finite temperature
1232  */
1233  virtual int calc_temp_e(fermion &n, fermion &p, double T,
1234  thermo &th);
1235  };
1236 
1237 #ifndef DOXYGEN_NO_O2NS
1238 }
1239 #endif
1240 
1241 #endif
1242 
root_cern def_sat_root
The default solver for calculating the saturation density.
Definition: eos_had_base.h:790
virtual double fesym_skew(double nb, double delta=0.0)
The skewness of the symmetry energy in .
virtual void set_sat_deriv(deriv_base<> &de)
Set deriv_base object to use to find saturation properties.
virtual void set_fermion_eval_thermo(fermion_eval_thermo &f)
Computing finite-temperature integrals.
virtual double feta_prime(double nb)
The derivative of the strength parameter for quartic terms in the symmetry energy.
double calc_ed(double nn, double np)
Compute the energy density as a function of the nucleon densities.
fermion_eval_thermo * fet
Fermion thermodynamics (default is def_fet)
Definition: eos_had_base.h:969
Equation of state base class.
Definition: eos_base.h:39
double comp
Compression modulus in .
Definition: eos_had_base.h:340
double calc_np_p(double mun, double mup)
Compute the proton density at fixed chemical potential.
virtual double fkprime(double nb, double delta=0.0)
Calculate skewness of nuclear matter in using calc_e()
deriv_gsl def_deriv
The default object for derivatives.
Definition: eos_had_base.h:768
const double barn
virtual double f_effm_neut(double nb, double delta=0.0)
Neutron (reduced) effective mass.
double t1_fun(double barn)
Compute t1 for gradient_qij().
void const_pf_derivs(double nb, double pf, double &dednb_pf, double &dPdnb_pf)
Compute derivatives at constant proton fraction.
void check_mu(fermion &n, fermion &p, thermo &th, double &mun_deriv, double &mup_deriv, double &mun_err, double &mup_err)
Check the chemical potentials by computing the derivatives numerically.
virtual double fmsom(double nb, double delta=0.0)
Calculate reduced neutron effective mass using calc_e()
double calc_mun_e(double nn, double np)
Compute the neutron chemical potential at fixed density.
double calc_pressure_nb(double nb, double delta=0.0)
Compute the pressure as a function of baryon density at fixed isospin asymmetry.
mroot * eos_mroot
The EOS solver.
Definition: eos_had_base.h:910
virtual double fcomp(double nb, double delta=0.0)
Calculate the incompressibility in using calc_e()
virtual double f_effm_scalar(double nb, double delta=0.0)
Scalar effective mass.
virtual double fn0(double delta, double &leoa)
Calculate saturation density using calc_e()
double esym
Symmetry energy in .
Definition: eos_had_base.h:343
virtual void f_inv_number_suscept(double mun, double mup, double &dednn, double &dednp, double &dedpp)
Compute the &#39;inverse&#39; number susceptibilities as a function of the densities, .
virtual double fesym_diff(double nb)
Calculate symmetry energy of matter as energy of neutron matter minus the energy of nuclear matter in...
void check_den(fermion &n, fermion &p, thermo &th, double &nn_deriv, double &np_deriv, double &nn_err, double &np_err)
Check the densities by computing the derivatives numerically.
A hadronic EOS at finite temperature based on a function of the densities [abstract base]...
A finite temperature hadronic EOS [abstract base].
Definition: eos_had_base.h:962
virtual const char * type()
Return string denoting type ("eos_had_base")
Definition: eos_had_base.h:878
deriv_base * sat_deriv
The derivative object for saturation properties.
Definition: eos_had_base.h:916
double n0
Saturation density in .
Definition: eos_had_base.h:346
double calc_press_over_den2(double nb, double delta=0.0)
Calculate pressure / baryon density squared in nuclear matter as a function of baryon density at fixe...
double calc_nn_p(double mun, double mup)
Compute the neutron density at fixed chemical potential.
root * sat_root
The solver to compute saturation properties.
Definition: eos_had_base.h:913
virtual int calc_p(fermion &n, fermion &p, thermo &th)=0
Equation of state as a function of the chemical potentials.
virtual void set_sat_root(root<> &mr)
Set class mroot object for use calculating saturation density.
A hadronic EOS at finite temperature based on a function of the chemical potentials [abstract base]...
void gradient_qij(fermion &n, fermion &p, thermo &th, double &qnn, double &qnp, double &qpp, double &dqnndnn, double &dqnndnp, double &dqnpdnn, double &dqnpdnp, double &dqppdnn, double &dqppdnp)
Calculate coefficients for gradient gradient part part of Hamiltonian.
double t2_fun(double barn)
Compute t2 for gradient_qij().
double calc_mup_e(double nn, double np)
Compute the proton chemical potential at fixed density.
fermion * proton
The proton object.
Definition: eos_had_base.h:925
deriv_gsl def_deriv2
The second default object for derivatives.
Definition: eos_had_base.h:775
double calc_edensity_delta(double delta, double nb)
Calculate energy density as a function of the isospin asymmetry at fixed baryon density.
A hadronic EOS based on a function of the densities [abstract base].
Definition: eos_had_base.h:932
virtual double f_effm_vector(double nb, double delta=1.0)
Vector effective mass.
double calc_musum_delta(double delta, double nb)
Compute the sum of the neutron and proton chemical potentials as a function of the isospin asymmetry...
deriv_base * sat_deriv2
The second derivative object for saturation properties.
Definition: eos_had_base.h:919
int nuc_matter_e(size_t nv, const ubvector &x, ubvector &y, double mun0, double mup0)
Solve for the densities given the chemical potentials.
Hadronic equation of state [abstract base].
Definition: eos_had_base.h:324
A hadronic EOS based on a function of the chemical potentials [abstract base].
Definition: eos_had_base.h:948
virtual double fcomp_err(double nb, double delta, double &unc)
Compute the incompressibility and its uncertainty.
double calc_edensity_nb(double nb, double delta=0.0)
Compute the energy density as a function of baryon density at fixed isospin asymmetry.
virtual int calc_e(fermion &n, fermion &p, thermo &th)=0
Equation of state as a function of density.
virtual void f_number_suscept(double mun, double mup, double &dPdnn, double &dPdnp, double &dPdpp)
Compute the number susceptibilities as a function of the chemical potentials, .
fermion def_proton
The defaut proton.
Definition: eos_had_base.h:759
virtual double feta(double nb)
The strength parameter for quartic terms in the symmetry energy.
virtual double fesym(double nb, double delta=0.0)
Calculate symmetry energy of matter in using calc_dmu_delta() .
double kprime
Skewness in .
Definition: eos_had_base.h:352
virtual double fesym_err(double nb, double delta, double &unc)
Calculate symmetry energy of matter and its uncertainty in .
virtual double f_effm_prot(double nb, double delta=0.0)
Proton (reduced) effective mass.
virtual double fesym_slope(double nb, double delta=0.0)
The symmetry energy slope parameter in .
int nuc_matter_p(size_t nv, const ubvector &x, ubvector &y, double nn0, double np0)
Solve for the chemical potentials given the densities.
double calc_pr(double nn, double np)
Compute the pressure as a function of the nucleon chemical potentials.
fermion * neutron
The neutron object.
Definition: eos_had_base.h:922
mroot_hybrids def_mroot
The default solver.
Definition: eos_had_base.h:782
virtual double feoa(double nb, double delta=0.0)
Calculate the energy per baryon in using calc_e()
double msom
Effective mass (neutron)
Definition: eos_had_base.h:349
double calc_dmu_delta(double delta, double nb)
Compute the difference between neutron and proton chemical potentials as a function of the isospin as...
virtual void saturation()
Calculates some of the EOS properties at the saturation density.
double eoa
Binding energy (without the rest mass) in .
Definition: eos_had_base.h:332
fermion_eff def_fet
Default fermion thermodynamics object.
virtual double fesym_curve(double nb, double delta=0.0)
The curvature of the symmetry energy in .
virtual void set_sat_deriv2(deriv_base<> &de)
Set the second deriv_base object to use to find saturation properties.
virtual void set_mroot(mroot<> &mr)
Set class mroot object for use in calculating chemical potentials from densities. ...
virtual void set_n_and_p(fermion &n, fermion &p)
Set neutron and proton.
fermion def_neutron
The defaut neutron.
Definition: eos_had_base.h:751

Documentation generated with Doxygen. Provided under the GNU Free Documentation License (see License Information).