# Solution of the Tolman-Oppenheimer-Volkov Equations¶

O2scl

The class tov_solve provides a solution to the Tolman-Oppenheimer-Volkov (TOV) equations given an equation of state (EOS), provided as an object of type eos_tov. These classes are particularly useful for static neutron star structure: given any equation of state one can calculate the mass vs. radius curve and the properties of any star of a given mass.

The EOS can be specified as an object of type eos_tov. The documentation of this parent class contains more information. The class eos_tov_interp is used most frequently. It uses linear interpolation to interpolate a user-specified table object. A faster lower-level EOS interpolation is performed by eos_tov_vectors. The Buchdahl EOS is given in eos_tov_buchdahl, a single polytrope EOS is given in eos_tov_polytrope, a linear EOS is given in eos_tov_linear, and an EOS with a polynomial form for the speed of sound is given in eos_cs2_poly

In units where $$c=1$$, the most general static and spherically symmetric metric is of the form

$ds^2 = - e^{2 \Phi(r)} d t^2 + e^{2 \Lambda(r)} d r^2 + r^2 d \theta^2 + r^2 \sin^2 \theta~d \phi^2$

where $$\theta$$ is the polar angle and $$\phi$$ is the azimuthal angle. Often we will not write explicitly the radial dependence for many of the quantities defined below, i.e. $$\Phi \equiv \Phi(r)$$.

This leads to the TOV equation (i.e. Einstein’s equations for a static and spherically symmetric star)

$\frac{d P}{d r} = - \frac{G \varepsilon m}{r^2} \left( 1+\frac{P}{\varepsilon}\right) \left( 1+\frac{4 \pi P r^3}{m} \right) \left( 1-\frac{2 G m}{r}\right)^{-1}$

where $$r$$ is the radial coordinate, $$m(r)$$ is the gravitational mass enclosed within a radius $$r$$, and $$\varepsilon(r)$$ and $$P(r)$$ are the energy density and pressure at $$r$$, and $$G$$ is the gravitational constant. The mass enclosed, $$m(r)$$, is related to the energy density through

$\frac{d m}{d r} = 4 \pi r^2 \varepsilon$

and these two differential equations can be solved simultaneously given an equation of state, $$P(\varepsilon)$$. The total gravitational mass is given by

$M = \int_0^R 4 \pi r^2 \varepsilon d r$

The boundary conditions are $$m(r=0)=0$$ and the condition $$P(r=R)=0$$ for some fixed radius $$R$$. These boundary conditions give a one-dimensional family solutions to the TOV equations as a function of the central pressure. Each central pressure implies a gravitational mass, $$M$$, and radius, $$R$$, and thus defines a mass-radius curve.

The metric function $$\Lambda$$ is

$e^{2 \Lambda} = \left( 1-\frac{2 G m}{r}\right)^{-1}$

The other metric function, $$\Phi(r)$$ is sometimes referred to as the gravitational potential. In vacuum above the star, it is

$e^{2 \Phi} = \left( 1-\frac{2 G M}{r}\right)$

and inside the star it is determined by

$\frac{d \Phi}{d r} = - \frac{1}{\varepsilon} \frac{ d P}{d r} \left(1+\frac{P}{\varepsilon}\right)^{-1} = \frac{G m}{r^2} \left( 1+\frac{4 \pi P r^3}{m} \right) \left( 1-\frac{2 G m}{r}\right)^{-1}$

Alternatively, that this can be rewritten as

$-d \Phi = \frac{d P}{P+\varepsilon} \, .$

In this form, $$\Phi$$ has no explicit dependence on $$r$$ so it can be computed (up to a constant) directly from the EOS.

If the neutron star is at zero temperature and there is only one conserved charge, (i.e. baryon number), then

$-d \Phi = \frac{d P}{\mu n} = \frac{d \mu}{\mu}$

and this implies that $$\mu e^{\Phi}$$ is everywhere constant in the star. If one defines the “enthalpy” by

$d h = \frac{dP}{P + \varepsilon}$

then

$-d \Phi = dh$

and thus $$\mu \propto e^{h}$$ or $$h = \ln \mu + C$$. This is the enthalpy used by the nstar_rot class.

Keep in mind that this enthalpy is determined by integrating the quantities in the stellar profile (which may be, for example, in beta-equilibrium). Thus, this is not equal the usual thermodynamic enthalpy which is

$H(P,S,N) = E + P V = T S + \sum_i \mu_i N_i$

or in differential form

$d H = T dS + V dP + \sum_i \mu_i d N_i \, .$

The proper boundary condition for the gravitational potential is

$\Phi(r=R) = \frac{1}{2} \ln \left(1-\frac{2 G M}{R} \right)$

which ensures that $$\Phi(r)$$ matches the metric above in vacuum. Since the expression for $$d\Phi/dr$$ is independent of $$\Phi$$, the differential equation can be solved for an arbitrary value of $$\Phi(r=0)$$ and then shifted afterwards to obtain the correct boundary condition at $$r=R$$ .

The surface gravity is defined to be

$g = \frac{G m}{r^2}\left(1-\frac{2 G m}{r}\right)^{-1/2}$

which is computed in units of inverse kilometers, and the redshift is defined to be

$z = \left(1-\frac{2 G m}{r}\right)^{-1/2} - 1$

which is unitless.

The baryonic mass is typically defined by

$M_B = \int_0^R 4 \pi r^2 n_B m_B \left(1-\frac{2 G m}{r}\right)^{-1/2} d r$

where $$n_B(r)$$ is the baryon number density at radius $$r$$ and $$m_B$$ is the mass one baryon (taken to be the mass of the proton by default and stored in o2scl::tov_solve::baryon_mass). If the EOS specifies the baryon density (i.e. if o2scl::eos_tov::baryon_column is true), then tov_solve will compute the associated baryonic mass for you.