Relativistic Hydrostatic Structure Equations
and Analytic Multilayer Stellar Model

Before moving on to detailed applications of the proposed relativistic hydrostatic structure equations, it is convenient to investigate how to determine local temperature from them. From the definition of $\check{E}_{r}$ and (2.5), a concise expression of $\check{E}_{r}$ with respect to thermodynamic variables is derived as

$$\check{E}_{r}=\frac{\frac{c^{2}r^{2}\nu^{\prime}}{G}-8\pi r^{3}P}{2(1+r\nu^{\prime})},$$

(2.12)

with $\nu^{\prime}=-2T^{\prime}/T$.²²2 The leading expression in (2.12), $\check{E}_{r}\approx\frac{c^{2}r^{2}\nu^{\prime}}{2G}$, may be reminiscent of the conventional mass expression of a star derived by incorporating the polytropic relation $P\propto\varrho^{1+1/n}$ and the Lane-Emden equation $\frac{1}{\xi^{2}}{d\over d\xi}(\xi^{2}{d\theta\over d\xi})=-\theta^{n}$ with $\xi\propto r$ as $m\propto\xi_{1}^{2}|\theta(\xi_{1})^{\prime}|$, where $\theta(\xi_{1})=0$ [32]. However, the expression (2.12) is much more general in the sense to hold without imposing any equation of state and to include the effects not only of general relativity but also of nonzero temperature. Differentiating both sides with respect to $r$ and using (2.3) and (2.7), one finds

$\displaystyle\nu^{\prime\prime}+\frac{2\nu^{\prime}}{r}+\nu^{\prime 2}-\frac{4\pi G}{c^{2}}\left((r\nu^{\prime}+1)(r\nu^{\prime}+2)\rho-\left(2r^{2}\nu^{\prime\prime}+r\nu^{\prime}(r\nu^{\prime}-3)-6\right)P\right)=0.$

(2.13)

This can be rewritten into a general coordinate invariant form as

$\displaystyle\nabla^{2}\nu=\frac{8\pi G}{c^{2}}(g^{\omega\mu}+2u^{\omega}u^{\mu})T_{\omega\mu},$

(2.14)

which is the general relativistic extension of the Poisson equation. To prove this, first compute the Laplacian acting on the scalar field $\nu$ as

	$\displaystyle\nabla^{2}\nu=$	$\displaystyle e^{-\lambda}(\nu^{\prime\prime}+(\frac{2}{r}+\frac{\nu^{\prime}-\lambda^{\prime}}{2})\nu^{\prime})$
	$\displaystyle=$	$\displaystyle e^{-\lambda}\nu^{\prime\prime}+\frac{\nu^{\prime}+4\pi Gr/c^{2}(P(r\nu^{\prime}+3)-(r\nu^{\prime}+1)\rho)+2/r}{r\nu^{\prime}+1}\nu^{\prime},$		(2.15)

where at the 1st equation the formula $\nabla^{2}\nu=\frac{1}{\sqrt{|g|}}\partial_{\omega}\sqrt{|g|}g^{\omega\mu}\partial_{\mu}\nu$ was used and at the 2nd one the term containing $\lambda^{\prime}$ was computed by employing (2.4) and (2.5) as

$\displaystyle e^{-\lambda}\frac{\nu^{\prime}-\lambda^{\prime}}{2}=4\pi Gr\frac{P-\rho}{c^{2}}+\frac{2G\check{E}_{r}}{r^{2}c^{2}},$

(2.16)

and (2.12) was used to substitute $\check{E}_{r}$. The coefficient of $\nu^{\prime\prime}$ in (2.15) is $e^{-\lambda}=\frac{1+8\pi Gr^{2}P/c^{2}}{r\nu^{\prime}+1}$, while that in (2.13) is $1+8\pi Gr^{2}P/c^{2}=e^{-\lambda}(r\nu^{\prime}+1)$. Thus multiplying $r\nu^{\prime}+1$ for both sides in (2.13) and subtracting it from (2.15), one can remove $\nu^{\prime\prime}$ in (2.15) and the rest is [22]

$$\nabla^{2}\nu=8\pi G(\rho+3P)/c^{2}.$$

(2.17)

This matches (2.14) and completes the proof. It is important to stress that the form of the relativistic Poisson equation (2.14) or (2.17) has been derived here without using any approximation and holds non-perturbatively in the Newton constant.

The generalized relativistic Poisson equation can be solved order by order in the Newton constant. The leading order is

$$\nu^{\prime\prime}+\frac{2\nu^{\prime}}{r}+\nu^{\prime 2}=0.$$

(2.18)

This has a trivial solution that $\nu$ is constant, while there is a non-trivial one as $\nu^{\prime}=\frac{1}{r(-1+C_{0}r)}=:\nu_{\star}^{\prime}$, where $C_{0}$ is an integration constant. This non-trivial solution corresponds to the exact one of the original differential equation (2.13) in the matter-empty region, $P=\rho=0$. The integration constant $C_{0}$ is fixed by plugging this back into (2.12)

$\displaystyle\check{E}_{r}=\frac{c^{2}-8\pi Gr^{2}(-1+C_{0}r)P}{2C_{0}G}.$

(2.19)

Then evaluate this at a point $r=R_{\star}$ so far away from the center of the system that the contribution of the pressure to $\check{E}_{r}$ is negligible, $\check{E}_{r}\approx\frac{c^{2}}{2C_{0}G}$ with $r\geq R_{\star}$. In this far region, $\check{E}_{r}$ is approximated as a constant and known as the gravitational mass by dividing it by $c^{2}$, which is denoted by $M_{\star}$. Then the integration constant is fixed as $C_{0}=\frac{1}{2M_{\star}G}$ and $\nu_{\star}$ is given by $\nu_{\star}=\log\left(\frac{1}{r}-\frac{1}{2M_{\star}G}\right)+\tilde{C}_{0}$, where $\tilde{C}_{0}$ is another integration constant. This integration constant originates in the arbitrariness of the base point of the gravitational potential, which can be fixed by choosing an arbitrary referencing point. Here it is fixed by requesting it to vanish at the evaluation point: $\nu(R_{\star})=0$. Then the gravitational potential in the matter-empty region is determined as

$$\nu_{\star}=\log\left(\frac{1-\frac{2M_{\star}G}{R_{\star}}}{1-\frac{2M_{\star}G}{r}}\right).$$

(2.20)

This exact gravitational potential in the matter-empty region is expanded with respect to the Newton constant as

$$\nu_{\star}=2M_{\star}G\left(-\frac{1}{r}+\frac{1}{R_{\star}}\right)+\cdots.$$

(2.21)

This leading term describes the Newton’s law of universal gravitation for an object of the mass $M_{\star}$, which is the origin of the name of $M_{\star}$ mentioned above. As a result the evaluation point $r=R_{\star}$ should be chosen to satisfy the pressure at the point $P(R_{\star})$ to satisfy an inequality

$$c^{2}\gg 4\pi R_{\star}^{3}\frac{P(R_{\star})}{M_{\star}},$$

(2.22)

with $R_{\star}\gg r_{\rm H}$ assumed. Accordingly the behavior of temperature in the matter-empty region is exactly determined as [23]

$\displaystyle T=T(R_{\star})e^{-\frac{\nu}{2}}=T(R_{\star})\sqrt{\frac{1-\frac{r_{\rm H}}{R_{\star}}}{1-\frac{r_{\rm H}}{r}}}=:T_{\rm vac},$

(2.23)

where $r_{\rm H}:={2M_{\star}G}$ is the radius of the horizon of the (Schwarzschild) black hole whose mass is identical to the stellar gravitational mass. The local temperature is formally divergent at the location of the horizon radius. It can be confirmed that (2.23) satisfies the differential equation $T_{\rm vac}\left(T_{\rm vac}^{\prime\prime}+\frac{2}{r}T_{\rm vac}^{\prime}\right)-3T_{\rm vac}^{\prime 2}=0$, which is identical to the leading order of the relativistic Poisson equation written in term of the local temperature as

		$\displaystyle T\left(T^{\prime\prime}+\frac{2}{r}T^{\prime}\right)-3T^{\prime 2}$
	$\displaystyle=$	$\displaystyle\frac{4\pi G}{c^{2}}\left(\rho\left(-2r^{2}T^{\prime 2}+3rTT^{\prime}-T^{2}\right)+P\left(4r^{2}T^{\prime 2}-rT\left(2rT^{\prime\prime}-3T^{\prime}\right)-3T^{2}\right)\right).$		(2.24)

This is the relativistic steady-state heat conduction equation for a relativistic hydrostatic equilibrium system with rotational symmetry.

The next-to-leading order of the general relativistic Poisson equation in the Newton constant is determined by expanding the differential equation (2.13) as $\nu=\nu_{0}+G\nu_{1}+{\cal O}(G^{2})$ with $\nu_{0}$ constant and taking the leading term as

$$\nu_{1}^{\prime\prime}+\frac{2}{r}\nu_{1}^{\prime}=\frac{8\pi}{c^{2}}(\rho_{0}+3P_{0}),$$

(2.25)

where $\rho_{0},P_{0}$ are the leading order of the energy density and that of the pressure, respectively. In terms of the variable of temperature, the next-to-leading order of the relativistic Poisson equation is given by

$$T_{1}^{\prime\prime}+\frac{2}{r}T_{1}^{\prime}=-\frac{4\pi T_{0}}{c^{2}}(\rho_{0}+3P_{0}),$$

(2.26)

where $T=T_{0}+GT_{1}+\cdots$ with $T_{0}$ constant. Inputing detailed information on fluid and solving the differential equation (2.25), (2.26) gives the subleading correction of the gravitational potential, the local temperature, respectively.

Comments are given below. Firstly, a perturbative solution for (2.13) or (2.24) obtained this way is useful in a weak-gravity regime with energy density and pressure sufficiently small such as the vicinity of a stellar surface, while in a strong-gravity region such as the interior of a heavy stellar core, more useful information will be obtained by solving the original relativistic Poisson equation directly in a non-perturbative fashion with respect to the Newton constant. Secondly, astrophysical observation is basically done in a region far away from the core, so $r\gg r_{\rm H}$. Therefore, the condition $r\gg 2G\check{E}_{r}/c^{2}$ is assumed for practical applications in what follows. In particular, this condition is satisfied in the Newtonian regime mentioned above.

To the end of the construction of a model of luminous stars, first consider a hydrostatic equilibrium system of a single ideal gas of particles with conserved particle number current, which are called baryonic particles in this paper. The equation of state of an ideal gas is

$$P=\hat{n}k_{\text{B}}T,$$

(3.1)

where $\hat{n}v$ gives the ordinary particle number density and $k_{\text{B}}$ is the Boltzmann constant, while the conservation of the particle number current requires that $\nabla_{\mu}J^{\mu}=0$, where $J^{\mu}=\hat{n}u^{\mu}$. The left hand side is computed in a radially moving frame as $\nabla_{\mu}J^{\mu}=u^{r}(\hat{n}^{\prime}+\hat{n}(\log(\sqrt{|g|}u^{r}))^{\prime})$. A relativistic fluid equation derived in [23] rewrites this as $\nabla_{\mu}J^{\mu}=u^{r}(\hat{n}^{\prime}-\hat{n}\frac{\rho^{\prime}}{\rho+P})$. Therefore the conservation of the number current leads to

$$\hat{n}^{\prime}=\hat{n}\frac{\rho^{\prime}}{\rho+P}.$$

(3.2)

Substituting (3.1) into this so as to eliminate $\hat{n}$ and using (2.9), one finds $P^{\prime}/P=\rho^{\prime}/\rho$. This can be solved as

$$P=w\rho,$$

(3.3)

where $w$ is an integration constant. The fluid for pressure to be proportional to energy density such as (3.3) is called a simple fluid in what follows for convenience. To fix the integration constant, recall that the energy density of an ideal gas is given by $u=nC_{V}T$, where $C_{V}$ is the heat capacity at constant volume and $n=\hat{n}v$. This can be rewritten as $\rho=\hat{n}C_{V}T$. By substituting this and (3.1) into (3.3), the integration constant can be fixed as

$$w=\frac{k_{\text{B}}}{C_{V}}=\gamma-1,$$

(3.4)

where $\gamma=C_{P}/C_{V}$ is the heat capacity ratio. Here $C_{P}$ is the heat capacity at constant pressure and given by the Mayer’s relation as $C_{P}=C_{V}+k_{\text{B}}$, which is guaranteed to hold by the laws of thermodynamics confirmed above. Plugging (3.3) into (2.9) with (3.4), one obtains

$$\frac{T^{\prime}}{T}=(1-\frac{1}{\gamma}){P^{\prime}\over P}.$$

(3.5)

This is exactly the saturation point of the inequality of the Schwarzschild stability condition for a hydrostatic equilibrium system [33] and known as an equation of state characteristic to a convection zone, in which the energy transport occurs mainly by convection. On the other hand, in order for fluid in hydrostatic equilibrium to satisfy (3.5), it needs to be a simple fluid satisfying (3.3). As a result, the necessary and sufficient condition for fluid to satisfy the equation (3.5) is that there is a proportionality relationship between pressure and energy density of fluid in the region. Note that the equation (3.5) can be easily solved as $P\propto T^{\frac{\gamma}{\gamma-1}}.$ Thus $\rho=aT^{\frac{\gamma}{\gamma-1}}$, where $a$ is a constant.

This new derivation of the temperature gradient equation in the convection zone leads to a new important prediction for a general property of stars. That is, in a region with the Newtonian approximation valid, pressure is negligible compared to energy density so that the parameter $w$ almost vanishes [16]. Combining this with the observational fact that the Newtonian approximation is valid near the surface of most stars, one concludes that the heat capacity ratio in any stellar convection zone with the Newtonian approximation valid is almost one:³³3 This is consistent with a result in fluid dynamics on flat spacetime that an isentropic fluid satisfies a polytropic relation $P\propto\rho^{\gamma}$.

$$\gamma\sim 1.$$

(3.6)

In particular, it is known that there exists considerably large convection zone underneath the solar surface [34, 35]. Therefore, it follows that fluid near the solar surface has the heat capacity ratio nearly equal to $1$. This is a robust prediction newly made in this paper. The deviation of the heat capacity ratio from one will be described by stellar global observables as (3.28) later.

Next, consider a system of an ideal gas consisting of non-relativistic particles for later use. The equation of state is given by $P=P_{\text{g}},\rho=\rho_{\text{g}}$, where

$\displaystyle P_{\text{g}}=\hat{n}_{\text{g}}k_{\text{B}}T,~{}~{}~{}\rho_{\text{g}}=(\bar{m}c^{2}+\frac{3}{2}k_{\text{B}}T)\hat{n}_{\text{g}},$

(3.11)

with $\bar{m}$ the averaged mass of the non-relativistic particles. This averaged mass is described as $\bar{m}=\mu_{\rm e}m_{\rm u}$, where $\mu_{\rm e}$ is the mean molecular weight per free electron and $m_{\rm u}$ is the atomic mass unit defined by $1{\rm g/mol}$. Note that $\mu_{\rm e}$ can be described by using the weight fraction of hydrogen denoted by $X$ as $\mu_{\rm e}=2/(1+X)$. Then consider to satisfy the temperature gradient equation (2.9), which is in the current setup given by $(P_{\text{g}}+\rho_{\text{g}})T^{\prime}=P^{\prime}_{\text{g}}T$. Solving this differential equation, one can fix the form of the particle number density as

$$\hat{n}_{\text{g}}=C_{\text{g}}T^{3/2}e^{-\frac{\bar{m}c^{2}}{k_{\text{B}}T}},$$

(3.12)

where $C_{\text{g}}$ is an integration constant. This result is indeed consistent with results of statistical physics, and the integration constant can be fixed to match the Maxwell distribution as $C_{\text{g}}=\bar{g}(\frac{2\pi\bar{m}k_{\text{B}}}{h^{2}})^{\frac{3}{2}}$, where $\bar{g}$ corresponds to the mean internal degrees of freedom of a non-relativistic particle.

As previously, the radial dependence of temperature is determined by solving the steady-state heat conduction equation. Differently from the previous case of simple fluid, it is very difficult to solve it exactly, so here it is solved up to the next-to-leading order. As shown in section 2.1, the leading order solution for (2.18) is a constant, $\nu_{0}=\mathrm{const}$, so is the temperature, $T_{0}=\mathrm{const}$, and thus the above thermodynamic quantities as well:

$$P_{0}=\hat{n}_{g0}k_{\text{B}}T_{0},~{}~{}~{}\rho_{0}=(\bar{m}c^{2}+\frac{3}{2}k_{\text{B}}T_{0})\hat{n}_{g0},~{}~{}~{}\hat{n}_{g0}=\bar{g}(\frac{2\pi\bar{m}k_{\text{B}}}{h^{2}})^{\frac{3}{2}}T_{0}^{3/2}e^{-\frac{\bar{m}c^{2}}{k_{\text{B}}T_{0}}}.$$

(3.13)

Therefore the right-hand side of the next-to-leading order of the relativistic Poisson equation (2.25) becomes a constant, so it can be solved as

$$\nu_{1}=\frac{4\pi r^{2}(\rho_{0}+3P_{0})}{3c^{2}}-\frac{C_{1}}{r}+\tilde{C}_{1},$$

(3.14)

where $C_{1},\tilde{C}_{1}$ are integration constants. Plugging this back into (2.12) leads to

$$\check{E}_{r}=\frac{8\pi c^{2}r^{4}\rho_{0}-3c^{4}C_{1}r}{16\pi Gr^{3}(\rho_{0}+3P_{0})+6c^{2}(r+C_{1}G)}.$$

(3.15)

Evaluating this at the surface of the inner layer, $r=R_{*}$, one finds

$$c^{2}M_{*}=\frac{8\pi c^{2}R_{*}^{4}\rho_{0}-3c^{4}C_{1}R_{*}}{16\pi GR_{*}^{3}(\rho_{0}+3P_{0})+6c^{2}(R_{*}+C_{1}G)},$$

(3.16)

where $M_{*}=\check{E}_{r}(R_{*})/c^{2}.$ From this, the integration constant $C_{1}$ is fixed as

	$\displaystyle C_{1}=$	$\displaystyle\frac{2R_{*}M_{*}\left(8\pi GR_{*}^{2}(\rho_{0}+3P_{0})/c^{2}+3\right)-8\pi R_{*}^{4}\rho_{0}/c^{2}}{3R_{*}-6GM_{*}}$
	$\displaystyle\approx$	$\displaystyle 2M_{*}-\frac{8}{3}\pi R_{*}^{3}\rho_{0}/c^{2},$		(3.17)

while $\tilde{C}_{1}$ is determined to request $\nu_{1}$ to vanish at $r=R_{\star}$ as introduced in section 2.1, outside which matter contribution is negligible so that (2.22) is satisfied:

$$\tilde{C}_{1}=\frac{2M_{*}}{R_{\star}}-\frac{4\pi\left(3P_{0}R_{\star}^{3}+\rho_{0}\left(R_{\star}^{3}+2R_{*}^{3}\right)\right)}{3c^{2}R_{\star}}.$$

(3.18)

Substituting this into (3.14) yields the solution of the gravitational potential at the sub-leading order as $\nu_{1}\approx\nu_{*}$, where

$\displaystyle\nu_{*}=2(\frac{1}{R_{\star}}-\frac{1}{r})\left(M_{*}+\frac{6\pi P_{0}rR_{\star}(r+R_{\star})+2\pi\rho_{0}\left(r^{2}R_{\star}+rR_{\star}^{2}-2R_{*}^{3}\right)}{3c^{2}}\right).$

(3.19)

This solution is valid in an annuls region with $R_{*}\leq r\leq R_{\star}$ and should connect smoothly to the one in the matter-empty region (2.21) at $r=R_{\star}$, which imposes a condition to satisfy a relation between the stellar gravitational mass $M_{\star}$ and $M_{*}$ as

$$M_{\star}=M_{*}+4\pi R_{\star}^{3}\frac{P_{0}}{c^{2}}+\frac{4}{3}\pi\left(R_{\star}^{3}-R_{*}^{3}\right)\frac{\rho_{0}}{c^{2}},$$

(3.20)

at the leading order of the Newton constant. Then the local temperature is determined up to the next-to-leading order in $G$ as

$$T=T(R_{\star})(1-G\nu_{*}/2+\cdots),$$

(3.21)

since $T_{0}=T(R_{\star})$. In particular, evaluating at $r=R_{*}$, one can compute the ratio of the temperatures at $r=R_{\star}$ and $r=R_{*}$ as

$\displaystyle\frac{T(R_{*})}{T(R_{\star})}=1+G(\frac{1}{R_{*}}-\frac{1}{R_{\star}})\left(M_{*}+\frac{6\pi P_{0}R_{*}R_{\star}(R_{*}+R_{\star})+2\pi\rho_{0}\left(R_{*}^{2}R_{\star}+R_{*}R_{\star}^{2}-2R_{*}^{3}\right)}{3c^{2}}\right).$

(3.22)

Note that the solutions (3.19) and (3.21) are applicable for a case with $\rho_{0},P_{0}$ general constants.

For most observed stars, the thermal kinetic energy $k_{\text{B}}T$ is much smaller than the rest energy of the constituent particle $\bar{m}c^{2}$ near surface. For instance, the temperature of the Sun at the core is estimated as of order $10^{7}{\rm K}\approx 1{\rm KeV}/k_{\text{B}}$, which is even smaller than the order of the electron mass. This means that if such a star would consist mainly of the ideal gas of non-relativistic particles, then it would collapse and not exist stably. Thus an ideal gas of non-relativistic particles will not be appropriate as main constituent material of stable stars.

To the end of the construction of a model of a star with multilayer structure, finally consider the mixture of two ideal gases consisting of baryonic particles and non-relativistic ones as an example. The ideal gas of non-relativistic particles here describes atmosphere consisting of neutral particles and ionized ones, while that of baryonic particles does main stellar constituent material. The total pressure and the total energy density are given by

$\displaystyle\rho=\rho_{\rm b}+\rho_{\text{g}},~{}~{}~{}P=P_{\rm b}+P_{\text{g}},$

(3.23)

where the subscript $b$ is used for the baryonic component. The equation of state for an ideal gas of baryons is as usual given by

$$P_{\rm b}=\hat{n}_{\rm b}k_{\text{B}}T,$$

(3.24)

while the baryonic particle number conservation implies

$$\hat{n}_{\rm b}^{\prime}=\hat{n}_{\rm b}\frac{\rho^{\prime}}{\rho+P},$$

(3.25)

as explained in section 3.1. Substituting (3.24) into (3.25) so as to remove $\hat{n}_{\rm b}$ and employing (2.9), one finds $P_{\rm b}^{\prime}(P_{\text{g}}+\rho)=P_{\rm b}(P_{\text{g}}+\rho)^{\prime}$. This can be solved as $P_{\rm b}=w_{\rm b}(P_{\text{g}}+\rho)$, where $w_{\rm b}$ is an integration constant. $\rho_{\rm b}$ can be determined in terms of $T$ by solving the temperature gradient equation (2.9).

For a realistic situation to a typical star, the local temperature away from the core is much smaller than the rest energy of the non-relativistic particle per the Boltzmann factor, $T\ll T_{\text{g}}$, where $T_{\text{g}}:=\frac{\bar{m}c^{2}}{k_{\text{B}}}$. In this situation, $\rho_{\text{g}}\approx\bar{m}c^{2}\hat{n}_{\text{g}}$, where $\hat{n}_{\text{g}}$ is given by (3.12), and $P_{\text{g}}\ll\rho_{\text{g}}$, so that $P=(1+w_{\rm b})P_{\text{g}}+w_{\rm b}\rho\approx w_{\rm b}\rho$. This implies that in the realistic situation baryonic particles are dominant and the system can be approximately regarded as a simple fluid near and below the surface, so that the parameter $w_{\rm b}$ is related to the heat capacity ratio of the stellar material $\gamma_{\rm b}$ as $w_{\rm b}=\gamma_{\rm b}-1$, as investigated in section 3.1.⁵⁵5 This consequence is indeed preferable as a property of luminous stars with taking into account that there exists considerably large convection zone underneath the solar surface [34, 35]. The total energy density was determined to satisfy (2.9) as $\rho\approx a_{\rm b}T^{\frac{\gamma_{\rm b}}{\gamma_{\rm b}-1}}$, where $a_{\rm b}$ is an integration constant. To fix it, introduce the boundary of a layer of the stellar material whose radius is denoted by $R_{*}$ and impose the energy density of main stellar constituent material to vanish at the boundary of the layer.

$$\rho_{\rm b}(R_{*})=0.$$

(3.26)

This fixes the integration constant as $a_{\rm b}\approx\rho(R_{*})/T(R_{*})^{\frac{\gamma_{\rm b}}{\gamma_{\rm b}-1}}$, where $\rho(R_{*})=\rho_{\text{g}}(R_{*})\approx\bar{m}c^{2}\hat{n}_{\text{g}}(R_{*})$ and $T(R_{*})$ is the temperature at the surface of the lower layer related to the one at the surface of the upper layer as (3.22). On the other hand, the local temperature inside the star is approximated as the one to obey the power law as (3.9), $T\approx T(R_{*})(\frac{R_{*}}{r})^{2(1-\frac{1}{\gamma_{\rm b}})}$, with $T(R_{*})=\left(\frac{c^{2}(\gamma_{\rm b}-1)}{2\pi a_{\rm b}G\left(\gamma_{\rm b}^{2}+4\gamma_{\rm b}-4\right)}\right)^{1-\frac{1}{\gamma_{\rm b}}}\frac{1}{R_{*}^{2(1-\frac{1}{\gamma_{\rm b}})}}.$ Substituting the expression of $a_{\rm b}$ into this, the radius of the baryonic layer is determined as

$$R_{*}\approx\sqrt{\frac{c^{2}(\gamma_{\rm b}-1)}{2\pi\rho(R_{*})G\left(\gamma_{\rm b}^{2}+4\gamma_{\rm b}-4\right)}}.$$

(3.27)

The physical implication of this result is twofold. One is that a layer with bigger density goes deeper inside the star, which is physically preferable and indeed observed as the onion structure of a star. The other is that from this expression the deviation of the heat capacity ratio from one can be computed as

$$\gamma_{\rm b}\approx 1+2\pi\frac{\rho(R_{*})}{c^{2}}GR_{*}^{2},$$

(3.28)

where $\gamma_{\rm b}\sim 1$ was used.⁶⁶6 For the case of the Sun, it is evaluated as $\gamma_{\rm b}-1\approx 1\times 10^{-6}$, where the following data were used: $\rho(R_{*})/c^{2}\to m_{\rm u}/(\frac{4}{3}\pi a_{\text{B}}^{3})\approx 2.7\times 10^{3}{\rm kg/m^{3}},R_{*}\to\frac{R_{\star}}{2}=\frac{R_{\odot}}{2}=3.5\times 10^{8}\rm m$, with $a_{\text{B}}$ the Bohr radius, $a_{\text{B}}=5.3\times 10^{-11}\rm m$. This rough estimation implies that the deviation of the heat capacity ratio from one is very small. The total energy density and the baryonic component thereof are determined as

$\displaystyle\rho\approx\rho(R_{*})\left(\frac{R_{*}}{r}\right)^{2},~{}~{}\rho_{\rm b}\approx\rho\times\left(1-(\frac{R_{*}}{r})^{1-\frac{3}{\gamma_{\rm b}}}e^{T_{\text{g}}(\frac{1}{T(R_{*})}-\frac{1}{T})}\right).$

(3.29)

Then the baryonic component of pressure is given by $P_{\rm b}\approx w_{\rm b}\rho$, and the baryonic number density is $\hat{n}_{\rm b}=P_{\rm b}/k_{\text{B}}T$.

These results of thermodynamic observables are valid in a parametric region with $r_{\rm H}\ll r\leq R_{*},T\ll T_{\text{g}}$. Indeed, these expressions of thermodynamic variables develop singularity at the center of the system. This implies the existence of a core layer around the center in which thermodynamic behavior changes from the power law. Such a core deep inside a star is expected to consist of degenerate matter [21], so it may be very roughly approximated as a hydrostatic equilibrium state with uniform energy density. The pressure and temperature can be analytically computed as follows [41, 23]

$\displaystyle P=$

$\displaystyle\frac{\sqrt{r_{\rm c}^{2}-r^{2}}-\sqrt{r_{\rm c}^{2}-R_{\rm c}^{2}}}{3\sqrt{r_{\rm c}^{2}-R_{\rm c}^{2}}-\sqrt{r_{\rm c}^{2}-r^{2}}}\rho_{\rm c},~{}~{}T=\frac{2\sqrt{r_{\rm c}^{2}-R_{\rm c}^{2}}}{3\sqrt{r_{\rm c}^{2}-R_{\rm c}^{2}}-\sqrt{r_{\rm c}^{2}-r^{2}}}T(R_{\rm c}),$

(3.30)

where $\rho_{\rm c}$ is the constant energy density, $R_{\rm c}$ is the radius of the surface of the core, $r_{\rm c}=\sqrt{\frac{3}{8\pi G\rho_{\rm c}}}$.