EOS generator eos_nuclei

The main code for generating tables is an executable called eos_nuclei. It has a command-line interface with a help feature, so try eos_nuclei -help.

Class eos_nuclei

class eos_nuclei : public eos

Solve for the EOS including nuclei.

Todo

Class eos_nuclei:

  • Rename n_nB2 to n_nB, etc.

  • Make child of eos_sn_base

  • Use c loaded instead of testing n_nB2==0

  • Move select_high_T to the parent

  • 1/15/22: I’m not sure if fix_cc() is really useful or not anymore?

  • 1/15/22: I’m not sure the non-derivative virial solver is needed anymore?

  • Future: Allow different form for the NS fit

  • Future: Allow user to specify where data files are located or to manually specify the nuclear model?

  • Future: Allow RMF rather than just Skyrme.

Stability tensors

o2scl::tensor_grid dmundYe
o2scl::tensor_grid dmundnB
o2scl::tensor_grid dmupdYe
o2scl::tensor_grid dsdT
o2scl::tensor_grid dsdnB
o2scl::tensor_grid dsdYe
o2scl::tensor_grid egv[4]
o2scl::tensor_grid tg_cs2
o2scl::tensor_grid tg_cs2_hom
o2scl::tensor_grid tg_Fint_old
o2scl::tensor_grid tg_F_old
size_t n_stability_fail
o2scl::part_funcs pfuncs

Partition functions for the nuclei.

o2scl::boson_rel relb

Particle object for bosons.

o2scl::boson_eff effb

Particle object for bosons.

bool inc_hrg

If true, include a hadron resonance gas (default false)

int solve_hrg(size_t nv, const ubvector &x, ubvector &y, double nB, double Ye, double T)

Solve for charge neutrality and fixed baryon fraction with a hadron resonance gas.

This is used in nuc_matter() .

Nuclear masses and their fits

o2scl::nucmass_mnmsk m95

Theoretical nuclear masses.

o2scl::nucmass_hfb_sp hfb

HFB masses for spin predictions.

o2scl::nucmass_ame ame

Experimental nuclear masses.

o2scl::nucmass_frdm frdm

Theoretical nuclear masses far from stability.

o2scl::nucmass_fit nm_fit

Fit theory masses.

int fit_frdm(std::vector<std::string> &sv, bool itive_com)

Fit the FRDM mass model.

<no parameters>

Fit the FRDM mass model.

Nucleus objects

std::vector<o2scl::nucleus> nuclei
o2scl::nucleus *nuc_alpha
o2scl::nucleus *nuc_deut
o2scl::nucleus *nuc_trit
o2scl::nucleus *nuc_he3
o2scl::nucleus *nuc_li4
o2scl::nucleus *nuc_heavy

Mathematical algorithm objects

o2scl::inte_qag_gsl iqg

Integrator.

o2scl::mroot_hybrids mh

Solver.

o2scl::root_brent_gsl rbg

Bracketing solver.

The nuclear partition function

partition_func part_func

Object which integrates partition functions.

ubvector Sneut

Neutron separation energies (in MeV)

ubvector Sprot

Proton separation energies (in MeV)

ubvector vomega

Partition function.

ubvector vomega_prime

Derivative of partition function with respect to temperature, in units of \( \mathrm{fm} \).

Parameters modifiable by the CLI user

double file_update_time

The time (in seconds) between output file updates for generate_table() (default 1800)

int file_update_iters

The number of iterations between file updates (default 100000)

int fd_A_max

The maximum value of A for a fixed distribution when alg_mode is 4 (default 600)

bool extend_frdm

If true, attempt to extend FRDM beyond the drip lines (default false).

double max_ratio

The maximum value of \( N/Z \) or \( Z/N \) (default 7)

double mh_tol_rel

Relative tolerance for the solver in the eos_fixed_dist() function (default \( 10^{-6} \))

std::string ext_guess

Filename containing separate table to use as a guess for the generate-table command (default is the empty string)

std::string nucleon_func

Function for delta Z and delta N in the single nucleus approximation.

This is the function used to specify the range of Z and N which is varied to find the smallest free energy

double max_time

Maximum time, in seconds, for the generate-table command. The default is 0.0 which is interpreted as no maximum time.

bool show_all_nuclei

If true, show all nuclei considered at every point (default false)

This applies to the point-nuclei command and the eos_vary_ZN() function.

bool rnuc_less_rws

If true, ensure that the nuclear radius is less than the Wigner-Seitz radius (default true)

std::map<std::string, std::string> vdet_units

Units of objects in the vdet arrays.

bool recompute

If true, recompute all points, irrespective of the value of the convergence flag (default false)

This setting is used in point-nuclei and generate-table.

bool verify_only

If true, verify points only and do not attempt to solve (default false)

This setting is principally used by the verify command which automatically flips this from false to true and then back.

std::string edge_list

If not empty, then recompute points where Z, A, log_xn, or log_xp reach a maximum or minimum (default is empty)

int alg_mode

Algorithm mode (default 4)

0 for SNA, 1 for old SNA method, 2 for vary dist., 3 for old vary dist., and 4 for fixed dist.

int fixed_dist_alg

Algorithm for eos_fixed_dist()

Modify the algorithm for the eos_fixed_dist() function. The 1s digit is the number of solves minus 1, the 10s digit is the number of brackets divided by 10, the 100s digit is the number of minimizes, and the 1000s digit is the number of random guesses to try divided by 1000. The default is 1111. Other good options are 1319, 1919, 1999, and 9999.

int six_neighbors

If true, when computing a point, perform the calculation using some of the neighboring points as an initial guess (default 0)

Values greater than 0 use the point at the next smallest density, values greater than 1 use the point at the next largest density, values greater than 2 use points at the next largest and next smallest temperature, and values greater than 4 use the next largest and smallest electron fraction.

int function_verbose

A new function verbose parameter.

Verbose for individual functions (default value 11111).

\t1s digit: fixed_ZN()

\t10s digit: vary_ZN()

\t100s digit: fixed_dist()

\t1000s digit: vary_dist()

\t10000s digit:

store_point().

bool propagate_points

If true, use previously computed points (or guesses) as an initial guess to compute adjacent points (default true)

bool survey_eqs

If true, survey the nB and Ye equations near a failed point (default false)

bool derivs_computed

If true, output all of the data necessary for a full EOS.

If true, include Eint, Pint, Sint, mun, and mup. If include_eg is additionally true, then add F, E, P, and S.

bool baryons_only

Always true, included for consistency with o2scl::eos_sn_base.

bool with_leptons

If true, include electrons and photons.

Requires that derivs_computed is also true

void store_hrg(double mun, double mup, double nn, double np, double T, o2scl::table_units<> &tab)

Desc.

Other internal objects

ubvector Ec

Coulomb energy (in \( \mathrm{fm}^{-1} \) )

This quantity is computed in solve_nuclei() and then used later in eos_fixed_dist().

o2scl::vec_index vi

Dictionary for mapping buffers to physical quantities.

bool loaded

True if an EOS is currently loaded.

std::vector<double> fd_rand_ranges

Ranges for randomly selected ranges in eos_fixed_dist()

o2scl::slack_messenger slack

Object for sending slack messages.

std::string Ye_list

The list of electron fractions to consider for the generate-table command.

Can be several comma-separated ranges e.g. “1-3,5-7,59-60”. Zero-indexed.

Other internal physics objects

eos_had_skyrme_ext skyrme_ext

Extended Skyrme model for finite-temperature corrections.

eos_had_lim_holt lim_holt

Extended Skyrme model for finite-temperature corrections.

Main EOS table storage

o2scl::tensor_grid tg_log_xn
o2scl::tensor_grid tg_log_xp
o2scl::tensor_grid tg_flag
o2scl::tensor_grid tg_F
o2scl::tensor_grid tg_E
o2scl::tensor_grid tg_P
o2scl::tensor_grid tg_S
o2scl::tensor_grid tg_Fint
o2scl::tensor_grid tg_Eint
o2scl::tensor_grid tg_Pint
o2scl::tensor_grid tg_Sint
o2scl::tensor_grid tg_mun
o2scl::tensor_grid tg_mup
o2scl::tensor_grid tg_mue
o2scl::tensor_grid tg_Z
o2scl::tensor_grid tg_A
o2scl::tensor_grid tg_Xn
o2scl::tensor_grid tg_Xp
o2scl::tensor_grid tg_Xalpha
o2scl::tensor_grid tg_Xnuclei
o2scl::tensor_grid tg_Ymu
o2scl::tensor_grid tg_Xd
o2scl::tensor_grid tg_Xt
o2scl::tensor_grid tg_XHe3
o2scl::tensor_grid tg_XLi4
o2scl::tensor_grid tg_A_min
o2scl::tensor_grid tg_A_max
o2scl::tensor_grid tg_NmZ_min
o2scl::tensor_grid tg_NmZ_max
bool include_detail

If true, include details in EOS file (default false)

Detail storage

o2scl::tensor_grid tg_zn
o2scl::tensor_grid tg_zp
o2scl::tensor_grid tg_F1
o2scl::tensor_grid tg_F2
o2scl::tensor_grid tg_F3
o2scl::tensor_grid tg_F4
o2scl::tensor_grid tg_Un
o2scl::tensor_grid tg_Up
o2scl::tensor_grid tg_msn
o2scl::tensor_grid tg_msp
o2scl::tensor_grid tg_g
o2scl::tensor_grid tg_dgdT
o2scl::tensor_grid tg_sigma
o2scl::tensor_grid tg_omega
o2scl::tensor_grid tg_rho

Other parameter objects

o2scl::cli::parameter_bool p_inc_hrg
o2scl::cli::parameter_bool p_survey_eqs
o2scl::cli::parameter_bool p_extend_frdm
o2scl::cli::parameter_bool p_show_all_nuclei
o2scl::cli::parameter_int p_fd_A_max
o2scl::cli::parameter_bool p_recompute
o2scl::cli::parameter_bool p_verify_only
o2scl::cli::parameter_string p_edge_list
o2scl::cli::parameter_string p_ext_guess
o2scl::cli::parameter_int p_six_neighbors
o2scl::cli::parameter_bool p_full_results
o2scl::cli::parameter_bool p_rnuc_less_rws
o2scl::cli::parameter_bool p_include_eg
o2scl::cli::parameter_bool p_propagate_points
o2scl::cli::parameter_bool p_include_detail
o2scl::cli::parameter_bool p_strange_axis
o2scl::cli::parameter_double p_mh_tol_rel
o2scl::cli::parameter_double p_max_time
o2scl::cli::parameter_string p_nucleon_func
o2scl::cli::parameter_int p_alg_mode
o2scl::cli::parameter_int p_fixed_dist_alg
o2scl::cli::parameter_int p_function_verbose
o2scl::cli::parameter_string p_Ye_list
o2scl::cli::parameter_double p_max_ratio
o2scl::cli::parameter_double p_file_update_time
o2scl::cli::parameter_int p_file_update_iters

MPI message values

static const int message_continue = 0
static const int message_done = 1

Flag values

static const int iflag_empty = 0

Point is empty.

static const int iflag_in_progress_empty = -1

Point is in progress, and empty.

static const int iflag_in_progress_with_guess = -5

Point is in progress, and has initial guess.

static const int iflag_guess = 5

Point has initial guess.

static const int iflag_done = 10

Point is finished.

Constructor and destructor

eos_nuclei()
virtual ~eos_nuclei()

Functions for the main algorithm

double solve_nuclei_ld(double x2, size_t nv, const ubvector &x, double nb, double ye, double T, int ix, double &mun_gas, double &mup_gas, o2scl::thermo &th_gas)

Use only one of the two equations for a function for a root bracketing algorithm.

double solve_nuclei_min(size_t nv, const ubvector &x, double nb, double ye, double T, double &mun_gas, double &mup_gas, o2scl::thermo &th_gas)

Solve for the nuclei via minimization.

int solve_nuclei(size_t nv, const ubvector &x, ubvector &y, double nb, double ye, double T, int loc_verbose, double &mun_gas, double &mup_gas, o2scl::thermo &th_gas, std::map<std::string, double> &vdet)

Construct equations to solve for a fixed baryon density and electron fraction.

int eos_fixed_ZN(double nb, double ye, double T, double &log_xn, double &log_xp, size_t nuc_Z1, size_t nuc_N1, o2scl::thermo &thx, double &mun_full, double &mup_full)

Determine the EOS presuming a fixed single heavy nucleus and solving for the log (base 10) of the free neutron and proton abundances.

int store_point(size_t i_nB, size_t i_Ye, size_t i_T, double nB, double Ye, double T, o2scl::thermo &th, double log_xn, double log_xp, double Zbar, double Nbar, double mun_full, double mup_full, ubvector &X, double A_min, double A_max, double NmZ_min, double NmZ_max, double flag, std::map<std::string, double> &vdet)

Store data in the tensor objects.

int eos_vary_ZN(double nb, double ye, double T, double &log_xn, double &log_xp, size_t &nuc_Z1, size_t &nuc_N1, o2scl::thermo &thx, double &mun_full, double &mup_full, bool nu_nuclei = false)

Determine the EOS allowing the Z and N of the nucleus to vary.

int eos_fixed_dist(double nB, double Ye, double T, double &log_xn, double &log_xp, o2scl::thermo &thx, double &mun_full, double &mup_full, int &A_min, int &A_max, int &NmZ_min, int &NmZ_max, std::map<std::string, double> &vdet, bool dist_changed, bool no_nuclei)

Determine the EOS presuming a distribution of nuclei with fixed limits in A and \( N-Z \).

int eos_fixed_dist_fix_table(double nB, double Ye, double T, double &log_xn, double &log_xp, o2scl::thermo &thx, double &mun_full, double &mup_full, int &A_min, int &A_max, int &NmZ_min, int &NmZ_max, bool dist_changed, bool no_nuclei)

Determine the EOS presuming a distribution of nuclei with fixed limits in A and \( N-Z \) but used to fix table only.

int nuc_matter(double nB, double Ye, double T, double &log_xn, double &log_xp, double &Zbar, double &Nbar, o2scl::thermo &thx, double &mun_full, double &mup_full, int &A_min, int &A_max, int &NmZ_min, int &NmZ_max, std::map<std::string, double> &vdet)

Compute the EOS presuming homogeneous nuclear matter.

int eos_vary_dist(double nB, double Ye, double T, double &log_xn, double &log_xp, double &Zbar, double &Nbar, o2scl::thermo &thx, double &mun_full, double &mup_full, int &A_min, int &A_max, int &NmZ_min, int &NmZ_max, std::map<std::string, double> &vdet, bool dist_changed, bool no_nuclei)

Determine the EOS presuming a distribution of nuclei and optimizing the limits in A and \( N-Z \).

int generate_table(std::vector<std::string> &sv, bool itive_com)

Generate an EOS table.

[out file]

This command is the full MPI calculation of the EOS table, given a model and using the current grid. If no output file name is specified, the results are placed in a file called eos_nuclei.o2 .

int save_compose(std::vector<std::string> &sv, bool itive_com)

Save in CompOSE format.

<output file prefix>

Save an EOS table in the CompOSE format

EOS post-processing functions

int eos_deriv(std::vector<std::string> &sv, bool itive_com)

Compute derivatives numerically.

<no parameters>

This command requires that an EOS has been loaded. It uses Steffen’s monotonic interpolation to compute the neutron and proton chemical potentials, and the interacting part of the internal energy, pressure, and entropy. The value of derivs_computed is set to true.

int interp_point(std::vector<std::string> &sv, bool itive_com)

Interpolate the EOS to fix the table near a point.

<nB> <Ye> <T MeV> <window> <st.o2>

This function requires that an EOS with leptons has been loaded.

int interp_internal(size_t i_fix, size_t j_fix, size_t k_fix, size_t window, interpm_krige_eos &ike)

The internal interpolation function.

Use interpolation to fix the table in the neighborhood of size window around (i_fix,j_fix,k_fix) using interpolation object ike.

int interp_fix_table(std::vector<std::string> &sv, bool itive_com)

Use interpolation to fix an entire table.

<input stability file> <window> <output table> <output stability file>

Under development.

int add_eg(std::vector<std::string> &sv, bool itive_com)

Add electrons and photons.

<no parameters>

This command is typically used after derivatives are computed with the “eos-deriv” command. This command requires a table has been loaded with derivs_computed equal to true, but does not require that a model has been selected.

int eg_table(std::vector<std::string> &sv, bool itive_com)

Construct an electrons and photon table.

<output file>

Construct a file consisting only of the electron and photon EOS. The grid is determined by the current grid specification. Muons are included if eos::include_muons is set to true. The output file has five tensor grid objects, F, E, P, S, and mue. If muons are included, then the file also includes Ymu. The electron (and muon) masses are also written to the table.

This command works independent of whether or not a table was loaded or if a model was selected. However, if a table is loaded when this command is invoked, the current table is cleared.

int beta_table(std::vector<std::string> &sv, bool itive_com)

Construct a table in beta equilibrium.

<filename>

Experimental. From a previously created table (which has both derivatives and leptons), construct a table over baryon density and temperature at beta equilibrium and store in the file named filename. Rank 2 tensors A, E, Eint, F, Fint, P, Pint, S, Sint, XHe3, XLi4, Xalpha, Xd, Xn, Xnuclei, Xp, Xt, Z, log_xn, log_xp, mue, mun, and mup are stored, as well as a tensor called Ye which stores the electron fraction at that density and temperature.

int edit_data(std::vector<std::string> &sv, bool itive_com)

Edit an EOS table.

<select func.> [tensor to modify] [value func.]

The edit-data command counts the number of (nB,Ye,T) points matching the criteria specified in <select func.>. If the remaining two arguments are given, then the values of [tensor to modify] for the selected points are changed to the result of the function [value func.].

int merge_tables(std::vector<std::string> &sv, bool itive_com)

Merge two output tables to create a third.

<input file 1> <input file 2> <output file>

Tables can only be merged if their grids and settings match. If the Fint table is anomalously small or large or not-finite, then this function calls the error handler. Otherwise, for each point in (nB,Ye,T), there are four reasons for which a point is copied from the second table to the first: (i) they both have flag=10 but the second has a smaller Fint, (ii) the second has flag=10 but the first does not, (iii) they both have flags less than 10 but the second has a non-zero flag with a smaller Fint, or (iv) the second table has a non-zero flag and the first does not. After the merge, the number of points modified is reported.

int compare_tables(std::vector<std::string> &sv, bool itive_com)

Compare two output tables.

<input file 1> <input file 2> [quantity]

Compare two EOS tables. If the optional argument “)+

is unspecified, then all quantities are compared. If [quantity] “+ is specified, then only that particular quantitiy is compared. “+

Only points for which flag=10 in both tables are compared. “+ If derivs_computed is true, then Pint, mun, and “+

mup are available for comparisons. If with_leptons is “+ true, then “+

F, E, P, and S, are also available for comparisons. Any current “+ EOS data stored is cleared before the comparison. If the “+ nB, Ye, or T grids do not match, then no comparison is performed.

int stats(std::vector<std::string> &sv, bool itive_com)

Output convergence statistics and simple checks.

<no parameters>

If an EOS is loaded, this function counts the number of points with each flag value, checks that the nuclear fractions add up to 1, checks that the free energy internal energy, and entropy are consistent, and checks the thermodynamic identity.

int point_nuclei(std::vector<std::string> &sv, bool itive_com)

Compute and/or show EOS results at one (n_B,Y_e,T) point.

<n_B> <Y_e> <T (MeV)> [log(xn) log(xp) Z N] [alg_mode 2-4: log(xn) log(xp) A_min A_max NmZ_min NmZ_max] [fname]

If an EOS is loaded, then the n_B, Y_e, and T values are modified to ensure that they lie on a grid point. If an initial guess is specified on the command line, it is used even if there is a good guess already in the table. If the flag is not 10 or if recompute is true, then the EOS is recomputed. If an EOS is loaded or the recompute was successful, then the results are output to the screen. If the point was successful it is stored in the current tables. If show_all_nuclei is true, then a file named dist.o2 is created which holds the full nuclear distribution.

int test_random(std::vector<std::string> &sv, bool itive_com)

Test an EOS at random points in (nB,Ye,T)

<n_tests> [“lg”]

This function tests the EOS at randomly chosen points in (nB,Ye,T) space. If the new calculation and the stored result disagree, then the new result is stored in the table. A table must be loaded and the full model must be specified (either with select-model or from the table).

If the additional argument “lg” (for liquid-gas) is given, then the random points are all selected at densities between 0.01 and 0.15 1/fm^3.

File I/O functions

int load(std::vector<std::string> &sv, bool itive_com)

Load an EOS table.

<filename>

Loads an EOS table in filename to memory. In the case where MPI is used, only one MPI rank reads the table at a time.

int output(std::vector<std::string> &sv, bool itive_com)

Output an EOS table to a file.

<filename>

Write an EOS table to a file. This command does not currently support multiple MPI ranks.

int write_results(std::string fname)

Write results to an HDF5 file.

Idea for Future:

Eventually replace this with eos_sn_base::output()?

int read_results(std::string fname)

Read results from an HDF5 file.

Idea for Future:

Eventually replace this with eos_sn_base::load()?

int write_nuclei(std::vector<std::string> &sv, bool itive_com)

Write the nuclear masses to an HDF5 file.

<output file>

This command creates a o2scl::table object and outputs that table to the specified file. The table contains a list of all nuclei used in the EOS, together with their spin degeneracy, mass, binding energy, neutron and proton separation energies, and flags which indicate where the masses and spins are obtained.

void load_nuclei()

Load nuclear masses.

This function is called in main().

void write_nuclei_intl(std::string fname)

Write nuclear masses to a file.

This function is called by write_nuclei().

Miscellaneous functions

virtual void setup_cli_nuclei(o2scl::cli &cli)

Setup the command-line interface.

void new_table()

Initialize tensors for a new EOS table.

int check_virial(std::vector<std::string> &sv, bool itive_com)

Check the virial EOS.

<no parameters>

This function checks the solver by using it to compute the EOS over a wide range of densities and temperatures This function is particularly good for checking to make sure that the virial part of the EOS does not lead to any discontinuities. For example, o2graph -read check_virial.o2 zn -Ye-slice 0.3 -set logz 1 -den-plot slice -show.

This function creates a file ‘check_virial.o2’ with four tensor_grid objects which store the neutron and proton fugacities. This file can be plotted with, e.g. o2graph -set logz 1 -read check_virial.o2 zn -set logx 1 -set logy 1 -set colbar 1 -to-table3d 0 2 slice 0.01 -den-plot slice -show.

int increase_density(std::vector<std::string> &sv, bool itive_com)

Use low densities to improve results at high densities.

<nB low> <nB high> <Ye low> <Ye high> <T low> <T high> <output

file>

This function computes the EOS at higher densities using initial guess from lower densities. It is particularly helpful in getting the phase transition between nuclei and nuclear matter correct. The outermost loop is temperature, the second loop is electron fraction and the inner loop is density. This function requires a table has been loaded and the EOS is specified. It has no parallelization support.

int fix_cc(std::vector<std::string> &sv, bool itive_com)

Increase nB to optimize the phase transition.

<output file>

Help

int verify(std::vector<std::string> &sv, bool itive_com)

Verify the EOS.

“random” <n_tests> <output file>, “random_lg” <n_tests> <output file>, “all” <output file>, “all_lg” <output

file>, or “point” <output file> <nB> <Ye> <T>

Verify the EOS, recompute if a point fails and the write final results to the specified output file. An table must be loaded and the full model must be specified (either with select-model or from the table).

This function only verifies that the baryon density and electron fraction equations are solved to within the current tolerance and does not attempt to solve them again. The test-random function is different, it actually re-solves the equations to show the answer is correct. Thus, this function requires a bit less running time at each point.

The first argument is a mode parameter which determines which points will be verified. If it is “random”, then random points will be verified. If it is “random_lg”, then random points with baryon densities from 0.01 to 0.16 fm^{-3} will be verified. If is is “all”, then all points will be verified. If it is “all_lg” then all points with baryon densities from 0.01 to 0.15 fm^{-3} will be verified. If it is “point”, then only the specified point will be verified. If any point fails, the verification procedure stops immediately and no further points are tested.

The baryon density, electron fraction, and temperature are output to the screen, as well as the integer from eos_vary_dist(), which is 0 for success.

This function does not appear to require electrons and photons but only works with alg_mode=4.

int mcarlo_nuclei(std::vector<std::string> &sv, bool itive_com)

Monte Carlo results with nuclei.

Params.

Help.

int mcarlo_nuclei2(std::vector<std::string> &sv, bool itive_com)

Monte Carlo results with nuclei (v2)

<nB> <Ye> <T> <N> <filename>

Help

int mcarlo_beta(std::vector<std::string> &sv, bool itive_com)

Monte Carlo neutrino opacity in beta equilibrium.

<filename> [n_point]

Help

int mcarlo_neutron(std::vector<std::string> &sv, bool itive_com)

Monte Carlo neutrino opacity in pure neutron matter.

<filename> [n_point]

Help

void compute_X(double nB, ubvector &X)

Compute the baryon number fractions and put them in X.

int select_high_T_internal(int option)

Select high-temperature EOS.

Note

This function is called by the constructor and thus cannot be virtual

int stability(std::vector<std::string> &sv, bool itive_com)

Compute the second derivatives and the stability matrix.

<output file> or <nB> <Ye> <T> or <inB low> <inB high> <iYe low> <iYe high> <iT low> <iT high>

This command computes the second derivatives, speed of sound, and stability matrix of the current EOS table. A table with electrons must be loaded and the full model must be specified (either with select-model or from the table).

If one output file argument is specified, the stability command creates an output file with several additional data objects for the second derivatives of the free energy, the eigenvalues of the curvature matrix, and the squared speed of sound.

Otherwise, if a density, electron fraction, and temperature are specified, the stability command compares the heterogeneous and homogeneous matter sound speeds at the grid point nearest to the specified values.

This command currently requires that a model has been selected so it can compute the speed of sound of homogeneous matter at acausal points for comparison.

int select_high_T(std::vector<std::string> &sv, bool itive_com)

Select the high-temperature Skyrme EOS.

<index>

Select 0 for the original DSH fit, 1 for NRAPR, 2 for Sk chi 414, 3 for Skchi450, 4 for Skchi500, 5 for ?, “+ and 6 for Sk chi m* (the default).

double f_min_search(size_t nvar, const ubvector &x, double nb, double ye, double T)

Compute eos with nuclei by searching minimum.

int init_function(size_t dim, const ubvector &x, ubvector &y)

Initialization for differential evolution approach.

int maxwell(std::vector<std::string> &sv, bool itive_com)

Maxwell construction.

Params.

Help.

int max_fun(size_t nv, const ubvector &x, ubvector &y, o2scl::interp_vec<std::vector<double>, ubvector> &itp_P, o2scl::interp_vec<std::vector<double>, ubvector> &itp_mun)

Desc.