Class nstar_rot (o2scl)

O2scl : Class List

class nstar_rot

Rotating neutron star class based on RNS v1.1d from N. Stergioulas et al.

Several changes have been made to the original code. The code using Numerical Recipes has been removed and replaced with an equivalent based on GSL and . The overall interface has been changed and some code has been updated with C++ equivalents.

Initial guess

The original RNS code suggests that the initial guess is typically a star with a smaller angular momentum.

References

The original RNS v1.1d can be obtained from http://www.gravity.phys.uwm.edu/rns/, and you may find Nick Stergioulas’s web page http://www.astro.auth.gr/~niksterg/, or Sharon Morsink’s page http://fermi.phys.ualberta.ca/~morsink/ useful.

See [Bonazzola73], [Bonazzola94], [Cook92], [Cook94], [Friedman88], [Gourgoulhon94], [Komatsu89], [Laarakkers99], [Nozawa98], [Stergioulas95], and [Stergioulas03].

Todo

In class nstar_rot:

  • Better documentation is needed everywhere.

  • Test the resize() function

  • It appears that KAPPA and KSCALE contant an arbitrary constant, try changing it and see if we get identical results. Try to ensure that the values are consistent between the eos_nstar_rot class and the nstar_rot class.

  • Variables r_is_gp, p_center, h_center, and others only occur in spherical_star(), integrate(), and make_center(), and can be moved to function parameters or otherwise reorganized.

  • Directly compare spherical_star() output with tov_solve results

Idea for Future:

Consider moving the int_z() algorithm to vector_derint.h

Remove the unit-indexed arrays everywhere.

Try moving some of the storage to the heap?

Some of the arrays seem larger than necessary.

The function o2scl::nstar_rot::new_search() is inefficient because it has to handle the boundary conditions separately. This could be improved.

Make the solvers more robust. The ang_vel() and ang_vel_alt() functions appear particularly unstable.

Equation of State

The thermodynamic enthalpy is

\[ H = E + P V \]
and the specific enthalpy is often written \( h = H/M \) where \( M \) is some mass scale. In GR this is often chosen to be a generic “baryon mass”, \( m_B \) so that
\[ h = \frac{\varepsilon + P}{\rho} \]
where \( \varepsilon \equiv E/V \) and \( \rho \equiv m_B/V \). The pseudo-enthalpy is typically defined by
\[ d\hat{h} = dP/(P+\varepsilon) \]
(sometimes denoted as \( H \) or \( h \) and called h in the code, but referred to here as \( \hat{h} \) to avoid confusion with the enthalpy and specific enthalpy defined above). Additionally, this quantity is sometimes defined with an additional factor of \( c^2 \). At \( T=0 \) and presuming only one conserved charge, one can use the Gibbs-Duhem relation, \( dP = (\varepsilon + P)/\mu d \mu \) to write
\[ d \hat{h} = \frac{d \mu}{\mu} \]
\[ \hat{h} = \ln \mu + C \]
The constant \( C \) is often chosen to be \( C \equiv - \ln \mu(P=0) \), so that
\[ \hat{h} = \ln \frac{\mu}{\mu(P=0)} \]
Typically, \( \mu(P=0) \) is around 931 MeV. The eos_nstar_rot_interp class takes an EOS tabulated in powers of \( \mathrm{fm}^{-1} \) and recasts it into a form which can be used by nstar_rot.

Coordinate system

The space is mapped to coordinates \( (s,\mu) \) where

\[ r = r_{\mathrm{eq}} \left( \frac{s}{1-s} \right) \]
and \( r_{\mathrm{eq}} \) is the equatorial radius, thus
\[ s = \frac{r}{r+r_{\mathrm{eq}}} \]
(This is Eq. 48 in Stergioulas’ thesis)

Draft documentation

For spherical stars, the isotropic radius \( r_{\mathrm{is}} \) is defined by

\[ \frac{d r}{d r_{\mathrm{is}}} = \left(\frac{r}{r_{\mathrm{is}}}\right) \left(1 - 2 \frac{m}{r}\right)^{1/2} \]

Quadrupole moments

Quadrupole moments computed using the method in [Laarakkers99].

Axisymmetric Instability

[Friedman88] shows that a secular axisymmetric instability sets in when the mass becomes maximum along a sequence of constant angular momentum. Equivalently, [Cook92] shows that the instability occurs when the angular momentum becomes minimum along a sequence of constant rest mass.

A GR virial theorem for a stationary and axisymmetric system was found in [Bonazzola73]. A more general two-dimensional virial identity was found in [Bonazzola94]. The three-dimensional virial identity found in [Gourgoulhon94] is a generalization of the Newtonial virial theorem.

Using the stationary and axisymmetric metric ( \( G = c = 1 \) )

\[ ds^2 = - e^{\gamma+\rho} dt^2 + e^{2 \alpha} \left( dr^2 + r^2 d\theta^2 \right) + e^{\gamma-\rho} r^2 \sin^2 \theta ( d \phi - \omega dt) ^2 \]
one solves for the four metric functions \( \rho(r,\theta) \), \( \gamma(r,\theta) \), \( \alpha(r,\theta) \) and \( \omega(r,\theta) \) .

It is assumed that matter is a perfect fluid, and the stress-energy tensor is

\[ T^{\mu \nu} = \left( \rho_0 + \rho_i + P \right) u^{\mu} u^{\nu} + P g^{\mu \nu} \]

Einstein’s field equations imply four field equations for a specified rotation law,

\[ u^{t} u_{\phi} = F(\Omega) \]
for some function \( F(\omega) \) .

Using Eq. (27) in [Cook92], one can write

\[ \rho(s,\mu) = - e^{-\gamma/2} \sum_{n=0}^{\infty} P_{2n}(\mu) \int_0^{1}~ds^{\prime} \int_0^1~d \mu f_{\rho}(n,s,s^{\prime}) P_{2n}{\mu^{\prime}} \tilde{S}(s^{\prime},\mu^{\prime}) \]
where the function \( f_{\rho} \) is defined by
\[ f_{\rho} \equiv \Theta(s^{\prime}-s) \left(\frac{1-s}{s}\right)^{2 n+1} \left[\frac{s^{\prime 2n}}{(1-s^{\prime})^{2n+2}}\right] + \Theta(s^{\prime}-s) \left(\frac{1-s}{s}\right)^{2 n+1} \left[\frac{s^{\prime 2n}}{(1-s^{\prime})^{2n+2}}\right] \]
This function is stored in f_rho . Similar definitions are made for f_gamma and f_omega .

The Keplerial orbit at the equator is

\[ \Omega_K = \frac{\omega^{\prime}}{2 \psi^{\prime}} ... \]

(eq. 31 in [Stergioulas03])

Note

This class is still experimental.

EOS member variables

double eq_radius_tol_rel

Relative accuracy for the equatorial radius, \( r_e \) (default \( 10^{-5} \))

Used in iterate() .

double alt_tol_rel

Accuracy for equatorial radius using alternate solvers (default \( 10^{-9} \))

int verbose

Verbosity parameter.

bool eos_set

If true, then an EOS has been set.

bool scaled_polytrope

If true, then use a polytrope and rescale.

eos_nstar_rot *eosp

Pointer to the user-specified EOS.

nstar_rot()
void output_table(o2scl::table3d &t)

Create an output table3d object from results.

This function creates slices named “ed pr h vsq rho gamma omega alpha” on the grid of \( s \) and \( \mu \).

void calc_masses_J(ubmatrix &rho_0)

Compute masses and angular momentum.

Output

double e_center

Central energy density (in units of \( 10^{15} \mathrm{g}/\mathrm{cm}^3 \))

double r_ratio

Ratio of polar to equatorial radius.

double r_e

Coordinate equatorial radius.

Quantities computed by nstar_rot::comp() (in order)

double r_p

Radius at pole.

double s_p

The value of the s-coordinate at the pole.

double s_e

The value of the s-coordinate at the equator.

double velocity_equator

The velocity at the equator.

double R_e

Circumferential radius in cm (i.e. the radius defined such that \( 2 \pi R_e \) is the proper circumference)

double Mass_p

Proper mass (in \( \mathrm{g} \) )

double Mass

Gravitational mass (in \( \mathrm{g} \) )

double Mass_0

Baryonic mass (in \( \mathrm{g} \) )

double J

Angular momentum (in \( \mathrm{g}~\mathrm{cm}^2 / \mathrm{s} \) )

double Omega

Angular velocity (in \( \mathrm{radians}/\mathrm{s} \) )

double T

Total rotational kinetic energy.

double I

Moment of inertia.

double W

Gravitational binding energy.

double Z_p

Polar redshift.

double Z_f

Forward equatorial redshift.

double Z_b

Backward equatorial redshift.

double Omega_K

Kepler rotation frequency (in 1/s)

double eccentricity

The eccentricity.

ubvector v_plus

Desc.

ubvector v_minus

Desc.

double vel_plus

Desc.

double vel_minus

Desc.

double h_plus

Height from surface of last stable co-rotating circular orbit in equatorial plane.

If this is zero then all orbits are stable.

double h_minus

Height from surface of last stable counter-rotating circular orbit in equatorial plane.

If this is zero then all orbits are stable.

double Omega_plus

Desc.

double u_phi

Desc.

double Omega_p

Angular velocity of a particle in a circular orbit at the equator.

double grv2

Desc.

double grv2_new

Desc.

double grv3

Desc.

double om_over_Om

Ratio of potential \( \omega \) to angular velocity \( \Omega \).

double mass_quadrupole

Mass quadrupole moment.

Settings

double cf

The convergence factor (default 1.0)

Internal constants

double C

Speed of light in vacuum (in CGS units)

double G

Gravitational constant (in CGS units)

double MSUN

Mass of sun (in g)

double KAPPA

Square of length scale in CGS units, \( \kappa \equiv 10^{-15} c^2/G \).

double MB

The mass of one baryon (in g)

double KSCALE

The value \( \kappa G c^{-4} \).

double PI

The constant \( \pi \).

void constants_rns()

Use the values of the constants from the original RNS code.

void constants_o2scl()

Use the O2scl constants (the default)

Grid quantities set in make_grid()

ubvector s_gp

The quantity \( s \), which runs from 0 to SMAX (which defaults to 0.9999) in steps of DS=SMAX/(SDIV-1)

ubvector s_1_s

\( s (1-s) \)

ubvector one_s

\( 1-s \)

ubvector mu

The quantity \( \mu \) which runs from 0 to 1 in steps of DM=1/(MDIV-1)

ubvector one_m2

\( 1-\mu^2 \)

ubvector theta

\( \theta \) defined by \( \mathrm{acos}~\mu \)

ubvector sin_theta

\( \sin \theta = \sqrt{1-\mu^2}\)

Grid values computed in integrate() for spherical_star()

ubvector r_gp

Isotropic radius.

ubvector r_is_gp

Radial coordinate.

ubvector lambda_gp

Metric function \( \lambda \).

ubvector nu_gp

Metric function \( \nu \).

ubvector m_gp

Enclosed gravitational mass.

ubvector e_d_gp

Energy density.

ubmatrix dgds

Desc.

ubmatrix dgdm

Desc.

Metric functions

ubmatrix rho

potential \( \rho \)

ubmatrix gamma

potential \( \gamma \)

ubmatrix omega

potential \( \omega \)

ubmatrix alpha

potential \( \alpha \)

Initial guess computed by the comp() function

double r_e_guess

Guess for the equatorial radius.

ubmatrix rho_guess

Guess for \( \rho \).

ubmatrix gamma_guess

Guess for \( \gamma \).

ubmatrix omega_guess

Guess for \( \alpha \).

ubmatrix alpha_guess

Guess for \( \omega \).

EOS quantities

ubmatrix energy

Energy density \( \epsilon \).

ubmatrix pressure

Pressure.

ubmatrix enthalpy

Enthalpy.

Other quantities defined over the full two-dimensional grid

ubmatrix velocity_sq

Proper velocity squared.

ubmatrix da_dm

Derivative of \( \alpha \) with respect to \( \mu \).

Quantities defined for fixed values of mu

ubvector gamma_mu_1

\( \gamma(s) \) at \( \mu=1 \)

ubvector gamma_mu_0

\( \gamma(s) \) at \( \mu=0 \)

ubvector rho_mu_1

\( \rho(s) \) at \( \mu=1 \)

ubvector rho_mu_0

\( \rho(s) \) at \( \mu=0 \)

ubvector omega_mu_0

\( \omega(s) \) at \( \mu=0 \)

double gamma_pole_h

The value of \( \hat{\gamma} \) at the pole.

double gamma_center_h

The value of \( \hat{\gamma} \) at the center.

double gamma_equator_h

The value of \( \hat{\gamma} \) at the equator.

double rho_pole_h

The value of \( \hat{\rho} \) at the pole.

double rho_center_h

The value of \( \hat{\rho} \) at the center.

double rho_equator_h

The value of \( \hat{\rho} \) at the equator.

double omega_equator_h

The value of \( \hat{\omega} \) at the equator.

double Omega_h

Angular velocity, \( \hat{\omega} \).

double p_center

Central pressure.

double h_center

Central enthalpy.

Helper functions for Green’s function expansion

tensor3 f_rho

\( f_{\rho}(s,n,s') \)

tensor3 f_gamma

\( f_{\gamma}(s,n,s') \)

tensor3 f_omega

\( f_{\omega}(s,n,s') \)

Legendre polynomials

ubmatrix P_2n

Legendre polynomial \( P_{2n}(\mu) \).

ubmatrix P1_2n_1

Associated Legendre polynomial \( P^1_{2n-1}(\mu) \).

ubmatrix D1_rho

Integrated term over m in eqn for \( \rho \).

ubmatrix D1_gamma

Integrated term over m in eqn for \( \gamma \).

ubmatrix D1_omega

Integ. term over m in eqn for \( \omega \).

ubmatrix D2_rho

Integrated term over s in eqn for \( \rho \).

ubmatrix D2_gamma

Integrated term over s in eqn for \( \gamma \).

ubmatrix D2_omega

Integ. term over s in eqn for \( \omega \).

ubmatrix S_gamma

source term in eqn for \( \gamma \)

ubmatrix S_rho

source term in eqn for \( \rho \)

ubmatrix S_omega

source term in eqn for \( \omega \)

double tol_abs

The tolerance for the functions with the prefix “fix” (default \( 10^{-4} \) )

Thermodyanmic quantities near the surface

double p_surface

Pressure at the surface.

double e_surface

Energy density at the surface.

double enthalpy_min

Minimum specific enthalpy.

Polytrope parameters

double n_P

Polytropic index.

double Gamma_P

Polytropic exponent.

Interpolation functions

int n_nearest

Cache for interpolation.

int new_search(int n, ubvector &x, double val)

Search in array x of length n for value val.

double interp(ubvector &xp, ubvector &yp, int np, double xb)

Driver for the interpolation routine.

First we find the tab. point nearest to xb, then we interpolate using four points around xb.

Used by int_z(), e_at_p(), p_at_e(), p_at_h(), h_at_p(), n0_at_e(), comp_omega(), comp_M_J(), comp(), spherical_star(), iterate().

double interp_4_k(ubvector &xp, ubvector &yp, int np, double xb, int k)

Driver for the interpolation routine.

Four point interpolation at a given offset the index of the first point k.

Used in comp() .

double int_z(ubvector &f, int m)

Integrate f[mu] from m-1 to m.

This implements a 8-point closed Newton-Cotes formula.

Used in comp() .

Basic Usage

inline void set_eos(eos_nstar_rot &eos)

Set the EOS.

inline void polytrope_eos(double index)

Use a polytropic EOS with a specified index.

int fix_cent_eden_axis_rat(double cent_eden, double axis_rat, bool use_guess = false)

Construct a configuration with a fixed central energy density and a fixed axis ratio.

The central energy density should be in \( \mathrm{g}/\mathrm{cm}^3 \) and the axis ratio is unitless. This is fastest of the high-level interface functions as it doesn’t require an additional solver.

int fix_cent_eden_grav_mass(double cent_eden, double grav_mass)

Construct a configuration with a fixed central energy density and a fixed gravitational mass.

The central energy density should be in \( \mathrm{g}/\mathrm{cm}^3 \) and the gravitational mass should be in solar masses.

int fix_cent_eden_bar_mass(double cent_eden, double bar_mass)

Construct a configuration with a fixed central energy density and a fixed baryonic mass.

The central energy density should be in \( \mathrm{g}/\mathrm{cm}^3 \) and the baryonic mass should be in solar masses.

int fix_cent_eden_with_kepler(double cent_eden)

Construct a configuration with a fixed central energy density and the Keplerian rotation rate.

The central energy density should be in \( \mathrm{g}/\mathrm{cm}^3 \) .

int fix_cent_eden_with_kepler_alt(double cent_eden, bool use_guess = false)

Experimental alternate form for fix_cent_eden_with_kepler()

int fix_cent_eden_grav_mass_alt(double cent_eden, double grav_mass, bool use_guess = false)

Experimental alternate form for fix_cent_eden_grav_mass()

int fix_cent_eden_bar_mass_alt(double cent_eden, double bar_mass, bool use_guess = false)

Experimental alternate form for fix_cent_eden_bar_mass()

int fix_cent_eden_ang_vel_alt(double cent_eden, double ang_vel, bool use_guess = false)

Experimental alternate form for fix_cent_eden_ang_vel()

int fix_cent_eden_ang_mom_alt(double cent_eden, double ang_mom, bool use_guess = false)

Experimental alternate form for fix_cent_eden_ang_mom()

int fix_cent_eden_non_rot(double cent_eden)

Construct a non-rotating configuration with a fixed central energy density.

The central energy density should be in \( \mathrm{g}/\mathrm{cm}^3 \) .

int fix_cent_eden_ang_vel(double cent_eden, double ang_vel)

Construct a configuration with a fixed central energy density and a fixed angular velocity.

The central energy density should be in \( \mathrm{g}/\mathrm{cm}^3 \) and the angular velocity should be in \( \mathrm{rad}/\mathrm{s} \). The final angular velocity (possibly slightly different than ang_vel is stored in Omega .

Note

In the original RNS code, the ang_vel argument is different because it was rescaled by a factor of \( 10^{4} \).

int fix_cent_eden_ang_mom(double cent_eden, double ang_mom)

Construct a configuration with a fixed central energy density and a fixed angular momentum.

The central energy density should be in \( \mathrm{g}/\mathrm{cm}^3 \). The angular momentum should be in units of \( G M_{\odot}^2/C \) .

Testing functions

These functions compare this class with hard-coded results obtained with the RNS code.

void test1(o2scl::test_mgr &t)

Test determining configuration with fixed central energy density and fixed radius ratio with EOS C.

void test2(o2scl::test_mgr &t)

Test configuration rotating and Keplerian frequency with a fixed central energy density and EOS C.

void test3(o2scl::test_mgr &t)

Test fixed central energy density and fixed gravitational mass with EOS C.

void test4(o2scl::test_mgr &t)

Test fixed central energy density and fixed baryonic mass with EOS C.

void test5(o2scl::test_mgr &t)

Test fixed central energy density and fixed angular velocity with EOS C.

void test6(o2scl::test_mgr &t)

Test fixed central energy density and fixed angular momentum with EOS C.

void test7(o2scl::test_mgr &t)

Test a series of non-rotating stars on a energy density grid with EOS C.

void test8(o2scl::test_mgr &t)

Test Keplerian frequency for a polytrope.

EOS functions

double e_at_p(double pp)

Compute \( \varepsilon(P) \)

Used in dm_dr_is(), dp_dr_is(), integrate() and iterate().

double p_at_e(double ee)

Compute \( P(\varepsilon) \)

Used in make_center() and integrate().

double p_at_h(double hh)

Pressure at fixed enthalpy.

Used in iterate().

double h_at_p(double pp)

Enthalpy at fixed pressure.

Used in make_center() and integrate().

double n0_at_e(double ee)

Baryon density at fixed energy density.

Used in comp_M_J() and comp() .

Derivatives on the grid

double s_deriv(ubvector &f, int s)

Returns the derivative w.r.t. s of an array f[SDIV+1].

double m_deriv(ubvector &f, int m)

Returns the derivative w.r.t. mu of an array f[MDIV+1].

double deriv_s(ubmatrix &f, int s, int m)

Returns the derivative w.r.t. s

double deriv_m(ubmatrix &f, int s, int m)

Returns the derivative w.r.t. mu.

double deriv_sm(ubmatrix &f, int s, int m)

Returns the derivative w.r.t. s and mu.

Initialization functions

double legendre(int n, double x)

Returns the Legendre polynomial of degree n, evaluated at x.

This uses the recurrence relation and is used in comp_f_P() which is called by the constructor.

void comp_f_P()

Compute two-point functions.

This function computes the 2-point functions \( f^m_{2n}(r,r') \) used to integrate the potentials \( \rho, \gamma \) and \( \omega \) (See Komatsu et al. 1989 for details). Since the grid points are fixed, we can compute the functions f_rho, f_gamma, f_omega, P_2n, and P1_2n_1 once at the beginning.

See Eqs. 27-29 of [Cook92] and Eqs. 33-35 of [Komatsu89]. This function is called by the constructor.

void make_grid()

Create computational mesh.

Create the computational mesh for \( s=r/(r+r_e) \) (where \( r_e \) is the coordinate equatorial radius) and \( \mu = \cos \theta \) using

\[ s[i]=\mathrm{SMAX}\left(\frac{i-1}{\mathrm{SDIV}-1}\right) \]
\[ \mu[j]=\left(\frac{i-1}{\mathrm{MDIV}-1}\right) \]
When \( r=0 \), \( s=0 \), when \( r=r_e \), \( s=1/2 \), and when \( r = \infty \), \( s=1 \) . Inverting the relationship between \( r \) and \( s \) gives \( r = r_e s / (1-s) \) .

Points in the mu-direction are stored in the array mu[i]. Points in the s-direction are stored in the array s_gp[j].

This function sets s_gp, s_1_s, one_s, mu, one_m2, theta and sin_theta . All of these arrays are unit-indexed. It is called by the constructor.

void make_center(double e_center)

Compute central pressure and enthalpy from central energy density.

For polytropic EOSs, this also computes rho0_center .

Post-processing functions

void comp_omega()

Compute Omega and Omega_K.

void comp_M_J()

Compute rest mass and angular momentum.

void comp()

Compute various quantities.

The main post-processing function

For computing spherical stars

void spherical_star()

Computes a spherically symmetric star.

The metric is

\[ ds^2 = -e^{2\nu}dt^2 + e^{2\lambda} dr^2 + r^2 d\theta^2 + r^2 sin^2\theta d\phi^2 \]
where \( r \) is an isotropic radial coordinate (corresponding to r_is in the code).

Todo:

AWS: 7/20/20: Better document out how this metric definition leads to \( \gamma=\nu+\lambda \) and \( \rho = \nu - \lambda \) and the relationship between r and r_is .

This function computes r_e_guess, R_e, Mass, and Z_p .

double dm_dr_is(double r_is, double r, double m, double p)

Derivative of gravitational mass with respect to isotropic radius.

double dp_dr_is(double r_is, double r, double m, double p)

Derivative of pressure with respect to isotropic radius.

double dr_dr_is(double r_is, double r, double m)

Derivative of radius with respect to isotropic radius.

void integrate(int i_check, double &r_final, double &m_final, double &r_is_final)

Integrate one of the differential equations for spherical stars.

Desc

int iterate(double r_ratio, double tol_rel)

Main iteration function.

Public Types

typedef boost::numeric::ublas::vector<double> ubvector
typedef boost::numeric::ublas::range ub_range
typedef boost::numeric::ublas::vector_range<boost::numeric::ublas::vector<double>> ubvector_range
typedef boost::numeric::ublas::matrix<double> ubmatrix

Public Functions

void resize(int MDIV_new, int SDIV_new, int LMAX_new, int RDIV_new)

Resize the grid.

inline int set_solver(o2scl::mroot<> &m)

Set new solver.

Public Members

int MDIV

The number of grid points in the \( \mu \) direction.

int SDIV

The number of grid points in the \( s \) direction.

int LMAX

The number of Legendre polynomials.

o2scl::mroot_hybrids def_mroot

Default solver.

Protected Functions

int solve_kepler(size_t nv, const ubvector &x, ubvector &y)

Solve for the Keplerian velocity.

int solve_grav_mass(size_t nv, const ubvector &x, ubvector &y, double grav_mass)

Solve for the gravitational mass.

int solve_bar_mass(size_t nv, const ubvector &x, ubvector &y, double bar_mass)

Solve for the gravitational mass.

int solve_ang_vel(size_t nv, const ubvector &x, ubvector &y, double ang_vel)

Solve for the gravitational mass.

int solve_ang_mom(size_t nv, const ubvector &x, ubvector &y, double ang_mom)

Solve for the gravitational mass.

Protected Attributes

o2scl::mroot *mrootp

Solver.

o2scl::root_bkt_cern<polytrope_solve> rbc

The polytrope solver.

o2scl::search_vec<ubvector_range> sv_ub

Array search object.

int RDIV

The number of grid points in integration of TOV equations for spherical stars.

double SMAX

Maximum value of s-coordinate (default 0.9999)

double DS

Spacing in \( s \) direction, \( \mathrm{SMAX}/(\mathrm{SDIV}-1) \).

double DM

Spacing in \( \mu \) direction, \( 1/(\mathrm{MDIV}-1) \).

double RMIN

Minimum radius for spherical stars (default \( 10^{-15} \))

class polytrope_solve

Subclass of nstar_rot which specifies the function to invert a polytropic EOS.

Public Functions

inline polytrope_solve(double Gamma_P, double ee)

Create a function object with specified polytropic index and ?

inline double operator()(double rho0)

The function.

Protected Attributes

double _Gamma_P

The polytropic index.

double _ee

The energy density.