All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
eos_tov.h
Go to the documentation of this file.
1 /*
2  -------------------------------------------------------------------
3 
4  Copyright (C) 2006-2014, 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_tov.h
24  \brief File defining \ref o2scl::eos_tov
25 */
26 #ifndef O2SCL_TOV_EOS_H
27 #define O2SCL_TOV_EOS_H
28 
29 #include <cmath>
30 #include <iostream>
31 #include <fstream>
32 
33 #include <boost/numeric/ublas/vector.hpp>
34 
35 #include <o2scl/constants.h>
36 #include <o2scl/lib_settings.h>
37 #include <o2scl/interp.h>
38 #include <o2scl/table_units.h>
39 #include <o2scl/vector_derint.h>
40 
41 #ifndef DOXYGEN_NO_O2NS
42 namespace o2scl {
43 #endif
44 
45  /** \brief A EOS base class for the TOV solver
46  */
47  class eos_tov {
48 
49  protected:
50 
51  /** \brief Set to true if the baryon density is provided in the
52  EOS (default false)
53  */
55 
56  friend class tov_solve;
57 
58  public:
59 
60  eos_tov() {
61  baryon_column=false;
62  verbose=1;
63  }
64 
65  virtual ~eos_tov() {}
66 
67  /// Control for output (default 1)
68  int verbose;
69 
70  /** \brief Given the pressure, produce the energy and number densities
71 
72  The arguments \c P and \c e should always be in \f$
73  M_{\odot}/\mathrm{km}^3 \f$ . The argument for \c nb should be
74  in \f$ \mathrm{fm}^{-3} \f$ .
75 
76  If \ref baryon_column is false, then \c nb is unmodified.
77  */
78  virtual void get_eden(double P, double &e, double &nb)=0;
79 
80  /** \brief Given the pressure, produce all the remaining quantities
81 
82  The argument \c P should always be in \f$
83  M_{\odot}/\mathrm{km}^3 \f$ .
84  */
85  virtual void get_aux(double P, size_t &np, std::vector<double> &auxp) {
86  np=0;
87  return;
88  }
89 
90  /** \brief Fill a list with strings for the names of the remaining
91  quanities
92  */
93  virtual void get_names_units(size_t &np,
94  std::vector<std::string> &pnames,
95  std::vector<std::string> &punits) {
96  np=0;
97  return;
98  }
99 
100  /** \brief Check that the baryon density is consistent
101  with the EOS
102  */
103  void check_nb(double &avg_abs_dev, double &max_abs_dev) {
104  if (!baryon_column) {
105  O2SCL_ERR2("Variable 'baryon_column' false in",
106  "eos_tov::check_nb().",exc_einval);
107  }
108  std::vector<double> edv, prv, nbv, dedn;
109  for (double pres=0.1;pres<3.0;pres*=1.001) {
110  double eps, nb;
111  get_eden(pres,eps,nb);
112  edv.push_back(eps);
113  prv.push_back(pres);
114  nbv.push_back(nb);
115  }
116  dedn.resize(edv.size());
117  vector_deriv_interp(edv.size(),nbv,edv,dedn,itp_linear);
118  avg_abs_dev=0.0;
119  max_abs_dev=0.0;
120  for(size_t i=0;i<edv.size();i++) {
121  double abs_dev=(fabs(edv[i]+prv[i]-dedn[i]*nbv[i])/
122  fabs(dedn[i]*nbv[i]));
123  if (abs_dev>max_abs_dev) max_abs_dev=abs_dev;
124  avg_abs_dev+=abs_dev;
125  }
126  avg_abs_dev/=((double)edv.size());
127  return;
128  }
129 
130  };
131 
132  /** \brief The Buchdahl EOS for the TOV solver
133 
134  Given the EOS
135  \f[
136  \varepsilon = 12 \sqrt{p_{*} P}- 5 P
137  \f]
138  the TOV equation has an analytical solution
139  \f[
140  R=(1-\beta) \sqrt{\frac{\pi}{288 p_{*} G (1-2 \beta)}}
141  \f]
142  where \f$ \beta = GM/R \f$.
143 
144  The baryon chemical potential is
145  \f[
146  \mu = \mu_1 \left[ \frac{\left(9 p_{*}-P_1\right) \left(3+t_1\right)
147  \left(3-t_2\right)}{\left(9 p_{*}-P\right)\left(3-t_1\right)
148  \left(3+t_2\right)}\right]^{1/4}
149  \f]
150  where \f$ t_1 = \sqrt{P}/\sqrt{p_{*}} \f$ and \f$ t_2 =
151  \sqrt{P_1/p_{*}} \f$ . The baryon density can then be obtained
152  directly from the thermodynamic identity. In the case that one
153  assumes \f$ \mu_1 = m_n \f$ and \f$ P_1 = 0 \f$, the baryon
154  density can be simplified to
155  \f[
156  n m_n = 12 \sqrt{P p_{*}} \left( 1-\frac{1}{3} \sqrt{P/p_{*}}
157  \right)^{3/2}
158  \f]
159  c.f. Eq. 10 in \ref Lattimer01.
160 
161  The central pressure and energy density are
162  \f[
163  P_c = 36 p_{*} \beta^2
164  \f]
165  \f[
166  {\varepsilon}_c = 72 p_{*} \beta (1-5 \beta/2)
167  \f]
168 
169  Physical solutions are obtained only for \f$ P< 25 p_{*}/144 \f$
170  (ensuring that the argument to the square root is positive)
171  and \f$ \beta<1/6 \f$ (ensuring that the EOS is not acausal).
172 
173  Based on \ref Lattimer01 .
174 
175  \future Figure out what to do with the buchfun() function
176  */
177  class eos_tov_buchdahl : public eos_tov {
178 
179  protected:
180 
181  /** \brief The baryon density at \c ed1
182  */
183  double nb1;
184 
185  /** \brief The energy density for which the baryon density is known
186  */
187  double ed1;
188 
189  /** \brief The pressure at \c ed1
190  */
191  double pr1;
192 
193  public:
194 
195  eos_tov_buchdahl() {
196  Pstar=3.2e-5;
197  }
198 
199  virtual ~eos_tov_buchdahl() {}
200 
201  /** \brief The parameter with units of pressure in units of solar
202  masses per km cubed (default value \f$ 3.2 \times 10^{-5} \f$
203  )
204  */
205  double Pstar;
206 
207  /** \brief Set the baryon density
208  */
209  void set_baryon_density(double nb, double ed) {
210  if (nb<0.0 || ed<0.0) {
211  O2SCL_ERR2("Negative densities not supported in ",
212  "eos_tov_polytrope::set_coeff_index().",exc_einval);
213  }
214  baryon_column=true;
215  nb1=nb;
216  ed1=ed;
217  if (36.0*Pstar*Pstar-5.0*Pstar*ed<0.0) {
218  O2SCL_ERR2("Values of 'Pstar' and 'ed' incommensurate in ",
219  "eos_tov_buchdahl::set_baryon_density().",exc_einval);
220  }
221  pr1=0.04*(72.0*Pstar-5.0*ed+
222  12.0*sqrt(36.0*Pstar*Pstar-5.0*Pstar*ed));
223  return;
224  }
225 
226  /** \brief Given the pressure, produce the energy and number densities
227 
228  If the baryon density is not specified, it should be set to
229  zero or \ref baryon_column should be set to false
230  */
231  virtual void get_eden(double P, double &e, double &nb) {
232  e=12.0*sqrt(Pstar*P)-5.0*P;
233  if (baryon_column) {
234  double mu1=(pr1+ed1)/nb1;
235  //double t1=P/sqrt(P*Pstar);
236  //double t2=pr1/sqrt(pr1*Pstar);
237  double t1=sqrt(P/Pstar);
238  double t2=sqrt(pr1/Pstar);
239  double mu=mu1*pow((-pr1+9.0*Pstar)*(3.0+t1)*(3.0-t2)/
240  (-P+9.0*Pstar)/(3.0-t1)/(3.0+t2),0.25);
241  nb=(P+e)/mu;
242  } else {
243  nb=0.0;
244  }
245  return;
246  }
247 
248  /// Given the pressure, produce all the remaining quantities
249  virtual void get_aux(double P, size_t &np, std::vector<double> &auxp) {
250  np=0;
251  return;
252  }
253 
254  /** \brief Fill a list with strings for the names of the remaining
255  quanities
256  */
257  virtual void get_names(size_t &np, std::vector<std::string> &pnames) {
258  np=0;
259  return;
260  }
261 
262  protected:
263 
264  /** \brief Solve to compute profiles
265 
266  After solving the two equations
267  \f{eqnarray*}
268  r^{\prime} &=& r \left(1-\beta+u\right)^{-1}
269  \left(1 - 2 \beta\right) \nonumber \\
270  A^2 &=& 288 \pi p_{*} G \left( 1 - 2 \beta \right)^{-1}
271  \f}
272  for \f$ u = \beta/(A r^{\prime}) \sin A r^{\prime} \f$
273  and \f$ r^{\prime} \f$,
274  one can compute the pressure and energy density
275  profiles
276  \f{eqnarray*}
277  8 \pi P &=& A^2 u^2 \left(1 - 2 \beta \right)
278  \left(1 - \beta + u \right)^{-2}
279  \nonumber \\
280  8 \pi \varepsilon &=& 2 A^2 u \left(1 - 2 \beta\right)
281  \left(1 - \beta - 3 u/2\right) \left(1 - \beta + u \right)^{-2}
282  \nonumber \\
283  \f}
284 
285  */
286  int solve_u_rp_fun(size_t bv, const std::vector<double> &bx,
287  std::vector<double> &by) {
288  double u, rp;
289  u=bx[0];
290  rp=bx[1];
291  //by[0]=rp*(1.0-beta+u)/(1.0-2.0*beta)-buchrad;
292  //by[1]=beta/biga/rp*sin(biga*rp);
293  return 0;
294  }
295 
296  };
297 
298  /** \brief Standard polytropic EOS \f$ P = K \varepsilon^{1+1/n} \f$
299 
300  The quantity \f$ K \f$ must be in units of
301  \f$ \left(M_{\odot}/km^3\right)^{-1/n} \f$ .
302  */
303  class eos_tov_polytrope : public eos_tov {
304 
305  protected:
306 
307  /** \brief The baryon density at \c ed1
308  */
309  double nb1;
310 
311  /** \brief The energy density for which the baryon density is known
312  */
313  double ed1;
314 
315  /** \brief The pressure at \c ed1
316  */
317  double pr1;
318 
319  /** \brief Coefficient (default 1.0)
320  */
321  double K;
322 
323  /// Index (default 3.0)
324  double n;
325 
326  public:
327 
329  K=1.0;
330  n=3.0;
331  }
332 
333  virtual ~eos_tov_polytrope() {}
334 
335  /** \brief Desc
336  */
337  void set_coeff_index(double coeff, double index) {
338  if (coeff<0.0 || index<0.0) {
339  O2SCL_ERR2("Negative coefficients and indices not supported in ",
340  "eos_tov_polytrope::set_coeff_index().",exc_einval);
341  }
342  K=coeff;
343  n=index;
344  if (baryon_column) {
345  pr1=K*pow(ed1,1.0+1.0/n);
346  }
347  }
348 
349  /** \brief Set the baryon density
350  */
351  void set_baryon_density(double nb, double ed) {
352  if (nb<0.0 || ed<0.0) {
353  O2SCL_ERR2("Negative densities not supported in ",
354  "eos_tov_polytrope::set_coeff_index().",exc_einval);
355  }
356  baryon_column=true;
357  nb1=nb;
358  ed1=ed;
359  pr1=K*pow(ed1,1.0+1.0/n);
360  return;
361  }
362 
363  /** \brief Given the pressure, produce the energy and number densities
364  */
365  virtual void get_eden(double P, double &e, double &nb) {
366  e=pow(P/K,n/(1.0+n));
367  if (baryon_column) {
368  nb=nb1*pow(e/ed1,1.0+n)/pow((e+P)/(ed1+pr1),n);
369  } else {
370  nb=0.0;
371  }
372  return;
373  }
374 
375  };
376 
377  /** \brief Linear EOS \f$ P = c_s^2 (\varepsilon-\varepsilon_0) \f$
378 
379  \note Experimental
380  */
381  class eos_tov_linear : public eos_tov {
382 
383  protected:
384 
385  /** \brief The baryon density at \c ed1
386  */
387  double nb1;
388 
389  /** \brief The energy density for which the baryon density is known
390  */
391  double ed1;
392 
393  /** \brief The pressure at \c ed1
394  */
395  double pr1;
396 
397  /** \brief Coefficient (default 1.0)
398  */
399  double cs2;
400 
401  /// Index (default 0.0)
402  double eps0;
403 
404  public:
405 
406  eos_tov_linear() {
407  cs2=1.0;
408  eps0=0.0;
409  }
410 
411  virtual ~eos_tov_linear() {}
412 
413  /** \brief Desc
414  */
415  void set_cs2_eps0(double cs2_, double eps0_) {
416  eps0=eps0_;
417  cs2=cs2_;
418  return;
419  }
420 
421  /** \brief Set the baryon density
422  */
423  void set_baryon_density(double nb, double ed) {
424  if (nb<0.0 || ed<0.0) {
425  O2SCL_ERR2("Negative densities not supported in ",
426  "eos_tov_polytrope::set_coeff_index().",exc_einval);
427  }
428  baryon_column=true;
429  nb1=nb;
430  ed1=ed;
431  pr1=cs2*(ed-eps0);
432  return;
433  }
434 
435  /** \brief Given the pressure, produce the energy and number densities
436  */
437  virtual void get_eden(double P, double &e, double &nb) {
438  e=P/cs2+eps0;
439  if (baryon_column) {
440  nb=nb1*pow(e+cs2*e-cs2*eps0,1.0/(1.0+cs2))*
441  pow(ed1+cs2*(-eps0+ed1),-1.0/(1.0+cs2));
442  } else {
443  nb=0.0;
444  }
445  return;
446  }
447 
448  /// Given the pressure, produce all the remaining quantities
449  virtual void get_aux(double P, size_t &np, std::vector<double> &auxp) {
450  np=0;
451  return;
452  }
453 
454  /** \brief Fill a list with strings for the names of the remaining
455  quanities
456  */
457  virtual void get_names(size_t &np, std::vector<std::string> &pnames) {
458  np=0;
459  return;
460  }
461 
462  };
463 
464  /** \brief An EOS for the TOV solver using simple linear
465  interpolation and an optional crust EOS
466 
467  The simplest usage of this class is simply to use \ref
468  read_table() to read a tabulated EOS stored in a \ref
469  table_units object and optionally specify a separate crust EOS.
470 
471  \note This stores a pointer to the user-specified table,
472  so if that pointer becomes invalid, the interpolation will
473  fail.
474 
475  Alternatively, the user can simply specify objects
476  of type <tt>std::vector<double></tt> which store the energy
477  density, pressure, and baryon density.
478 
479  There are two methods to handle the crust-core interface. The
480  default, <tt>smooth_trans</tt> uses the crust below pressure \f$
481  P_1 \f$ (equal to the value of \ref trans_pres divided by \ref
482  trans_width) and the core above pressure \f$ P_2 \f$ (the value
483  of \ref trans_pres times \ref trans_width) and then in between
484  uses
485  \f[
486  \varepsilon(P) = [1-\chi(P)] \varepsilon_{\mathrm{crust}}(P) +
487  \chi(P) \varepsilon_{\mathrm{core}}(P)
488  \f]
489  where the value \f$ \chi(P) \f$ is determined by
490  \f[
491  \chi(P) = (P-P_1)/(P_2-P_1) \, .
492  \f]
493  This method is a bit more faithful to the original EOS tables,
494  but the matching can result in pressures which decrease with
495  increasing energy density. Alternatively the <tt>match_line</tt>
496  method uses
497  \f$ \varepsilon_1=\varepsilon_{\mathrm{crust}}(P_1) \f$ and
498  \f$ \varepsilon_2=\varepsilon_{\mathrm{core}}(P_2) \f$ and
499  \f[
500  \varepsilon(P) = (\varepsilon_2 - \varepsilon_1) \chi
501  + \varepsilon_1 \, .
502  \f]
503  (using the same expression for \f$ \chi \f$ ). This method less
504  frequently results in decreasing pressures, but can deviate
505  further from the original tables.
506 
507  Internally, energy and pressure are stored in units of \f$
508  \mathrm{M}_{\odot}/\mathrm{km}^3 \f$ and baryon density is
509  stored in units of \f$ \mathrm{fm}^{-3} \f$ . The user-specified
510  EOS table is left as is, and unit conversion is performed as
511  needed in get_eden() and other functions from the units
512  specified in the input \ref table_units object.
513 
514  \todo It might be useful to exit more gracefully when non-finite
515  values are obtained in interpolation, analogous to the
516  err_nonconv mechanism elsewhere.
517 
518  \todo Create a sanity check where core_auxp is nonzero only if
519  core_table is also nonzero. Alternatively, this complication is
520  due to the fact that this class works in two ways: one where it
521  reads a table (and adds a crust), and one where it reads in
522  vectors (with no crust). Maybe these two modes of operation
523  should be separated into two classes.
524 
525  \todo Read different vector types than just std::vector<>.
526  */
527  class eos_tov_interp : public eos_tov {
528 
529  public:
530 
531  eos_tov_interp();
532 
533  virtual ~eos_tov_interp();
534 
535  /// \name Mode of transitioning between crust and core EOS
536  //@{
537  int transition_mode;
538  static const int smooth_trans=0;
539  static const int match_line=1;
540  //@}
541 
542  /// \name Basic usage
543  //@{
544  /// Specify the EOS through a table
545  void read_table(table_units<> &eosat, std::string s_cole,
546  std::string s_colp, std::string s_colnb="");
547 
548  /** \brief Read the EOS from a set of equal length
549  vectors for energy density, pressure, and baryon density
550  */
551  void read_vectors(size_t n_core, std::vector<double> &core_ed,
552  std::vector<double> &core_pr,
553  std::vector<double> &core_nb);
554 
555  /** \brief Read the EOS from a pair of equal length
556  vectors for energy density and pressure
557  */
558  void read_vectors(size_t n_core, std::vector<double> &core_ed,
559  std::vector<double> &core_pr);
560  //@}
561 
562  /// \name Crust EOS functions
563  //@{
564  /// Default crust EOS from \ref Negele73 and \ref Baym71
565  void default_low_dens_eos();
566 
567  /// Crust EOS from \ref Shen11b
568  void sho11_low_dens_eos();
569 
570  /** \brief Crust EOS from \ref Steiner12
571 
572  Current acceptable values for \c model are <tt>APR</tt>,
573  <tt>Gs</tt>, <tt>Rs</tt> and <tt>SLy4</tt>.
574  */
575  void s12_low_dens_eos(std::string model="SLy4",
576  bool external=false);
577 
578  /** \brief Crust EOS from Goriely, Chamel, and Pearson
579 
580  From \ref Goriely10, \ref Pearson11, and \ref Pearson12 .
581  */
582  void gcp10_low_dens_eos(std::string model="BSk20",
583  bool external=false);
584 
585  /** \brief Crust EOS from \ref Newton13 given L in MeV
586 
587  Current acceptable values for \c model are <tt>PNM</tt>
588  and <tt>J35</tt>.
589 
590  */
591  void ngl13_low_dens_eos(double L, std::string model="PNM",
592  bool external=false);
593 
594  /** \brief Crust EOS from \ref Newton13 given S and L in MeV
595  and a transition density
596 
597  Note that this function works only if \f$ 28 < S < 38 \f$ MeV,
598  \f$ 25 < L < 115 \f$ MeV, \f$ 0.01 < n_t < 0.15 \f$,
599  and \f$ L > 5 S-65~\mathrm{MeV} \f$
600  . If \c fname is a string of length 0 (the default),
601  then this function looks for a file named \c newton_SL.o2
602  in the \o2 data directory specified by
603  \ref o2scl::lib_settings_class::get_data_dir() .
604  */
605  void ngl13_low_dens_eos2(double S, double L, double nt,
606  std::string fname="");
607 
608  /// Compute with no crust EOS
610  use_crust=false;
611  return;
612  }
613  //@}
614 
615  /// \name Functions used by the tov_solve class
616  //@{
617  /** \brief Given the pressure, produce the energy and number densities
618 
619  The arguments \c P and \c e should always be in \f$
620  M_{\odot}/\mathrm{km}^3 \f$ . The argument for \c nb should be
621  in \f$ \mathrm{fm}^{-3} \f$ .
622 
623  If the baryon density is not specified, it should be set to
624  zero or \ref baryon_column should be set to false
625  */
626  virtual void get_eden(double pres, double &ed, double &nb);
627 
628  /** \brief Given the pressure, produce all the remaining quantities
629 
630  The argument \c P should always be in
631  \f$ M_{\odot}/\mathrm{km}^3 \f$ .
632  */
633  virtual void get_aux(double P, size_t &nv, std::vector<double> &auxp);
634 
635  /** \brief Fill a list with strings for the names of the remaining
636  quanities
637  */
638  virtual void get_names_units(size_t &np,
639  std::vector<std::string> &pnames,
640  std::vector<std::string> &punits);
641  //@}
642 
643  /// \name Other functions
644  //@{
645  /** \brief Get the energy density from the pressure in the
646  user-specified unit system
647  */
648  virtual void get_eden_user(double pres, double &ed, double &nb);
649 
650  /** \brief Get the transition pressure (in the user-specified
651  unit system) and width
652  */
653  void get_transition(double &ptrans, double &pwidth);
654 
655  /** \brief Set the transition pressure and "width"
656 
657  Sets the transition pressure and the width (specified as a
658  number greater than unity in \c pw) of the transition between
659  the two EOSs. The transition should be in the same units of
660  the user-specified EOS. The transition is done smoothly using
661  linear interpolation between \f$ P=\mathrm{ptrans}/pmathrm{pw}
662  \f$ and \f$ P=\mathrm{ptrans} \times pmathrm{pw} \f$.
663  */
664  void set_transition(double ptrans, double pw);
665  //@}
666 
667  /// \name User EOS
668  //@{
669  /// Energy densities from full EOS
670  std::vector<double> full_vece;
671  /// Pressures from full EOS
672  std::vector<double> full_vecp;
673  /// Baryon densities from full EOS
674  std::vector<double> full_vecnb;
675  /// Number of lines in full EOS
676  size_t full_nlines;
677  //@}
678 
679 #ifndef DOXYGEN_INTERNAL
680 
681  protected:
682 
683  void internal_read();
684 
685  /// \name Crust EOS variables
686  //@{
687  /// Set to true if we are using a crust EOS (default false)
688  bool use_crust;
689 
690  /// Energy densities
691  std::vector<double> crust_vece;
692  /// Pressures
693  std::vector<double> crust_vecp;
694  /// Baryon densities
695  std::vector<double> crust_vecnb;
696  /// Number of EOS entries
697  size_t crust_nlines;
698  //@}
699 
700  /// \name User EOS
701  //@{
702  /// True if core table has been specified
703  bool core_set;
704  /// Full user EOS table
706  /// Column for energy density in EOS file
707  size_t cole;
708  /// Column for pressure in EOS file
709  size_t colp;
710  /// Column for baryon density in EOS file
711  size_t colnb;
712  /// True if an EOS has been specified
713  bool eos_read;
714  /// Number of additional columns in the core EOS
715  size_t core_auxp;
716  //@}
717 
718  /// \name Interpolation objects
719  //@{
722  interp<std::vector<double> > gen_int;
723  //@}
724 
725  /// \name Unit conversion factors for core EOS
726  //@{
727  /// Unit conversion factor for energy density (default 1.0)
728  double efactor;
729  /// Unit conversion factor for pressure (default 1.0)
730  double pfactor;
731  /// Unit conversion factor for baryon density (default 1.0)
732  double nfactor;
733  //@}
734 
735  /// \name Properties of transition
736  //@{
737  /** \brief Transition pressure (in \f$ M_{\odot}/\mathrm{km}^3 \f$)
738  */
739  double trans_pres;
740  /// Transition width (unitless)
741  double trans_width;
742  //@}
743 
744 #endif
745 
746  };
747 
748 #ifndef DOXYGEN_NO_O2NS
749 }
750 #endif
751 
752 #endif
753 
754 
std::vector< double > full_vecnb
Baryon densities from full EOS.
Definition: eos_tov.h:674
void set_baryon_density(double nb, double ed)
Set the baryon density.
Definition: eos_tov.h:423
double pr1
The pressure at ed1.
Definition: eos_tov.h:395
virtual void get_eden(double P, double &e, double &nb)=0
Given the pressure, produce the energy and number densities.
Linear EOS .
Definition: eos_tov.h:381
void ngl13_low_dens_eos2(double S, double L, double nt, std::string fname="")
Crust EOS from Newton13 given S and L in MeV and a transition density.
virtual void get_eden(double P, double &e, double &nb)
Given the pressure, produce the energy and number densities.
Definition: eos_tov.h:365
void default_low_dens_eos()
Default crust EOS from Negele73 and Baym71.
std::vector< double > crust_vecnb
Baryon densities.
Definition: eos_tov.h:695
void set_baryon_density(double nb, double ed)
Set the baryon density.
Definition: eos_tov.h:209
std::vector< double > crust_vecp
Pressures.
Definition: eos_tov.h:693
void sho11_low_dens_eos()
Crust EOS from Shen11b.
size_t crust_nlines
Number of EOS entries.
Definition: eos_tov.h:697
void vector_deriv_interp(size_t n, ovec_t &v, vec2_t &dv, size_t interp_type=itp_cspline)
size_t full_nlines
Number of lines in full EOS.
Definition: eos_tov.h:676
std::vector< double > full_vecp
Pressures from full EOS.
Definition: eos_tov.h:672
double eps0
Index (default 0.0)
Definition: eos_tov.h:402
void ngl13_low_dens_eos(double L, std::string model="PNM", bool external=false)
Crust EOS from Newton13 given L in MeV.
virtual void get_names_units(size_t &np, std::vector< std::string > &pnames, std::vector< std::string > &punits)
Fill a list with strings for the names of the remaining quanities.
Definition: eos_tov.h:93
double ed1
The energy density for which the baryon density is known.
Definition: eos_tov.h:187
double pr1
The pressure at ed1.
Definition: eos_tov.h:317
double trans_pres
Transition pressure (in )
Definition: eos_tov.h:739
virtual void get_names_units(size_t &np, std::vector< std::string > &pnames, std::vector< std::string > &punits)
Fill a list with strings for the names of the remaining quanities.
bool eos_read
True if an EOS has been specified.
Definition: eos_tov.h:713
double K
Coefficient (default 1.0)
Definition: eos_tov.h:321
double ed1
The energy density for which the baryon density is known.
Definition: eos_tov.h:391
virtual void get_eden_user(double pres, double &ed, double &nb)
Get the energy density from the pressure in the user-specified unit system.
virtual void get_names(size_t &np, std::vector< std::string > &pnames)
Fill a list with strings for the names of the remaining quanities.
Definition: eos_tov.h:457
void read_table(table_units<> &eosat, std::string s_cole, std::string s_colp, std::string s_colnb="")
Specify the EOS through a table.
double pfactor
Unit conversion factor for pressure (default 1.0)
Definition: eos_tov.h:730
void set_baryon_density(double nb, double ed)
Set the baryon density.
Definition: eos_tov.h:351
double Pstar
The parameter with units of pressure in units of solar masses per km cubed (default value ) ...
Definition: eos_tov.h:205
void gcp10_low_dens_eos(std::string model="BSk20", bool external=false)
Crust EOS from Goriely, Chamel, and Pearson.
std::vector< double > crust_vece
Energy densities.
Definition: eos_tov.h:691
bool core_set
True if core table has been specified.
Definition: eos_tov.h:703
double nb1
The baryon density at ed1.
Definition: eos_tov.h:387
virtual void get_aux(double P, size_t &nv, std::vector< double > &auxp)
Given the pressure, produce all the remaining quantities.
double n
Index (default 3.0)
Definition: eos_tov.h:324
double ed1
The energy density for which the baryon density is known.
Definition: eos_tov.h:313
An EOS for the TOV solver using simple linear interpolation and an optional crust EOS...
Definition: eos_tov.h:527
#define O2SCL_ERR2(d, d2, n)
table_units * core_table
Full user EOS table.
Definition: eos_tov.h:705
A EOS base class for the TOV solver.
Definition: eos_tov.h:47
virtual void get_eden(double pres, double &ed, double &nb)
Given the pressure, produce the energy and number densities.
size_t core_auxp
Number of additional columns in the core EOS.
Definition: eos_tov.h:715
virtual void get_eden(double P, double &e, double &nb)
Given the pressure, produce the energy and number densities.
Definition: eos_tov.h:231
Class to solve the Tolman-Oppenheimer-Volkov equations.
Definition: tov_solve.h:273
int solve_u_rp_fun(size_t bv, const std::vector< double > &bx, std::vector< double > &by)
Solve to compute profiles.
Definition: eos_tov.h:286
size_t cole
Column for energy density in EOS file.
Definition: eos_tov.h:707
void set_coeff_index(double coeff, double index)
Desc.
Definition: eos_tov.h:337
void no_low_dens_eos()
Compute with no crust EOS.
Definition: eos_tov.h:609
Standard polytropic EOS .
Definition: eos_tov.h:303
virtual void get_aux(double P, size_t &np, std::vector< double > &auxp)
Given the pressure, produce all the remaining quantities.
Definition: eos_tov.h:85
void s12_low_dens_eos(std::string model="SLy4", bool external=false)
Crust EOS from Steiner12.
void set_cs2_eps0(double cs2_, double eps0_)
Desc.
Definition: eos_tov.h:415
double cs2
Coefficient (default 1.0)
Definition: eos_tov.h:399
double pr1
The pressure at ed1.
Definition: eos_tov.h:191
virtual void get_names(size_t &np, std::vector< std::string > &pnames)
Fill a list with strings for the names of the remaining quanities.
Definition: eos_tov.h:257
void get_transition(double &ptrans, double &pwidth)
Get the transition pressure (in the user-specified unit system) and width.
void check_nb(double &avg_abs_dev, double &max_abs_dev)
Check that the baryon density is consistent with the EOS.
Definition: eos_tov.h:103
double nb1
The baryon density at ed1.
Definition: eos_tov.h:183
void read_vectors(size_t n_core, std::vector< double > &core_ed, std::vector< double > &core_pr, std::vector< double > &core_nb)
Read the EOS from a set of equal length vectors for energy density, pressure, and baryon density...
double efactor
Unit conversion factor for energy density (default 1.0)
Definition: eos_tov.h:728
double nfactor
Unit conversion factor for baryon density (default 1.0)
Definition: eos_tov.h:732
void set_transition(double ptrans, double pw)
Set the transition pressure and "width".
size_t colnb
Column for baryon density in EOS file.
Definition: eos_tov.h:711
bool use_crust
Set to true if we are using a crust EOS (default false)
Definition: eos_tov.h:688
int verbose
Control for output (default 1)
Definition: eos_tov.h:68
bool baryon_column
Set to true if the baryon density is provided in the EOS (default false)
Definition: eos_tov.h:54
std::vector< double > full_vece
Energy densities from full EOS.
Definition: eos_tov.h:670
virtual void get_eden(double P, double &e, double &nb)
Given the pressure, produce the energy and number densities.
Definition: eos_tov.h:437
virtual void get_aux(double P, size_t &np, std::vector< double > &auxp)
Given the pressure, produce all the remaining quantities.
Definition: eos_tov.h:249
double nb1
The baryon density at ed1.
Definition: eos_tov.h:309
double trans_width
Transition width (unitless)
Definition: eos_tov.h:741
The Buchdahl EOS for the TOV solver.
Definition: eos_tov.h:177
size_t colp
Column for pressure in EOS file.
Definition: eos_tov.h:709
virtual void get_aux(double P, size_t &np, std::vector< double > &auxp)
Given the pressure, produce all the remaining quantities.
Definition: eos_tov.h:449

Documentation generated with Doxygen. Provided under the GNU Free Documentation License (see License Information).
Hosted at Get Object-oriented Scientific Computing
Lib at SourceForge.net. Fast, secure and Free Open Source software
downloads..