#  Relativistic Hydrostatic Structure Equations
and Analytic Multilayer Stellar Model

Shuichi Yokoyama,
∗ Department of Physical Sciences, College of Science and Engineering,
Ritsumeikan University, Shiga 525-8577, Japan syr18046[at]fc.ritsumei.ac.jp

The relativistic extension of the classic stellar structure equations is
investigated. It is pointed out that the Tolman-Oppenheimer-Volkov (TOV)
equation with the gradient equation for local gravitational mass can be made
complete as a closed set of differential equations by adding that for the
Tolman temperature with one equation of state, and the set is proposed as the
relativistic hydrostatic structure equations. The exact forms of the
relativistic Poisson equation and the steady-state heat conduction equation in
the curved spacetime are derived. The application to an ideal gas of particles
with the conserved particle number current leads to a strong prediction that
the heat capacity ratio almost becomes one in any Newtonian convection zone
such as the solar surface. The steady-state heat conduction equation is solved
exactly in the system and thermodynamic observables exhibit the power law
behavior, which implies the possibility for the system to be a new model of
stellar corona and a flaw in the earlier one obtained by using the non-
relativistic stellar structure equations. The mixture with another ideal gas
yields multilayer structure to a stellar model, in which classic stellar
structure equations are reproduced and analytic multilayer structure of
luminous stars including the Sun is revealed in suitable approximation.

##  1 Introduction

As also read from descriptions in myths and Bible, the brightness of stars has
been recognized as something indispensable and reverent from ancient time. The
Sun has given energy to living beings during the day, shining stars have
attracted people like jewels in the night. How they maintain energy and beauty
has been a mystery before physicists unveil it. To reveal the mystery,
Helmholtz and Kelvin proposed a hypothesis that the brightness of the Sun
originated in the gravitational potential energy due to the contraction of
solar fluid [1, 2]. (See [3] for related history.) This proposal was rejected
because only the gravitational contraction cannot produce sufficient energy to
maintain the solar visibility longer than the time-scale of terrestrial
minerals. Taking into account the progress of nuclear physics, Eddington
pointed out that what is necessary for the Sun to gain sufficient luminous
energy is certain subatomic reaction process [4], whose insight was
subsequently developed further by Gamov and Bethe [5, 6]. For a spherically
symmetric star, this Eddington’s insight is described by the following
equation

$\displaystyle{dL_{r}\over dr}=$ $\displaystyle 4\pi r^{2}\varrho\epsilon,$
(1.1)

where $L_{r}$ is the energy flowing outwards across a sphere of radius $r$
called luminosity, $\varrho$ the mass density, and $\epsilon$ the net energy
production rate per unit mass liberated from all reaction processes including
subatomic ones. If such produced energy is transported by radiation, then the
temperature gradient inside a star is induced as [4]

$\displaystyle{dT\over dr}=$ $\displaystyle-\frac{3\kappa\varrho}{16\sigma
T^{3}}\frac{L_{r}}{4\pi r^{2}},$ (1.2)

where $\sigma=\pi^{2}ck_{\text{B}}^{4}/(60(c\hbar)^{3})$ is the Stefan-
Boltzmann constant and $\kappa$ is the opacity of stellar fluid per unit mass,
while if energy is transported by convection, the temperature gradient is
replaced by

$\displaystyle{dT\over dr}=$ $\displaystyle-(1-\frac{1}{\gamma})\frac{G_{\rm
N}M_{r}\varrho T}{r^{2}P},$ (1.3)

where $\gamma$ is the specific heat ratio of stellar fluid, $P$ the local
pressure, $G_{\rm N}$ the Newton constant, $M_{r}$ defined by

$\displaystyle{dM_{r}\over dr}=$ $\displaystyle 4\pi r^{2}\varrho.$ (1.4)

On the other hand, the visibility of a luminous star is based on its stable
existence, which is achieved as a result of the balance between the
gravitational force and the internal pressure of stellar fluid. In a region
where the gravitational force can be approximated by the Newtonian one, the
hydrostatic equilibrium can be described by

$\displaystyle{dP\over dr}=$ $\displaystyle-\frac{G_{\rm
N}M_{r}}{r^{2}}\varrho.$ (1.5)

The above set of differential equations is basic to investigate stellar
structure and known as the stellar structure equations. Solving them combining
equations of state for $P,\epsilon,\kappa$ in terms of $\varrho,T$ with
suitable boundary conditions, one can determine four macroscopic variables
$\varrho,T,M_{r},L_{r}$ as functions of $r$. (See standard textbooks [7, 8, 9,
10, 11, 12] and references therein.) These differential equations also yield
information on interdependent relations between global observables of luminous
stars such as the mass-luminosity relation [13] and the luminosity-temperature
one depicted by the Hertzsprung-Russell diagram [5]. (See [14, 15] for their
recent data.)

The above traditional stellar structure equations have certainly played
important roles to extract information on interior of luminous stars and
interdependent relations of their global observables. However, they are
apparently built on the basis of non-relativistic physics, which restricts
them to be applicable only in the Newtonian regime such that density and
pressure are sufficiently small like in the neighborhood of a stellar surface
and the inside of a light star. (See [16] for instance.) Indeed, the
relativistic extension of the equation for hydrostatic equilibrium (1.5) has
already been investigated and established as the Tolman-Oppenheimer-Volkov
(TOV) equation [17], and the deviation between non-relativistic results and
relativistic ones increases as the mass density does [18, 19]. (See also
[20].) Since there is an evolutionary process from a main sequence star to a
highly dense one consisting of degenerate matter [21], it would be desirable
to extend the stellar structure equations fully to their relativistic ones so
as to be applicable in the non-Newtonian regime such as the interior of a
degenerate star and in order to keep track of a stellar evolutionary process.

With taking this into account, it is natural to ask how the traditional
structure equations should be extended to the general relativistic ones, and
whether any significant consequence can be drawn from this extension
particularly to physics of luminous stars less dense than compact stars. (See
textbooks, for instance [22], on the study of a compact star using the TOV
equation.) The purpose of this paper is to address these questions by
employing some latest results of the author on relativistic local
thermodynamics in relativistic hydrostatic equilibrium with spherical symmetry
[23]. In the work, entropy current and entropy density were constructed by the
method proposed in [24], and the constructed entropy density was shown to
satisfy the local Euler’s relation and the first law of thermodynamics
concurrently and non-perturbatively with respect to the Newton constant, in
which the local temperature is exactly coincident with the one given by Tolman
[25], which is called the Tolman temperature in this paper. On top of this,
the established local thermodynamics was applied to a hydrostatic equilibrium
system with uniform energy density and the relativistic stellar structure was
completely determined with the analytic expressions of all local thermodynamic
observables. A lessen derived from this simple application is that the TOV
equation and two gradient equations for local gravitational mass and the
Tolman temperature form a closed set of differential equations with one
equation of state given.

Based on these results, the author puts forward that the TOV equation and two
gradient equations for the gravitational mass and the Tolman temperature are
the relativistic structure equations for spherically symmetric hydrostatic
equilibrium upgraded from the traditional stellar structure equations. A
nontrivial point in this proposal is that the proposed relativistic structure
equations do not contain variables on rate such as luminosity and energy
production rate, so that they cannot be obtained simply by extending each
traditional non-relativistic stellar structure equation to the relativistic
one. Then the main issue of this paper is twofold: whether the proposed
structure equations can reproduce the classic stellar structure equations
leading to important stellar properties, and whether any new significant
consequence can be drawn from the proposed relativistic ones properly.

The goal of the paper is to answer these questions positively and provide
evidence for validity and efficiency of the proposal. This is not trivial at
all taking into account the fact that the proposed relativistic hydrostatic
structure equations have less number of equations and that of variables as
well. To the end, firstly, the setup of a relativistic hydrostatic equilibrium
system with spherical symmetry is fixed and the necessary and sufficient set
of the proposed structure equations is explicitly presented in section 2. It
is pointed out that the proposed temperature gradient equation is consistent
with a thermodynamic relation known in the ordinary thermodynamics and thus it
is expected to hold for any local thermodynamic equilibrium system. Then the
relativistic Poisson equation is derived non-perturbatively in the Newton
constant, and it is converted into the differential equation for the Tolman
temperature using its relation to the gravitational potential. This is the
steady-state heat conduction equation exactly holding in this curved spacetime
and plays a key role to determining the hydrostatic structure.

In section 3, the proposed hydrostatic structure equations are applied to the
construction of a model with multilayer structure of luminous stars including
the Sun. Through this application, advantages of the relativistic extension of
the structure equations become transparent. One of the advantages is seen in
the temperature gradient equation. In the conventional stellar structure
equations, the temperature gradient equation needs to be chosen as either
(1.2) or (1.3) suitably in accordance with the way of energy transportation,
or is newly computed by using the so-called mixing length theory [26, 27]. In
the proposed ones, the temperature gradient equation is unchanged regardless
of transport phenomena and the energy transport of fluid is described by its
equation of state. In particular, it is shown that a simple fluid whose
pressure is proportional to its energy density satisfies the conventional
temperature gradient equation in a convective zone (1.3). Such a simple fluid
can be realized by an ideal gas of particles with conserved particle number
current, which is called baryonic particles in this paper. This model predicts
that the heat capacity ratio is almost one in any Newtonian convection zone
such as the solar surface. A characteristic feature of this new hydrostatic
equilibrium model is that the steady-state heat conduction equation can be
solved exactly and thermodynamic observables are determined non-perturbatively
in the Newton constant. They exhibit the power law behavior with the power law
index related to the heat capacity ratio and imply that this model can be also
used to describe an ionized state of fluid and a plasma one in stellar corona.
An interesting feature in this model is that pressure is totally well-behaved
and vanishing at the asymptotically far region. This result conflicts with the
earlier one obtained by using the non-relativistic Newtonian gravity [28], and
it concludes that the Newtonian approximation is not valid any more in stellar
corona. After applying to an ideal gas of non-relativistic particles to
determine the local temperature by perturbation in section 3.2, an analytic
stellar model is investigated as the hydrostatic equilibrium of two types of
ideal gases of baryonic particles and non-relativistic ones in section 3.3. In
the model, the ideal gas of baryonic particles forms a layer of the main
stellar material and that of non-relativistic ones does a layer of atmosphere.
In addition, by considering a situation with the system coupling to radiation,
the luminosity can be related to the local temperature, which leads to the
relativistic extension of the conventional temperature gradient equation in a
radiation zone (1.2). A summary of the multilayer stellar model is given in
table 1.

Layer | Region | Temperature $T$ | Energy density $\rho$ | Matter
---|---|---|---|---
Empty | $R_{\star}\leq r<\infty$ | $\displaystyle{T_{\star}\sqrt{\frac{1-\frac{r_{\rm H}}{R_{\star}}}{1-\frac{r_{\rm H}}{r}}}}$ | $0$ | ($\gamma$)
Atmosphere | $R_{*}\leq r\leq R_{\star}$ | $\displaystyle{T_{\star}\left(1-\frac{G_{\rm N}\nu_{*}}{2c^{2}}+\cdots\right)}$ | $(\bar{m}c^{2}+\frac{3}{2}k_{\text{B}}T)\hat{n}_{\text{g}}$ | NR(+$\gamma$)
Interior | ${R_{\rm c}\leq r\leq R_{*}}$ | $\displaystyle{T(R_{*})\left(\frac{R_{*}}{r}\right)^{2(1-\frac{1}{\gamma_{\rm b}})}}$ | $\displaystyle{\rho(R_{*})\left(\frac{R_{*}}{r}\right)^{2}}$ | NR+B(+$\gamma$)
Core | $0<r\leq R_{\rm c}$ | $T(R_{\rm c})\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}}}$ | $\rho_{\rm c}$ | D
Table 1: A summary of a simple analytic multi-layer stellar model is shown. In
each layer, the local temperature and the energy density determined in the
main text are shown. In the column of ’Matter’, the symbols ’NR’, ’B’,
’$\gamma$’, and ’D’ mean ’Non-Relativistic particles’, ’Baryonic particles’,
’photons’, and ’Degenerate matter’, respectively. The reason why photon is
enclosed by bracket is that it is not involved in the stellar structure
equations. $r_{\rm H}=\frac{2G_{\rm N}M_{\star}}{c^{2}},r_{\rm
c}=\sqrt{\frac{3c^{2}}{8\pi G_{\rm N}\rho_{\rm c}}}$, $M_{\star}$ is the
stellar gravitational mass, $\gamma_{\rm b}$ is the heat capacity ratio of the
ideal gas of baryonic particles, while $\nu_{*}$ is given by (3.19),
$\hat{n}_{\text{g}}$ (3.12), $T_{\star}$ (3.41), $R_{\star}$ (3.42),
$T(R_{*})$ (3.22), $R_{*}$ (3.27), in the main text.

Employing obtained results, a concrete model of the Sun is constructed in
section 3.3.1.

Section 4 is devoted to summary and discussion including open problems and
future works.

##  2 Relativistic hydrostatic structure equations

In order to generalize the hydrostatic structure equations so as to include
relativistic effect, a hydrostatic equilibrium system needs to be formulated
by general theory of relativity. Such a system with spherical symmetry
formulated by general relativity was studied in [25, 17, 29], (see also
appendix in [30]) in which the line element generally takes the form such that

$\displaystyle g_{\mu\nu}\mathrm{d}x^{\mu}\mathrm{d}x^{\nu}=$
$\displaystyle-e^{\nu}\mathrm{(}cdt)^{2}+e^{\lambda}\mathrm{d}r^{2}+r^{2}(\mathrm{d}\theta^{2}+(\cos\theta)^{2}\mathrm{d}\phi^{2}),$
(2.1)

where $\nu,\lambda$ are functions of the radial coordinate $r$, and the
energy-stress tensor is of the form of a perfect fluid

$\displaystyle
T_{\mu\nu}=(\rho+P)\frac{u_{\mu}}{c}\frac{u_{\nu}}{c}+Pg_{\mu\nu},$ (2.2)

where $\rho$ is local energy density, $P$ is local pressure, $u^{\mu}$ is the
fluid four velocity. Note that $c^{2}\nu/2$ corresponds to the gravitational
potential probed by a point particle with unit mass in the weak gravity regime
[22], so $\nu$ itself will be also called the gravitational potential for
convenience. The balance between the gravitational attractive force and the
repulsive one of fluid is described by the Einstein equation, which reduces in
the comoving frame to

$\displaystyle P^{\prime}=$ $\displaystyle-\frac{(P+\rho)}{2}\nu^{\prime},$
(2.3) $\displaystyle\lambda^{\prime}=$ $\displaystyle re^{\lambda}8\pi
G\frac{\rho}{c^{2}}-\frac{e^{\lambda}-1}{r},$ (2.4)
$\displaystyle\nu^{\prime}=$ $\displaystyle re^{\lambda}8\pi
G\frac{P}{c^{2}}+\frac{e^{\lambda}-1}{r},$ (2.5)

where $G:=G_{\rm N}/c^{2}$. Eliminating $\lambda$ and $\nu$ leads to the TOV
equation [17]

$\displaystyle P^{\prime}=$
$\displaystyle-\frac{(P+\rho)}{c^{2}}\frac{G(\check{E}_{r}+4\pi
r^{3}P)}{r^{2}(1-2\frac{G\check{E}_{r}}{rc^{2}})},$ (2.6)
$\displaystyle\check{E}_{r}^{\prime}=$ $\displaystyle 4\pi r^{2}\rho,$ (2.7)

where $\check{E}_{r}$ is defined by $\check{E}_{r}:={rc^{2}\over
2G}(1-e^{-\lambda})$. Note that $\check{E}_{r}$ has the dimension of energy.

Employing the method given in [24] one can construct entropy current and
entropy density as a conserved current and a conserved charge density for this
system [23]. (See also [31].) A key to the construction is to find a vector
field $\xi^{\mu}$ to satisfy a differential equation such that
$T^{\mu}\\!_{\nu}\nabla_{\mu}\xi^{\nu}=0$. This vector field enables one to
construct a conserved current as $J^{\mu}=\sqrt{-g}T^{\mu}\\!_{\nu}\xi^{\nu}$
and a conserved charge as $Q=\int d^{3}xJ^{t}$ in general. In order to find
conserved entropy current for the fluid, look for such a vector field
$\xi^{\mu}$ to be proportional to the fluid velocity $u^{\mu}$. The entropy
density constructed in this way was shown to satisfy the local Euler’s
relation and the first law of thermodynamics concurrently and non-
perturbatively with respect to the Newton constant in the comoving frame [23].

$Ts=u+Pv,~{}~{}Ts^{\prime}=u^{\prime}+Pv^{\prime},$ (2.8)

where $v=r^{2}\cos\theta(1-{2G\check{E}_{r}\over rc^{2}})^{-{1\over 2}}$,
$u=\rho v$ is the internal energy density, $T\propto e^{-\frac{\nu}{2}}$ is
the local temperature given by Tolman [25]. In particular, it was shown that
the entropy density does not depend on the radial coordinate, which is a
natural consequence for a steady state with vanishing heat flux in the energy
stress tensor. Using (2.3), one can derive the temperature gradient equation
as [23]

$T^{\prime}={P^{\prime}\over P+\rho}T.$ (2.9)

Plugging (2.6) into this leads to

$T^{\prime}=-\frac{G(\check{E}_{r}+4\pi
r^{3}P)}{c^{2}r^{2}(1-2\frac{G\check{E}_{r}}{rc^{2}})}T.$ (2.10)

This temperature gradient equation (2.10) and the TOV equation (2.6) with
(2.7) form a closed set of differential equations including one equation of
state, and this set of three differential equations is proposed as the
relativistic extension of classic stellar structure equations.

Comments are in order. Firstly, the temperature gradient equation (2.9)
involves only thermodynamic variables and their gradients, which suggests that
a corresponding relation exists in global thermodynamics on flat spacetime.
The answer of such a thermodynamic relation is

$\bar{T}\left({\partial\bar{P}\over\partial\bar{T}}\right)_{\bar{V},\bar{N}}=\left({\partial\bar{U}\over\partial\bar{V}}\right)_{\bar{T},\bar{N}}+\bar{P},$
(2.11)

since this can be rewritten as
$({\partial\bar{P}\over\partial\bar{T}})_{\bar{V},\bar{N}}=\frac{\bar{T}}{\bar{P}+\bar{\rho}}$,
where $\bar{\rho}=({\partial\bar{U}\over\partial\bar{V}})_{\bar{T},\bar{N}}$
and the thermodynamic variables with overline stand for global quantities used
here only.111 It may not be so trivial to obtain the local thermodynamic
relation (2.9) from (2.11), because the energy density $\rho$ appearing in the
energy stress tensor and the internal energy density $u$ in the thermodynamic
relations are different on curved spacetime as shown above, which could make
it difficult to presume (2.9) from (2.11).  Therefore the relation (2.9) is
expected to hold for any local thermodynamic equilibrium system. Secondly, as
is clear from the above derivation, there are several equivalent expressions
for the temperature gradient equation. (2.10) may be useful for numerical
analysis, while (2.9) or (2.3) may be more useful for non-numerical one. Which
expression is best to be used will depend on contexts and ways of analysis.
Thirdly, although the constructed entropy density does not play any role in
the structure equations, the existence thereof has important implications. One
is that this ensures the validity of the definitions of the macroscopic
quantities such as the temperature and the internal energy due to the fact
that they satisfy the local Euler’s relation and the first law of
thermodynamics. Another is to assure the system to be in a local equilibrium,
which justifies to solve the hydrostatic structure equations by assuming the
distribution of an ideal statistical ensumble for local thermodynamic
quantities. This point together with (2.9) will become important to include
the temperature effect dynamically. In earlier study, the TOV equation (2.6)
and (2.7) were used to determine the behavior of pressure and density by
solving them combined with an equation of state with respect to pressure and
density only [22], but they can be used to determine the local temperature and
include the thermal effect incorporating the temperature gradient equation
(2.9). Fourthly, the proposed hydrostatic structure equations do not contain
the rate variables such as the luminosity or the corresponding equation for
the luminosity gradient (1.1). They are closed without the variable of
luminosity, which implies that the variables on rate are not fundamental for
hydrostatic equilibrium. This result is preferable to the fact of the
existence of a non-radiating star in hydrostatic equilibrium, while this in
turns implies that the proposed hydrostatic structure equations are not
applicable to a situation where there is so strong radiant energy flux to
involve hydrostatic structure. The restriction of applicability is discussed
in section 3.3.2. Finally, on the non-relativistic reduction, it is easily
seen that the TOV equation with (2.7) reduce to the non-relativistic
hydrostatic equation (1.5) with (1.4) in the Newtonian regime where
$P\ll\rho,P\ll\frac{\check{E}_{r}}{4\pi r^{3}},\frac{2G_{\rm
N}\check{E}_{r}}{c^{2}}\ll r$ so that $\rho\to\varrho c^{2},\check{E}_{r}\to
M_{r}c^{2}$. On the other hand, it needs some preparation of setup to
reproduce (1.2) and (1.3), so it will be confirmed later.

###  2.1 General relativistic Poisson equation and steady-state heat
conduction equation

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$.222 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.

##  3 Applications

In this section, the relativistic hydrostatic structure equations proposed in
the previous section are applied to constructing a model of luminous stars.
This application also demonstrates how classic results on a luminous star are
reproduced or modified with taking into account the effect of general
relativity and how to determine local thermodynamic observables inside a
luminous star.

###  3.1 Ideal gas of baryonic particles

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:333 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.

####  3.1.1 Power law behavior of thermodynamic observables

The radial dependence of temperature is determined by solving the relativistic
steady-state heat conduction equation (2.24), which reduces to

$\displaystyle T\left(T^{\prime\prime}+\frac{2}{r}T^{\prime}\right)-3T^{\prime
2}$ $\displaystyle=$ $\displaystyle\frac{4\pi
aG}{c^{2}}T^{\frac{\gamma}{\gamma-1}}\left(-6r^{2}T^{\prime
2}+2r^{2}TT^{\prime\prime}+2T^{2}+\gamma\left(4r^{2}T^{\prime
2}-rT\left(2rT^{\prime\prime}-3T^{\prime}\right)-3T^{2}\right)\right).$ (3.7)

This can be solved perturbatively in the Newton constant as illustrated in
section 2.1, while this expression implies the existence of an exact solution
obeying the power law, $T=Ar^{n}$, where $A$ is a nonzero constant and $n$ is
a negative number so that the temperature increases as it goes deeper into the
center. These parameters are determined by substituting back into the
differential equation, which simplifies to

$A^{2}(2n-1)r^{2n-2}\left(4\pi
aGA^{\frac{\gamma}{\gamma-1}}(n(\gamma-2)+3\gamma-2)r^{\frac{n}{\gamma-1}+n+2}+c^{2}n\right)=0.$
(3.8)

There exists a nontrivial solution if and only if $\frac{n}{\gamma-1}+n+2=4\pi
aGA^{\frac{\gamma}{\gamma-1}}(n(\gamma-2)+3\gamma-2)r^{\frac{n}{\gamma-1}+n+2}+c^{2}n=0$.
Thus the parameters for a nontrivial solution are determined as
$n=-2(\gamma-1)/\gamma,A=\left(\frac{c^{2}(\gamma-1)}{2\pi
aG\left(\gamma^{2}+4\gamma-4\right)}\right)^{\frac{\gamma-1}{\gamma}}$. As a
result, the relativistic steady-state heat conduction equation has been solved
exactly as

$T=\left(\frac{c^{2}(\gamma-1)}{2\pi
aG\left(\gamma^{2}+4\gamma-4\right)}\right)^{1-\frac{1}{\gamma}}\frac{1}{r^{2(1-\frac{1}{\gamma})}}.$
(3.9)

That is, the local temperature obeys the power law. Thus the density and the
pressure do so as well in such a way that $P\propto\rho\propto 1/r^{2}$.
Therefore it concludes that a hydrostatic equilibrium system of a simple fluid
exhibits the power law behavior of macroscopic observables.

It is important to stress that this solution cannot be obtained by
perturbation in the Newton constant, so that the power law behavior cannot be
seen by solving the system perturbatively. It is also interesting that the
parameters in this solution are given only by data of internal stellar fluid
and not by external data such as the gravitational mass. From this result, the
gravitational potential for this model can be also determined as
$\nu=4(1-\frac{1}{\gamma})\log(r/R_{\star}).$

####  3.1.2 New model for stellar corona

The above results of the power law behavior of thermodynamic observables imply
that this model is applicable to stellar corona, in which temperature is so
much higher than the other part of stellar superficial region that gas is
almost or completely ionized. This suggests that macroscopic behavior and
transport phenomena in stellar corona can be investigated by combining plasma
physics with stellar structure equations. For instance, temperature is
predicted to fall off in a power law in solar corona. It was pointed out by
Parker that combining the result of plasma physics for fully ionized gas with
the conventional non-relativistic hydrostatic structure equations leads to
non-vanishing pressure at infinity [28]. He speculated on this as a signal of
the impossibility for gas in the solar corona to reach hydrostatic equilibrium
and presumed that the non-vanishing pressure at asymptotic region gives rise
to gas streaming outward from the Sun called the solar wind, whose existence
was suggested by Biermann earlier [36, 37]. The solar wind was indeed observed
by a spacecraft launched some years later [38, 39].

This is a successful interplay between theory and experiment. However, it is
also pointed out that there still remains unsolved issues for the corona and
the solar wind by using the traditional stellar structure equations and
results of plasma physics [40]. In the application of the above relativistic
hydrostatic model of a simple fluid to stellar corona, however, there is no
problem for the system to stay in hydrostatic equilibrium. Indeed, not only
the local temperature but also the local pressure behave in the power law as
shown above, which leads to pressure vanishing in the asymptotic region.

This result conflicts with the earlier one obtained by using the non-
relativistic structure equations explained above. Then one might wonder which
result is indeed correct. The answer of the author is that the earlier model
contains a flaw in the argument. To explain this, one first needs to accept
the fact that the relativistic hydrostatic structure equations presented in
this paper, which are rigorously derived from the Einstein gravity, are more
fundamental than the non-relativistic ones derived from the Newtonian gravity.
Indeed, it has been confirmed that the relativistic hydrostatic structure
equations reduce to the non-relativistic ones by neglecting terms due to the
effect of general relativity. Such an approximation is valid only in the
Newtonian regime, where the pressure is negligible compared to other terms.
This implies that the Newtonian approximation cannot be justified for a
solution such that the pressure is non-vanishing at a far distance. More
concretely, one of the conditions to justify the non-relativistic reduction,
$4\pi r^{3}P\ll\check{E}_{r}$, is clearly violated for such a situation, and
the contribution from pressure in (2.6) is not negligible at all. Parker
derived the analytic expression of the number density for an ideal gas with
temperature behaving as the power law by solving the non-relativistic
hydrostatic equation (1.5), and found that it becomes divergent at infinity
[28], while, in the current system of an ideal gas with temperature obeying
the power law, the number density has been obtained as

$\hat{n}=\frac{P}{k_{\text{B}}T}=\frac{a(\gamma-1)}{k_{\text{B}}}\left(\frac{c^{2}(\gamma-1)}{2\pi
aG\left(\gamma^{2}+4\gamma-4\right)}\right)^{\frac{1}{\gamma}}\frac{1}{r^{\frac{2}{\gamma}}},$
(3.10)

which also obeys the power law and is vanishing in the asymptotic region. The
latter fully contains general relativistic effect.

As a result it concludes that the Newtonian approximation is not valid for a
hydrostatic equilibrium system in a plasma state such as stellar corona. Note
that this conclusion does not deny the existence of the solar wind, whose
existence must be accounted for by another mechanism. It is an interesting
future work to study whether a newly proposed model compatible with
relativistic structure equations can demonstrate any property on stellar
corona and make prediction on it.

####  3.1.3 Comment on application to radiative fluid

One might wonder that this result is applicable also to an ideal gas of photon
or the radiation dominant region, which is described by the equation of state
(3.3) with $w=\frac{1}{3}$, so $P_{\rm rad}=\frac{1}{3}\rho_{\rm rad}\propto
T^{4}$, $\rho_{\rm rad}=aT^{4}$. The integration constant $a$ is fixed by
considering a surface of a layer444 This layer may be called photosphere,
which is conventionally specified by a single surface at a certain optical
depth $\bar{\tau}=\frac{2}{3}$ defined by
$\bar{\tau}=\int_{R}^{\infty}\kappa\varrho dr$, where $R$ is a radius of a
luminous star [11].  at which the Stefan-Boltzmann’s law holds, so that
$a=4\sigma/c$. This leads to $\rho_{\rm rad}=\frac{4\sigma}{c}T^{4},$ which is
consistent with the Bose distribution in equilibrium. Then the behavior of the
local temperature would be given by (3.9) with the parameters specified above
for the radiation dominant region. However, there is an important caveat for
the application of the relativistic hydrostatic structure equations. That is,
they were derived in the comoving observer with the perfect fluid, and such an
observer does not exist for null fluid. Thus the application of the above
results to radiation is not supported strictly speaking, though it may be
useful merely as a rough estimation. Note that this caveat is also the case
for the non-relativistic stellar structure equations. The proper analysis to
include the effect of radiation to hydrostatic structure is left as future
work.

###  3.2 Ideal gas of non-relativistic particles

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.

###  3.3 Analytic multilayer stellar model

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.555 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.666 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}}}$.

Beyond the layer of main stellar material, $r>R_{*}$, the existing stellar
matter is only the ideal gas of non-relativistic particles. Therefore,
thermodynamic observables in this region are already computed in the previous
section 3.2.

Once energy density is determined, $M_{*}$ and $M_{\star}$ are computed from
(2.7) as

$M_{*}=\frac{1}{c^{2}}\int_{0}^{R_{*}}dr4\pi
r^{2}\rho,~{}~{}M_{\star}=M_{*}+\frac{1}{c^{2}}\int_{R_{*}}^{R_{\star}}dr4\pi
r^{2}\rho.$ (3.31)

The difference of these masses can be estimated as (3.20). If the core layer
is negligibly small and the energy density can be always approximated by
(3.29) inside the star, then $M_{*}$ is evaluated as

$M_{*}\approx\frac{\rho(R_{*})}{c^{2}}4\pi R_{*}^{3}$ (3.32)

This has a simple physical meaning that the total mass of a star can be
estimated just by the product of mass density at the surface of the baryonic
layer and its volume up to a numerical factor $3$. To confirm validity of this
expression, below estimate the solar mass by using it. To this end, assume
that the contribution of the atmosphere to the mass is negligible,
$M_{\star}\approx M_{*}$, and that the radius of the upper layer is set to be
identical to the solar one related to that of the lower layer as
$R_{\star}\equiv R_{\odot}=\eta R_{*}$, where $\eta$ is a numerical factor
greater than 1. On the other hand, the mass density at the surface of the
baryonic layer is estimated as
$\rho(R_{*})/c^{2}=\bar{m}\hat{n}_{\text{g}}(R_{*})=\mu_{\rm e}m_{\rm
u}/(\frac{4}{3}\pi a_{\text{B}})^{3}$, where the number density is
$\hat{n}_{\text{g}}(R_{*})\approx 1/(\frac{4}{3}\pi a_{\text{B}})^{3}$ with
$a_{\text{B}}$ the Bohr radius and $\mu_{\rm e}\approx 1.2$ since the weight
fraction of hydrogen is $X=0.6$ for the case of the Sun [42].777This
approximation means that the Sun is treated almost as a dense pack of
hydrogens at surface. From these data, the solar gravitational mass is
evaluated as

$M_{\star}\approx\frac{\mu_{\rm e}m_{\rm u}}{\frac{4}{3}\pi
a_{\text{B}}^{3}}\times
4\pi(\frac{R_{\odot}}{\eta})^{3}\approx\frac{1}{\eta^{3}}\times 10^{34}\rm g.$
(3.33)

This estimation yields the same order of the solar mass, $M_{\odot}=2.0\times
10^{33}\rm g$.

In order to construct a model of a luminous star, consider energy transported
by radiation measured by the luminosity $L_{r}$. In this paper, the luminosity
is literally defined as the energy of radiation only flowing outwards across
the sphere of radius $r$ per unit time here.888 In some contexts, the
luminosity is defined as the net energy flow transported not only by radiation
but also by convection in a more involved context such as the mixing length
theory. The reason to adopt the definition in this paper is because the net
energy flow will be evaluated by using entropy as $TdS$.  This is evaluated in
the background metric (2.1) as

$L_{r}=4\pi r^{2}e^{\lambda/2}F,$ (3.34)

where $F$ is the radiative energy flux summed all over the frequencies $F=\int
d\nu F_{\nu}$ and $F_{\nu}$ is the energy flux component of radiation with
frequency $\nu$. Now assume that the source of radiation and the absorption
body are absent. Then the energy flux component $F_{\nu}$ is related to its
pressure $P_{\nu}$ by ${dP_{\nu}\over d\tau_{\nu}}=\frac{F_{\nu}}{c},$ where
$\tau_{\nu}$ is the frequency-sensitive optical depth defined through the
absorption coefficient $\alpha_{\nu}$ as $d\tau_{\nu}=-\alpha_{\nu}dr.$
Therefore ${dP_{\nu}\over dr}=-\alpha_{\nu}\frac{F_{\nu}}{c}$. Summing all
over the modes labeled by the frequency $\nu$ leads to

$P_{\rm rad}^{\prime}=-\alpha_{\rm R}\frac{F}{c},$ (3.35)

where $\alpha_{\rm R}$ is an averaged absorption coefficient called the
Rosseland mean. The radiative component of pressure is given by $P_{\rm
rad}=\frac{4\sigma}{3c}T^{4}$, so $T^{\prime}=\frac{3c}{16\sigma T^{3}}P_{\rm
rad}^{\prime}$. Substituting (3.35) into this yields

$T^{\prime}=-\frac{3\alpha_{\rm R}}{16\sigma T^{3}}F.$ (3.36)

Then using (3.34) one finds the relation between the temperature gradient and
the local luminosity as

$\displaystyle T^{\prime}=-\frac{3\alpha_{\rm R}}{16\sigma
T^{3}}\frac{L_{r}}{4\pi r^{2}(1-{2G\check{E}_{r}\over rc^{2}})^{-{1\over
2}}}.$ (3.37)

This is the relativistic extension of the temperature gradient equation (1.2)
for the radiation dominant zone. As a byproduct, one can derive a relativistic
result on the luminosity by plugging (2.10) into (3.37) as

$\displaystyle L_{r}=\frac{64\pi\sigma}{3c^{2}\alpha_{\rm
R}}\frac{G(\check{E}_{r}+4\pi
r^{3}P)}{(1-2\frac{G\check{E}_{r}}{rc^{2}})^{3/2}}T^{4}.$ (3.38)

Indeed, near the stellar surface, (3.37) reduces to (1.2) with taking into
account the relation between the Rosseland mean and the opacity as
$\alpha_{\rm R}\approx\varrho\kappa$ with $\varrho$ the mass density of
stellar matter, as asserted before, while (3.38) reduces to

$\displaystyle L_{r}\approx\frac{64\pi\sigma}{3c^{2}\alpha_{\rm
R}}G\check{E}_{r}T^{4}.$ (3.39)

Note that these expressions for the luminosity are valid in the region with
matter present. In the matter empty region, the Rosseland mean and opacity
vanish, and from (3.36), the temperature is constant. This is consistent with
the result of the local temperature in the matter empty region (2.23)
effectively.

To extract physical implication, evaluate (3.39) in the atmospheric layer,
$R_{*}\leq r\leq R_{\star}$. Near the upper surface, $r\lesssim R_{\star}$,
the pressure is negligible as long as (2.22) is satisfied, and $\check{E}_{r}$
and $T$ are almost constant from (2.12) and (2.23), respectively, as seen in
section 2.1. Thus

$L_{r}\approx\frac{64\pi\sigma}{3\alpha_{\rm
R\star}}GM_{\star}T(R_{\star})^{4}.$ (3.40)

Here $\alpha_{R\star}$ is the Rosseland mean near the upper surface and
becomes a constant. This shows that the luminosity in this region is also
constant, which is consistent with the Stefan-Boltzmann law, $L_{r}\approx
L_{\star}$, where

$L_{\star}=4\pi R_{\star}^{2}\sigma T_{\star}^{4}.$ (3.41)

Here the radius $R_{\star}$, which was introduced so as to satisfy the
condition (2.22), is defined for the Stefan-Boltzmann law to hold at the
radius, and $T_{\star}$ is the effective surface temperature identified with
$T_{\star}=T(R_{\star})$. Then a novel relation is obtained as

$R_{\star}^{2}=\frac{16GM_{\star}}{3\alpha_{\rm R\star}}.$ (3.42)

This can be used to estimate the Rosseland mean at the measured point.999 The
expression (3.42) cannot be used in a region of stellar corona, since, as
argued in section 3.1.2, the Newtonian approximation is not valid any more in
a plasma state of hydrostatic equilibrium. In order to estimate the Rosseland
mean in stellar corona, one has to use an expression obtained without using
the Newtonian approximation such that $\alpha_{\rm
R\star}=\frac{16G(M_{\star}+4\pi
R_{\star}^{3}P(R_{\star})/c^{2})}{3R_{\star}^{2}}$.  It would be interesting
to compare the estimated result and detailed data of the opacity [43, 44, 11].

####  3.3.1 Analytic multilayer solar model

In order to build a concrete model of a star from the above multilayer stellar
model, internal data of the star is necessary. To illustrate how to build such
a concrete model, here construct a multilayer model of the Sun as an example.
It is known that the Sun has a convection zone up to around three quarters of
the total radius from the surface and the deeper rest is a radiation zone
except the core [34, 35]. Here, for simplicity, the region of solar corona is
neglected for simplicity, and the convective zone and the radiative one are
both realized by a baryonic layer consisting of an ideal gas of baryonic
particles with the heat capacity ratio varying in the layer. Since the total
solar mass is already estimated correctly as seen in (3.33), below investigate
behaviors of local observables.

The surface of the atmosphere starts from the radius $R_{\star}$, which is set
to the solar one, $R_{\star}:=R_{\odot}$. The layer of atmosphere continues
roughly up to the half of the solar radius, $R_{*}={1\over 2}R_{\star}$, from
which the baryonic layer starts. The baryonic layer is roughly divided into a
convective zone, which ends at the radius $\tilde{R}_{*}=\frac{3}{4}R_{*}$,
and a radiative one, which ends at $r=R_{\rm c}$ and connects to the core
smoothly. Here the radius of the core is chosen as $R_{\rm c}:=5r_{\rm H}$ by
hand for illustration. Accordingly the heat capacity ratio starts with $1$ at
the surface of the baryonic layer and increases as it goes deeper. It takes
the value $4/3$ at the transition point from the convective zone to the
radiative zone, and $5/3$ at the one from the baryonic layer to the core,
since the core is expected to be supported by repulsive force of degenerate
electron gas [21], whose heat capacity ratio is $5/3$ [45]. Note that this
realization of the radiation zone is different from the traditional one,
because the temperature gradient equation obeys (3.5) in an ideal gas of
baryonic particles, which is supposed to hold not in a radiative zone but in a
convective one. There are two reasons for this. One is that, as argued in
section 3.1.3, the presented relativistic hydrostatic structure equations
cannot be applied for null fluid. This does not allow one to solve the
structure equations combining with the equation (3.37). The other is that if
only the equation (3.37) is used as the stellar structure equations without
using (3.5), then one of the relativistic structure equations (2.9), which is
expected to always hold in local thermodynamic equilibrium, is not satisfied.
On these reasons, the radiative zone is realized by a simple fluid to satisfy
the equation (3.5) with $w=1/3$ in this model.

The behaviors of local observables inside the Sun are summarized in table 1
and figure 1.

(a) Solar multilayer structure.

(b) Behaviors of local observables.

Figure 1: (a) A schematic picture of solar multilayer structure and (b) the
behaviors of local thermodynamic observables are drawn. Three local
observables are plotted in the same graph even though they have different
units in order to describe their behaviors concisely.

Their characteristic values are easily calculated as follows. On the local
temperature, as an initial condition, $T_{\star}$ is set to the solar surface
temperature observed by an asymptotic observer in accord with the Stefan-
Boltzmann law, $T_{\star}:=T_{\odot}\approx 6\times 10^{3}{\rm K}$. The
temperature does not change so much in the layer of atmosphere, so that it
lasts at the surface of the baryonic layer, $T(R_{*})\approx T_{\odot}$. Then
it increases as it goes deeper. At the transition point from the convective
zone to the radiative one,
$T(\tilde{R}_{*})=T(R_{*})\left(\frac{R_{*}}{\tilde{R}_{*}}\right)^{1/2}\approx
7\times 10^{3}{\rm K}$, while at the point from the baryonic layer to the core
$T(R_{\rm c})=T(R_{*})\left(\frac{R_{*}}{R_{\rm c}}\right)^{5/4}\approx
2\times 10^{7}{\rm K}$. Note that in the current case, $r_{\rm c}\gg R_{\rm
c}$, so the temperature of the center of the system is almost equal to the one
at the surface of the core, $T(0)\approx T(R_{\rm c})$. The core temperature
estimated here is of the same order as the earlier one obtained by using the
classic stellar structure equation [34, 46].

On the density, an initial condition is given at the boundary between the
atmosphere and the interior, $r=R_{*}$, where the density is estimated as a
dense pack of hydrogens, $\rho(R_{*})/c^{2}\approx m_{\rm u}/(\frac{4}{3}\pi
a_{\text{B}}^{3})\approx 3\times 10^{3}{\rm kg/m^{3}}$. Then at the transition
point from the convective zone to the radiative one
$\rho(\tilde{R}_{*})/c^{2}=\rho(R_{*})/c^{2}\left(\frac{R_{*}}{\tilde{R}_{*}}\right)^{2}\approx
5\times 10^{3}{\rm kg/m^{3}}$, while at the core $\rho_{\rm
c}/c^{2}=\rho(R_{*})/c^{2}\left(\frac{R_{*}}{R_{\rm c}}\right)^{2}\approx
1\times 10^{12}{\rm kg/m^{3}}$. This density at the core is much bigger than
the earlier result around $1.5\times 10^{5}{\rm kg/m^{3}}$ [34, 46]. In fact,
the value is still bigger than the known average density of a typical white
dwarf around $3\times 10^{9}{\rm kg/m^{3}}$ but less than that of a neutron
star around $5\times 10^{14}{\rm kg/m^{3}}$ [41]. However, this result of high
density is in some sense reasonable because the core is assumed to be
supported by repulsive force of degenerate electrons so its density is
expected to be comparable to that of a white dwarf. Accordingly the pressure
in the model also becomes much higher than the earlier result. The behaviors
of local observables for this solar model are plotted in figure 1 (b) by
determining the heat capacity ratio in each zone by the linear interpolation
using the boundary values. They are qualitatively different from the earlier
ones [34, 46]. This qualitative difference presumably arises due to that of
the realization of the radiative zone.

It is important to stress that the values obtained here are very rough without
fine tuning. Therefore they can be possibly made more precise by tuning the
parameter with more precise data. Such a precise evaluation and further
discussion of physical interpretations of these differences are left to future
work.

####  3.3.2 Comment on the Eddington bound

In the above model of luminous stars, the radiation pressure is assumed to be
too small to involve hydrostatic structure [43]. This approximated treatment
is not used in the mixing length theory [26, 27], in which both the convective
energy flux and the radiant energy flux involve the traditional hydrostatic
structure equations. (See also [47].) This approximation will be valid as long
as the outward force caused by radiation is sufficiently small compared to the
gravitational one produced by matter inside the star:

$|{dP_{\rm rad}\over dr}|\ll|{dP\over dr}|.$ (3.43)

Using (3.34) and (3.38), one can compute the left-hand side as $|P_{\rm
rad}^{\prime}|=\frac{\alpha_{\rm R}L_{r}}{4\pi
r^{2}(1-2\frac{G\check{E}_{r}}{rc^{2}})^{-{1\over 2}}c}$, while the right-hand
is already computed as (2.6). Therefore the above inequality becomes

$L_{r}\ll(1+w_{\rm b})\rho\frac{4\pi G(\check{E}_{r}+4\pi r^{3}P)}{\alpha_{\rm
R}c(1-2\frac{G\check{E}_{r}}{rc^{2}})^{\frac{3}{2}}},$ (3.44)

where $P\approx w_{\rm b}\rho$ was used. Further estimation of the right-hand
side needs to specify the expression of the Rosseland mean. In a non-
relativistic hydrostatic system, the Rosseland mean takes the form of the
product of the mass density $\varrho$ and the opacity $\kappa$. In a general
relativistic hydrostatic system, the energy density $\rho$ is more fundamental
than the mass density $\varrho$, so the Rosseland mean is defined by using the
energy density as $\alpha_{\rm R}=\kappa\rho/c^{2}$. Then the above can be
rewritten as

$L_{r}\ll\gamma_{\rm b}\frac{4\pi cG(\check{E}_{r}+4\pi
r^{3}P)}{\kappa(1-2\frac{G\check{E}_{r}}{rc^{2}})^{\frac{3}{2}}}.$ (3.45)

The right-hand side corresponds to the relativistic extension of the Eddington
limit for this stellar model. Indeed, evaluating this around a stellar surface
with the Newtonian approximation valid, one finds

$L_{R_{\star}}\ll\frac{4\pi cG_{\rm N}M_{\star}}{\kappa}.$ (3.46)

##  4 Discussion

The relativistic extension of the classic stellar structure equations has been
investigated in the Einstein gravity. It has been pointed out that the TOV
equation and two gradient ones for local gravitational mass and the Tolman
temperature form a closed set of differential equations, and the set has been
put forward as the desired relativistic extension of the stellar structure
equations. The proposed relativistic hydrostatic structure equations have been
shown to be endorsed by laws of local thermodynamics.101010 This point should
not be overlooked for a precise model building of a celestial body. A certain
polytropic relation is often assumed to study a stellar structure in classic
papers such as [32] and standard textbooks on astrophysics, but it is not
clear or described whether it is indeed endorsed by laws of local
thermodynamics. If it is not, then such a model is a rough one which awaits
replacement by a thermodynamically correct one, even though it has contributed
enormously to the understanding of physics of an astronomical body.  From
them, the exact form of the relativistic Poisson equation has been derived,
and it has been converted into the steady-state heat conduction equation
holding non-perturbatively with respect to the Newton constant by using the
relation between the gravitational potential and the Tolman temperature.

A couple of applications have been presented. One is to a hydrostatic
equilibrium system consisting of a single ideal gas of particles with the
particle number current conserved called baryonic particles. In this system,
the proposed temperature gradient equation reduces precisely to the
conventional one known in the convection zone. It has been predicted that the
heat capacity ratio almost becomes one in the Newtonian convective regime such
as the neighborhood of the solar surface. The generalized steady-state heat
conduction equation has been solved exactly and thermodynamic observables
determined non-perturbatively in the Newton constant, which exhibit the power
law behavior such that $T\propto r^{-2(\gamma-1)\over\gamma},\rho\propto
p\propto r^{-2}$, with $\gamma$ the heat capacity ratio. This result implies
that this model is applicable to a plasma state and thus to stellar corona.
Another application is to a hydrostatic equilibrium system of an ideal gas of
non-relativistic particles, for which the local temperature has been
determined perturbatively in the Newton constant. Finally, by combining two
ideal gases of baryonic particles and non-relativistic ones, an analytic
multilayer structure of luminous stars has been revealed with a simple
example. In this model, the interior consists mainly of the ideal gas of
baryonic particles, while the atmosphere is constituted by non-relativistic
particles. This model also admits the layer of degenerate core and the stellar
corona by changing the constituents suitably. By coupling the system to
radiation, the traditional temperature gradient equation in a radiation zone
has been extended to its relativistic version. Inputing observational data
into the model, a concrete model of the Sun has been constructed and analytic
expressions of thermodynamic observables have been determined. Finally, the
Eddington bound has been obtained as a condition to justify the approximation
such that the contribution of radiation is too small to be involved in the
structure equations.

In the proposed hydrostatic structure equations, the relativistic extension
for the luminosity gradient (1.1) is not included. This is simply because
structure equations are closed without the variable of luminosity. This
implies that the variables on rate are not fundamental for relativistic
hydrostatic equilibrium. This conclusion is actually preferable with taking
into account the existence of a non-radiating star which can be described by a
relativistic hydrostatic equilibrium system. However, this does not mean that
the gradient equation for luminosity is not useful in a relativistic
hydrostatic equilibrium system. Its naive relativistic extension will be given
by

$\displaystyle{dL_{r}\over dr}=$ $\displaystyle\frac{4\pi
r^{2}}{(1-2\frac{G\check{E}_{r}}{rc^{2}})^{\frac{1}{2}}}\rho\epsilon,$ (4.1)

which could be useful to extract information of the radiant energy production
rate $\epsilon$ per unit energy for a hydrostatic stellar system with
radiation. The information of energy production rate becomes important to
build a model in detail [34], particularly to investigate an evolutionary
process of a star as a non-equilibrium open system including dissipation.
Further extension including such a hydrodynamical process is an interesting
future work.

If the relativistic hydrostatic structure equations with a simple fluid is
applied to stellar corona, then pressure vanishes at the far asymptotic
region. This result conflicts with the earlier one obtained by using the non-
relativistic stellar structure equations, and what is a flaw in the earlier
argument is ignoring a term containing general relativistic effect in the
stellar equilibrium equation which cannot be neglected in a situation with
pressure finite in the far asymptotic region. It would be interesting to
investigate this model further to solve problems known in solar corona [34,
35].

In the multilayer stellar model presented in this paper, it consists of two
ideals gases of baryonic particles and non-relativistic ones as a simple
example. It is also possible to pile up a similar layer or another one
consisting of a different type of fluid. Such a multilayer structure should be
arranged to build a model for each observed star. In the analysis, radiation
does not involve in the structure equations but just plays a role of a probe
of local temperature. This approximated treatment of radiation is inevitable
since the proposed hydrostatic structure equations, and thus the classic
stellar ones, cannot be applied to null fluid, strictly speaking, and it is
justified as long as the Eddington bound is sufficiently satisfied. However
the radiative contribution may not be negligible inside the star, in
particular, with so high temperature that an effective method of analysis is
only numerical calculation. In fact, it is expected that there exists a large
radiative zone in the interior of the Sun [34, 35]. In this paper, the
radiative zone has been realized by a simple fluid with $w=1/3$ as a
provisional treatment. It is important to extend the presented relativistic
hydrostatic structure equations so as to include the radiative contribution
and test the validity of the approximation used in this paper by estimating
the radiative contribution.

The application of the proposed relativistic hydrostatic structure equations
has been focused on the construction of a model of luminous stars. It would be
tempted to apply them to a different type of stars such as a compact star
consisting of degenerate matter and to investigate the finite temperature
effect. A caveat for the application of the relativistic hydrostatic structure
equations to such compact stars is to include the local chemical potential and
temperature so as to be consistent with the local thermodynamic relations. It
would be of a great interest to investigate how the thermodynamic relations
shown in [23] including chemical potential in such a system is derived and
whether such thermal interaction gives rise to any significant effect on the
system.

The author hopes to address these issues and will report the progress in the
near future.

## Acknowledgement

The author would like to thank Masaki Mori for a helpful comment on the first
version. This work is supported in part by the Grant-in-Aid of the Japanese
Ministry of Education, Sciences and Technology, Sports and Culture (MEXT) for
Scientific Research (No. JP22K03596).

## References

  * [1] H. L. F. von Helmholtz, “Lxiv. on the interaction of natural forces,” Philosophical Magazine Series 1 11 489–518.
  * [2] W. Kelvin, Popular Lectures and Addresses.  Nature series. Macmillan and Company, 1889.
  * [3] J. N. Bahcall, “How the sun shines,” J. Roy. Astron. Soc. Canada 94 (2000) 219, astro-ph/0009259.
  * [4] A. Eddington, The Internal Constitution of the Stars.  Cambridge University Press/Dover Publications/Cambridge University Press, 1930/1959/1988.
  * [5] G. Gamow, “Nuclear energy sources and stellar evolution,” Phys. Rev. 53 (Apr, 1938) 595–604.
  * [6] H. A. Bethe, “Energy production in stars,” Phys. Rev. 55 (1939) 434–456.
  * [7] M. Harwit, Astrophysical Concepts.  Astronomy and Astrophysics Library. Springer New York, 1988.
  * [8] A. Phillips, The Physics of Stars.  Manchester Physics Series. Wiley, 1994.
  * [9] R. Tayler, The Stars: Their Structure and Evolution.  Cambridge University Press, 1994.
  * [10] C. J. Hansen and S. D. Kawaler, Stellar Interiors.  Springer, 1994.
  * [11] R. Kippenhahn, A. Weigert, and A. Weiss, Stellar Structure and Evolution.  2013\.
  * [12] A. Choudhuri, Astrophysics for Physicists.  Astrophysics for Physicists. Cambridge University Press, 2010.
  * [13] G. P. Kuiper, “The Empirical Mass-Luminosity Relation.,” Astrophysical Journal 88 (Nov., 1938) 472.
  * [14] A. W. Mann, E. Gaidos, and M. Ansdell, “Spectro-thermometry of m dwarfs and their candidate planets: Too hot, too cool, or just right?,” ApJ 779 (2013) 188.
  * [15] M. Cuntz and Z. Wang, “The mass–luminosity relation for a refined set of late-k/m dwarfs,” Research Notes of the AAS 2 (jan, 2018) 19.
  * [16] C. W. Misner, K. S. Thorne, and J. A. Wheeler, Gravitation.  W. H. Freeman, San Francisco, 1973.
  * [17] J. R. Oppenheimer and G. M. Volkoff, “On massive neutron cores,” Phys. Rev. 55 (Feb, 1939) 374–381.
  * [18] G. S. Saakyan and M. A. Mnatsakanyan, “Neutron configurations in the generalized newtonian theory of gravitation.,” Astrophysics (Engl. Transl.) 4: 62-7(Summer 1968).
  * [19] G. S. Saakyan and Y. L. Vartanyan, “Basic Parameters of Baryon Configurations,” Soviet Astronomy 8 (Oct., 1964) 147.
  * [20] J. Hartle, Gravity: An Introduction to Einstein’s General Relativity.  Addison-Wesley, 2003.
  * [21] R. H. Fowler, “On Dense Matter,” Monthly Notices of the Royal Astronomical Society 87 (12, 1926) 114–122, https://academic.oup.com/mnras/article-pdf/87/2/114/3623303/mnras87-0114.pdf.
  * [22] Y. Zel’dovich and I. Novikov, Stars and Relativity.  Dover Books on Physics. Dover Publications, 2014.
  * [23] S. Yokoyama, “Local Thermodynamics and Entropy for Relativistic Hydrostatic Equilibrium,” 2304.06196.
  * [24] S. Aoki, T. Onogi, and S. Yokoyama, “Charge conservation, entropy current and gravitation,” Int. J. Mod. Phys. A 36 (2021), no. 29, 2150201, 2010.07660.
  * [25] R. C. Tolman, “On the weight of heat and thermal equilibrium in general relativity,” Phys. Rev. 35 (Apr, 1930) 904–924.
  * [26] L. Biermann, “Konvektion in rotierenden Sternen,” Z. Astrophys 25 (Jan., 1948) 135.
  * [27] E. Vitense, “Die wasserstoffkonvektionszone der sonne. mit 11 textabbildungen,” Z. Astrophys 32 (1953) 135.
  * [28] E. N. Parker, “Dynamics of the Interplanetary Gas and Magnetic Fields.,” Astrophysical Journal 128 (Nov., 1958) 664.
  * [29] O. Heckmann, “REVIEW: Relativity, Thermodynamics and Cosmology, by R. C. Tolman,” apj 82 (Dec., 1935) 435.
  * [30] S. W. Hawking and G. F. R. Ellis, The Large Scale Structure of Space-Time.  Cambridge Monographs on Mathematical Physics. Cambridge University Press, 2, 2011.
  * [31] S. Aoki, T. Onogi, and S. Yokoyama, “Conserved charges in general relativity,” Int. J. Mod. Phys. A 36 (2021), no. 10, 2150098, 2005.13233.
  * [32] S. Chandrasekhar, “The highly collapsed configurations of a stellar mass (Second paper),” Monthly Notices of the Royal Astronomical Society, Vol. 95, p.207-225 95 (Jan., 1935) 207–225.
  * [33] K. Schwarzschild, “On the equilibrium of the Sun’s atmosphere,” Nachrichten von der Königlichen Gesellschaft der Wissenschaften zu Göttingen. Math.-phys. Klasse 195 (jan, 1906) 41–53.
  * [34] J. N. Bahcall and R. K. Ulrich, “Solar models, neutrino experiments, and helioseismology,” Rev. Mod. Phys. 60 (Apr, 1988) 297–372.
  * [35] NASA, “Solar physics,” https://solarscience.msfc.nasa.gov/interior.shtml.
  * [36] L. Biermann, “Kometenschweife und solare Korpuskularstrahlung,” Z. Astrophys 29 (Jan., 1951) 274.
  * [37] L. Biermann, “Solar corpuscular radiation and the interplanetary gas,” The Observatory 77 (June, 1957) 109–110.
  * [38] K. I. Gringauz, V. V. Bezrokikh, V. D. Ozerov, and R. E. Rybchinskii, “A Study of the Interplanetary Ionized Gas, High-Energy Electrons and Corpuscular Radiation from the Sun by Means of the Three-Electrode Trap for Charged Particles on the Second Soviet Cosmic Rocket,” Soviet Physics Doklady 5 (Sept., 1960) 361.
  * [39] C. W. Snyder and M. Neugebauer, “Interplanetary solar-wind measurements by mariner ii,” in Proceedings of the Plasma Space Science Symposium, C. C. Chang and S. S. Huang, eds., pp. 67–90.  Springer Netherlands, Dordrecht, 1965.
  * [40] D. Verscharen, K. G. Klein, and B. A. Maruca, “The multi-scale nature of the solar wind,” 2019.
  * [41] N. Glendenning, Compact Stars: Nuclear Physics, Particle Physics and General Relativity.  Astronomy and Astrophysics Library. Springer New York, 2012.
  * [42] O. Pols, Stellar Structure and Evolution.  Astronomical Institute Utrecht, 2011.
  * [43] S.-M. Kung, “The opacity and the internal structure of the Sun,” Commmunications of the Konkoly Observatory Hungary 42 (Jan., 1957) 109\.
  * [44] S. Shapiro and S. Teukolsky, Black Holes, White Dwarfs, and Neutron Stars: The Physics of Compact Objects.  Wiley, 2008.
  * [45] D. Prialnik, An Introduction to the Theory of Stellar Structure and Evolution.  Cambridge University Press, 2000.
  * [46] B. Carroll and D. Ostlie, An Introduction to Modern Astrophysics.  Pearson Addison-Wesley, 2007.
  * [47] J. P. Cox, “Principles of stellar structure. volume 2. applications to stars.,”.


# Active microrheology of fluids with orientational order

J. S. Lintuvuori  $|$ , A. Würger<EMAIL_ADDRESS><EMAIL_ADDRESS>
(Received July 19, 2023, in final form November 9, 2023)

###### Abstract

We study the dynamics of a driven spherical colloidal particle moving in a
fluid with a broken rotational symmetry. Using a nematic liquid crystal as a
model, we demonstrate that when the applied force is not aligned along or
perpendicular to the orientational order, the colloidal velocity does not
align with the force, but forms an angle with respect to the pulling
direction. This leads to blue an anisotropic hydrodynamic drag tensor which
depends on the material parameters. In the case of nematic liquid crystal, we
give an analytical expression and discuss the resulting implications for
active microrheology experiments on fluids with broken rotational symmetry.

Key words: liquid crystals, colloids, microrheology, lattice Boltzmann methods

###### Abstract

Äîñëäæóòüñÿ äèíàìêà êåðîâàíî¿ ñôåðèíî¿ êîëî¿äíî¿ àñòèíêè, ùî ðóõàòüñÿ â ðäèí ç
ïîðóøåíîþ îáåðòàëüíîþ ñèìåòðþ. Âèêîðèñòîâóþè â ÿêîñò ìîäåë íåìàòèíèé ðäêèé
êðèñòàë, ïîêàçàíî, ùî êîëè ïðèêëàäåíà ñèëà íå ä âçäîâæ àáî ïåðïåíäèêóëÿðíî
îðíòàöéíîìó ïîðÿäêó, øâèäêñòü àñòèíêè íå êîëíåàðíîþ ç íàïðÿìêîì ñèëè, à óòâîðþ
êóò âäíîñíî íàïðÿìêó ðîçòÿãóâàííÿ. Öå ïðèçâîäèòü äî ‘‘ñèíüîãî’’ àíçîòðîïíîãî
ãäðîäèíàìíîãî òåíçîðà îïîðó, ÿêèé çàëåæèòü âä ïàðàìåòðâ ìàòåðàëó. Ó âèïàäêó
íåìàòèíîãî ðäêîãî êðèñòàëà ìè íàâîäèìî éîãî àíàëòèíèé âèðàç îáãîâîðþìî íàñëäêè
äëÿ àêòèâíèõ ìêðîðåîëîãíèõ åêñïåðèìåíòâ íà ðäèíàõ ç ïîðóøåíîþ îáåðòàëüíîþ
ñèìåòðþ.

Ключов слова: ðäê êðèñòàëè, êîëî¿äè, ìêðîðåîëîãÿ, ðàòêîâ ìåòîäè Áîëüöìàíà

##  1 Introduction

Understanding the flow properties of complex fluids is a challenging goal.
Various techniques have been developed to probe the hydrodynamic response of
the fluids to external stimuli. An example of this is provided by shear
rheology, which allows the measurement of the bulk viscosity $\eta$. Another
modern technique is called microrheology, which consists of measuring the
hydrodynamic drag experienced by a spherical colloidal particle moving in the
fluid due to external [1] or thermal [2] forces. In the case of an external
force, the resulting force-velocity curves can be analyzed by using Stokes’
law [3]: a particle of size $R$ moving at velocity $v$ in a fluid of viscosity
$\eta$ experiences an hydrodynamic drag proportional to these quantities,
$-\xi v$ with $\xi=6\piup\eta R$, which balances the external force $F$.

The situation is markedly different when the host fluid has a broken
rotational symmetry, such as the orientational order observed in liquid
crystals [4], sheared polymer solutions [5, 6, 7] and various biological
fluids [8, 9, 10]. Then, the viscous response is described by a fourth-order
tensor $\eta_{ijkl}$ acting on the stress components $\sigma_{kl}$ [4], and
the resulting hydrodynamic drag force on a spherical bead does not necessarily
align with the velocity $\mathbf{v}$ of the particle [11]. In other words, the
drag coefficient $\xi_{i}$ becomes a second-order tensor.

A typical example of fluids with orientational order is provided by
thermotropic liquid crystals in a nematic state. Here, the molecules are, on
average, aligned along a common direction, characterised by a vector field
$\hat{\mathbf{n}}$ with a head-to-tail symmetry. Recent experimental,
theoretical and simulations studies have quantified the drag coefficient
$\xi_{i}=6\piup\eta_{i}R$ of a spherical colloid, parallel and perpendicular
to $\hat{\mathbf{n}}$, and observed
$\frac{\xi_{\perp}}{\xi_{||}}=\frac{\eta_{\perp}}{\eta_{||}}\sim 2$ [12, 13,
11, 14, 15, 16, 17]; this relation also defines effective viscosities
$\eta_{\perp}$ and $\eta_{\parallel}$.

Here, we study a general case of a spherical particle (radius $R$) dragged
through a nematic liquid crystal (LC) with an angle $\theta$ between the
pulling force $\mathbf{F}$ and the nematic director $\hat{\mathbf{n}}$. We
show that when the pulling direction is between the limiting case of parallel
($\theta=0^{\circ}$) and perpendicular ($\theta=90^{\circ}$) to the nematic
direction, the force and velocity do not align, but a sliding motion is
observed. This is characterised by an angle $\alpha$ between the vectors
$\mathbf{F}$ and $\mathbf{v}$.

Interpreting this by using the general force-velocity relation via Stokes’ law
leads to an effective drag response (and thus viscosity) which depends both on
the $\theta$ and the ratio of the drag coefficients parallel and perpendicular
to the nematic director. We give analytical expressions for both the
generalized drag coefficient $\xi_{\mathrm{eff}}=-F/v$ and the sliding angle
$\alpha$, and validate them by lattice Boltzmann simulations.

##  2 Methods

The nematic order is described by a symmetric and traceless tensor order
parameter $\mathbf{Q}$, which time evolution follows a hydrodynamic equation
[18]

$(\partial_{t}+u_{\nu}\partial_{\nu})Q_{\alpha\beta}-S_{\alpha\beta}=\Gamma
H_{\alpha\beta},$ (2.1)

where the first part describes the advection and $S_{\alpha\beta}$ describes
the possible rotation/stretching of $\mathbf{Q}$ by the flow [18]. $\Gamma$ is
the rotational diffusion constant. The molecular field is

$H_{\alpha\beta}=-{\delta{\cal F}/\delta
Q_{\alpha\beta}}+(\delta_{\alpha\beta}/3){\mbox{\rm Tr}}({\delta{\cal
F}/\delta Q_{\alpha\beta}}),$ (2.2)

where ${\cal F}$ is a Landau-de Gennes free-energy whose density can be
expressed in terms of a symmetric and traceless order parameter tensor
$\mathbf{Q}$ as ${\cal
F}=F(Q_{\alpha\beta})+\tfrac{K}{2}(\partial_{\beta}Q_{\alpha\beta})^{2}$, with

$F(Q_{\alpha\beta})=A_{0}\left(1-\frac{\gamma}{3}\right)\frac{Q_{\alpha\beta}^{2}}{2}-\frac{\gamma}{3}Q_{\alpha\beta}Q_{\beta\gamma}Q_{\gamma\alpha}+\frac{\gamma}{4}(Q_{\alpha\beta}^{2})^{2},$
(2.3)

where Greek indices denote Cartesian coordinates and summation over repeated
indices is implied. $A_{0}$ is a free energy scale, $\gamma$ is a temperature-
like control parameter giving a order/disorder transition at $\gamma\sim 2.7$,
and $K$ is an elastic constant. The anchoring at the particle surface is
modelled by $f_{s}=W(Q_{\alpha\beta}-Q^{0}_{\alpha\beta})^{2}$, where $W$ is
the anchoring strength and $Q^{0}_{\alpha\beta}$ is the preferred alignment of
the nematic director at the particle surface.

Figure 1: (Colour online) (a) Schematic of the system. (b) The observed steady
state trajectories of the particles in the $x-z$ space (symbols) in the
absence of surface anchoring ($WR/K=0$). The solid lines show the direction of
the applied force. When $0^{\circ}<\theta<90^{\circ}$, the force and velocity
are not aligned, but the particles move on a straight trajectory with an angle
$\alpha$ with respect to the applied force.

The fluid velocity obeys the continuity $\partial_{\alpha}u_{\alpha}=0$, and
the Navier-Stokes equations. These are coupled to the LC via a stress tensor.
We employ a 3D lattice Boltzmann algorithm to solve the equations of motion
(for further details see e.g. [16, 19]).

The colloids are modelled as spherical particles, with a no-slip boundary
condition. The no-slip boundary condition at the fluid/solid interface is
realized by a standard method of bounce-back on links (BBL) [20, 21] and it
can be modified to take into account the moving particle surface [22].

Figure 2: (Colour online) An example of the active microrheological
measurement of particles with no surface anchoring ($WR/K=0$), corresponding
to the sample in figure 1. The particles were pulled in the nematic LC with a
force $F$ at an angle $\theta$ with respect to the nematic director
$\hat{\mathbf{n}}$, and the magnitude $v$ of the velocity $\mathbf{v}$ was
measured. The solid lines are linear fits $F=\xi v$.

The dynamics of our system is governed by the Reynolds (Re) and Ericksen (Er)
numbers, measuring the ratios between inertial and viscous forces as well as
viscous and elastic forces, respectively. Using our simulation parameters
111Parameters were (in simulation units): $A_{0}=0.1,~{}K\simeq
0.01,~{}\xi=0.7,~{}\gamma=3.0$, $q=1/2$ , $\Gamma=0.3$, $W=0$ or $0.01$ and a
particle radius $R=4.0$. These give $\tau=0$, and
$\gamma_{1}=\frac{2q^{2}}{\Gamma}=\frac{5}{3}$. A cubic simulation box
$16R\times 16R\times 16R$ with periodic boundary conditions, was used., we get
the upper limits of $\mathrm{Re}=\frac{vR}{\eta}\approx 0.03$ and
$\mathrm{Er}=\frac{\gamma_{1}vR}{K}\approx 1.5$, where $\gamma_{1}$ is the
rotational viscosity of the nematic LC.

##  3 Results

Figure 3: (Colour online) The difference between the sliding angle $\alpha$
and the pulling angle $\theta$ as a function of $\theta$. The symbols
correspond to LB simulations and the solid lines are theoretical predictions
(equation 3.3) calculated using the ratio $\frac{\xi_{\perp}}{\xi_{||}}$ from
figure 4.

We study the hydrodynamic drag of a spherical particle pulled through a
nematic liquid crystal with a force $\mathbf{F}$ on an angle $\theta$ with
respect to the nematic director $\hat{\mathbf{n}}$ (figure 1a). In addition to
the situation where particle surface does not influence the orientation of the
nearby liquid crystal molecules ($WR/K=0$), we also consider a strong
homeotropic or degenerate planar anchoring conditions ($WR/K=4$), where the
director adopts either normal or planar orientation with respect to the
particle surface, giving a rise to Saturn ring or boojum defects,
respectively.

When the force is aligned along ($\theta=0^{\circ}$) or perpendicular
($\theta=90^{\circ}$) to $\hat{\mathbf{n}}$, the particle moves along the
force with a speed $v$. This allows for a microrheological measurement of the
drag coefficients along and perpendicular to the nematic order, where
$\frac{\xi_{\perp}}{\xi_{||}}=\frac{\eta_{\perp}}{\eta_{||}}\sim 2$ is
observed (figure 2 and 4) in agreement with the previous studies [12, 13, 14,
15, 16, 17].

In the general case, when the external force makes an angle $\theta$ with the
nematic director $\mathbf{F}\cdot\hat{\mathbf{n}}=F\cos\theta$, the velocity
components parallel and perpendicular to the director satisfy the relations

$\xi_{||}v_{||}=F\cos\theta,\quad\xi_{\perp}v_{\perp}=F\sin\theta.$ (3.1)

Unless the drag coefficients are equal $\xi_{||}=\xi_{\perp}$, the velocity is
not parallel to the force, but forms an angle $\alpha$ with the pulling
direction (figure 1b), leading to a sliding motion towards the lower viscosity
direction ($x$-axis in our case).

The force-velocity relation is linear (figure 2), which is a minimum
requirement for an active microrheological measurement. From a simple
geometrical considerations, one finds

$\tan\left(\theta-\alpha\right)=\frac{v_{\perp}}{v_{||}}=\frac{\xi_{||}}{\xi_{\perp}}\tan\theta$
(3.2)

or

$\alpha=\theta-\arctan\left(\frac{\xi_{||}}{\xi_{\perp}}\tan\theta\right),$
(3.3)

which agree very well with the simulations (figure 3).

Figure 4: (Colour online) The observed normalized effective viscosity
$\frac{\eta_{\mathrm{eff}}}{\eta_{||}}$ as a function of $\theta$ from the
simulations (symbols). The solid lines are the theoretical predictions
(equation 3.4) in the text plotted using $\xi_{||}$ and $\xi_{\perp}$ from the
LB simulations.

The linear velocity $v=\sqrt{v^{2}_{||}+v^{2}_{\perp}}$ defines the effective
drag coefficient through

$\frac{1}{\xi_{\mathrm{eff}}(\theta)}=\frac{v}{F}=\sqrt{\frac{\cos^{2}\theta}{\xi^{2}_{\perp}}+\frac{\sin^{2}\theta}{\xi^{2}_{||}}},$
(3.4)

and thus, an angle dependent effective viscosity
$\eta_{\mathrm{eff}}=\xi_{\mathrm{eff}}/6\piup R$.

##  4 Conclusions

Using the lattice Boltzmann simulations, we carried out microrheological
experiments of colloidal particles in a nematic liquid crystal. Both the angle
$\theta$ and the surface properties of the particles were varied. In all the
cases, we observed a linear relation between the magnitude $F$ of the imposed
force and the observed speed $v$. The drag coefficients are trivially
available from linear fits $F=\xi v$ (see e.g. figure 2). The simulations
confirmed that the sliding motion towards the lower viscosity direction, leads
to a non-trivial flow response (figure 4). As can be expected from the sliding
angle (figure 3), the effect is most pronounced for colloids with a strong
homeotropic anchoring at their surface followed by particles with a planar
surface anchoring. This can be understood in terms of the strength of coupling
between the particle surface and the local nematic order. The effective drag
coefficients follow, to a very high degree, the theoretical predictions of
equation (3.4), for all the three surface anchorings considered (figure 4).

## References

  * [1] Squires T. M., Mason T. G., Annu. Rev. Fluid Mech., 2010, 42, 413, doi:10.1146/annurev-fluid-121108-145608.
  * [2] Manon T. G., Rheol. Acta, 2000, 39, 371, doi:10.1007/s003970000094.
  * [3] Landau L. D., Lifshitz E. M., Fluid Mechanics, Pergamon Press, Oxford, 1987.
  * [4] De Gennes P.-G., Prost J., The Physics of Liquid Crystals, 2nd Ed., Clarendon Press, Oxford, 1993.
  * [5] Larson R. G., The Structure and Rheology of Complex Fluids, Oxford University Press, New York, 1999.
  * [6] Pujolle-Robic C., Noirez L., Nature, 2001, 409, 167, doi:10.1038/35051537 .
  * [7] Parisi D., Seo J., Schaake R. P., Rhoades A. M., Colby R. H., Prog. Polym. Sci., 2021, 112, 101323, doi:10.1016/j.progpolymsci.2020.101323.
  * [8] Doostmohammadi A., Ladoux B., Trends Cell Biol., 2022, 32, 140, doi:10.1016/j.tcb.2021.09.012.
  * [9] Broedersz C. P., MacKintosh F. C., Rev. Mod. Phys., 2014, 86, 995, doi:10.1103/RevModPhys.86.995.
  * [10] Nishiguchi D., Nagai K. H., Chaté H., Sano M., Phys. Rev. E, 2017, 95, 020601, doi:10.1103/PhysRevE.95.020601.
  * [11] Stark H., Phys. Rep., 2001, 351, 387, doi:10.1016/S0370-1573(00)00144-7.
  * [12] Ruhwandl R. W., Terentjev E. M., Phys. Rev. E, 1996, 54, 5204, doi:10.1103/PhysRevE.54.5204.
  * [13] Stark H., Ventzki D., Phys. Rev. E, 2001, 64, 031711, doi:10.1103/PhysRevE.64.031711.
  * [14] Loudet J. C., Hanusse P., Poulin P., Science, 2004, 306, 1525, doi:10.1126/science.1102864.
  * [15] Gleeson H., Wood T. A., Dickinson M., Phil. Trans. R. Soc. A, 2006, 364, 2789, doi:10.1098/rsta.2006.1855.
  * [16] Lintuvuori J. S., Marenduzzo D., Stratford K., Cates M. E., J. Mater. Chem., 2010, 20, 10547, doi:10.1039/c0jm01824d.
  * [17] Mondiot F., Loudet J.-C., Mondain-Monval O., Snabre P., Vilquin P. A., Würger A., Phys. Rev. E, 2012, 86, 010401, doi:10.1103/PhysRevE.86.010401.
  * [18] Beris A. N., Edwards B. J., Thermodynamics of Flowing Systems, Oxford University Press, Oxford, 1994.
  * [19] Cates M. E., Henrich O., Marenduzzoa D., Stratfordb K., Soft Matter, 2009, 5, 3791, doi:10.1039/b908659p.
  * [20] Ladd A. J. C., J. Fluid Mech., 1994, 271, 285, doi:10.1017/S0022112094001771.
  * [21] Ladd A. J. C., J. Fluid Mech., 1994, 271, 311, doi:10.1017/S0022112094001783.
  * [22] Nguyen N.-Q., Ladd A. J. C., Phys. Rev. E, 2002, 66, 046708, doi:10.1103/PhysRevE.66.046708.

Ukrainian 0

Àêòèâíà ìêðîðåîëîãÿ ïëèíâ ç îðíòàöéíèì ïîðÿäêîì [É. Ñ. Ëíòóâóîð, À. Âþðãåð] É.
Ñ. Ëíòóâóîð, À. Âþðãåð

Óíâåðñèòåò Áîðäî, CNRS, UMR 5798, F-33400 Òàëåíñ, Ôðàíöÿ