# Structure of Stellar Remnants with Coupling to a Light Scalar

Christina Gao  and Albert Stebbins

###### Abstract

In this paper we study how a Yukawa coupling of the Standard Model fermions to
a light scalar field effects the stellar structure of cold stellar remnants
such as neutron stars. We elucidate the stellar structure phenomenology using
a simple model of a massive scalar coupled to a single dominant fermion with
no other interactions. For a broad scalar mass range
($10^{-10}\,\mathrm{eV}\ll m_{\phi}\ll 10^{3}\,\mathrm{eV}$ for neutron stars)
we show that the equation-of-state and stellar structure depends only the
effective coupling $\mathfrak{g}=\frac{g_{f}\,m_{f}}{m_{\phi}}$, where $g_{f}$
is the Yukawa coupling, $m_{f}$ the fermion mass and $m_{\phi}$ is the scalar
kinematic mass at nuclear densities. If $\mathfrak{g}>\mathcal{O}(1)$ the
Yukawa coupled matter exhibits various anomalous behaviors including
hydrodynamic instability, negative pressure, distinct phases (soft and hard)
of matter with sharp phase boundaries between them and with vacuum. These
anomalies can lead to stars consisting of only soft, only hard or hybrid of
soft and hard matter. These stars can have varying signs of the slope of the
mass-radius relation, anomalously large and small masses, gaps in allowed
radii, multiple radii for the same mass, very thin crusts and radiate
anomalously large amounts of energy when they form (in the form of neutrinos
for neutron stars). To the extent that these anomalies have not and/or will
not be observed limits the effective coupling to
$\mathfrak{g}<\mathcal{O}(1)$. We argue this phenomenology is generic to
realistic models of stars with Yukawa coupled matter.

##  1 Introduction

A light scalar field is one of the simplest extension to the Standard Model
(SM) but is only one of the plethora of proposed extensions referred to
collectively as Beyond the Standard Model (BSM) physics. A scalar can act as a
“mediator” between the SM and the “dark sector” of particles which is believed
to be associated with a putative dark matter particle, i.e. it couples
directly to SM and dark sector particles. The direct SM coupling could make a
mediator easier to detect than other dark sector particles. Such a mediator
will also give a non-SM (BSM) coupling of SM matter with itself which is the
topic of this paper. We make no reference to a dark sector or to dark matter.
Here we consider the effect of a linear scalar field coupling to SM particles
on the structure of stellar remnants.

Neutron stars provide a unique astrophysical laboratory for both SM and BSM
physics in the combination of properties

  1. 1.

like many astronomical objects, gravity determines their structure, which
provides a probe of BSM interactions of matter at the gravitational strength
in addition to or in place of the standard gravity (GR / general relativity),

  2. 2.

their small size ($\sim 10\,$km) provides a probe of these interactions on
smaller length scales than other astronomical objects

  3. 3.

their structure on these small scales is well constrained by a variety of very
different types of observations,

  4. 4.

they consist of bulk quantities of nuclear density matter (the highest density
in the present universe) providing a probe of SM and BSM interactions at these
densities as well as their density dependence,

  5. 5.

their SM equation of state (EoS) is (ideally) simple with essentially no free
parameters providing a foil against which to test BSM effects.

While there may be no free parameters in the SM EoS, there is currently a
great deal of uncertainty in the EoS of cold matter at the nuclear densities
due to difficulties in accounting for the SM hadronic interactions. The
uncertainty in the SM EoS leads to the uncertainty in neutron star structure
and thus limits observational sensitivities to BSM effects from the
observations of neutron star structure.

The simple EoS property derives from a short cooling time ($\sim 10\,$sec by
neutrinos) allowing the bulk111The state of the surface (crust) of a neutron
star may be impacted by the exterior environment and history. to cool well
below the Fermi energy of its constituent particles and therefore $T=0$ is a
good approximation. Before cooling, the bulk of the stellar matter was hot and
dense enough to approach nuclear statistical equilibrium (NSE) which means
that the composition is known and described by a single baryon chemical
potential222Conductivity by $e^{-}$, $p^{+}$, etc. leads to a negligible
charge separation.. At nuclear densities, neutron star consists primarily of
neutrons. In a given SM or BSM context, a neutron star is like a black hole in
that it is well described by only its mass and angular momentum. Unlike black
hole stellar remnants, the observed and expected rotation rate of neutron
stars is small enough to have little effect on their structure. While black
holes have “no hair”, neutron stars can and do; most importantly the embedded
magnetic fields. However these fields are weak enough to have even less effect
on the structure of the star than angular momentum. To $\lesssim 10\%$
accuracy, mass is the only parameter of neutron star structure, although
angular momentum333Neutron stars used to constrain BSM physics would have to
have known rotation rates. and to a lesser extent magnetic fields would become
important if we could determine its structure with greater precision. The
stars are thus nearly spherically symmetric and their structure is given by
their density profile or $M(R)$. This paper focuses on determining $M(R)$ when
normal matter is coupled to a scalar field.

White dwarfs share much of the simplicity of neutron stars. While white dwarfs
are not in NSE and their chemical composition will vary, there is relatively
little composition dependence of their EoS. In §7.1 we briefly discuss how our
findings for neutrons stars apply almost directly to white dwarfs.

###  1.1 Overview

In this paper we examine how the existence of a Yukawa coupling to a BSM
scalar field will effect the structure of stellar remnants. It is not assumed
that this scalar is present in any significant abundance in the universe, say,
as dark matter. The paper is organized in the following way. In §2 we describe
the coupling under consideration and the level of approximation with which we
treat it. In §3 we review stellar structure and the ingredient needed to
determine this structure which is the equation of state (EoS). In §4 we
formulate the EoS including a scalar coupled to an ideal gas of free fermions.
Specifically in §4.3 we show how the scalar is polarized to non-zero field
values in the presence of fermions. In §4.4 we identify three regimes for
computation of this polarization, heavy, light and ultralight which depends on
the mass of the scalar and the density of fermions. Specializing to the most
interesting light scalar regime in §4.5.1 we give analytic expressions for the
ideal gas +Yukawa, or IGY EoS which for a single fermion depends only on the
effective coupling $\mathfrak{g}_{f}=\frac{g_{f}\,m_{f}}{m_{\phi}}$ ($g_{f}$
is the Yukawa coupling, $m_{f}$ the fermion mass and $m_{\phi}$ the scalar
mass). In §5 we highlight anomalies in the IGY EoS which appear when
$\mathfrak{g}_{f}$ exceeds a critical value. These anomalies lead to different
matter phases depending on the matter density: a range of densities forbidden
in stars separates a low density soft matter phase from a higher density hard
matter phase. The soft phase can exist in a metastable state and all low
density matter may be metastable to conversion to hard matter if
$\mathfrak{g}_{f}$ is large enough. For a small range of $\mathfrak{g}_{f}$ a
stable phase boundary can exist between hard and soft matter (§5.2.1) and for
larger $\mathfrak{g}_{f}$ a stable boundary between hard matter and vacuum
(§5.2.2). IGY stars, which are spherical hydrostatic configurations of matter
with a single fermion (the neutron) IGY EoS provide a prototype for neutron
stars with a Yukawa coupling between a scalar and baryons. This simple model
is only considered a prototype as it has a number of deficiencies listed in
§8. This phenomenology is described in §6 showing their masses and radii in
§6.1-§6.4 and selected density profiles in §6.5. Above the critical coupling
stars may consist purely of hard matter or of soft matter or a hybrid of the
two; and some stars may be metastable. Very high density but very low mass
hard matter configurations, bound by the scalar field and not by gravity,
exist which we call hard drops (§6.2) which could constitute the dark matter.
In §6.6 we explore the energetics of collapse to an IGY star and the
implications for the energy radiated during formation. In §6.7 we list some of
the anomalous features of the IGY stellar structure relative to that of SM
stars when $\mathfrak{g}_{f}$ is large enough: the mass may increase or
decrease with radius; the maximum mass stable star may or may not have the
minimum radius; there may be a gap in allowed stellar radii; there are stable
configurations with nuclear densities but extremely small radii and masses;
stars much larger in size and mass than the SM neutron stars are possible; the
crust of an IGY star may be very thin; neutrino bursts from IGY stars may be
very large. We argue in §7 that the qualitative phenomenology found for IGY
stars will be common to most models of stellar remnants with matter coupled to
scalars. In §7.1 we illustrate this point with the example of white dwarfs. In
§9 we discuss implications of our results for detecting or excluding a Yukawa
coupling from observations of neutron stars and white dwarfs:
$\mathfrak{g}_{f}<\mathcal{O}(1)$ if no anomalous behavior is observed where
$f$ is fermion which dominates the pressure. More realistic modeling of
Yukawa-coupled neutron stars and white dwarfs and comparison to data are left
to future works which will set bounds on Yukawa couplings from astronomical
observations of stellar structure and compare these with bounds obtained by
other methods.

The EoS phenomenology developed here has much in common with that developed
for mass-varying neutrinos [1, 2].

###  1.2 Conventions and Notation

Neutron stars are mildly relativistic so our calculations should account for
space-time curvature. We assume Einstein metric gravity and use the
$\\{+,-,-,-\\}$ metric signature convention, Einstein summation convention and
;μ for covariant derivatives. We will use HEP natural units where $\hbar=c=1$
except where they appear explicitly.

##  2 A Scalar Linearly Coupled To Matter

We posit the existence of a light BSM real scalar field $\phi$ which has a
linear coupling to quarks, leptons, gluons and/or photons444Here for
simplicity the coupling is assumed diagonal between different quarks and
gluons.

$\displaystyle\mathcal{L}\supset\left(\sum_{q}g_{q}\,\bar{\psi}_{q}\,\psi_{q}+\sum_{\ell}g_{\ell}\,\bar{\psi}_{\ell}\,\psi_{\ell}+g_{\mathrm{g}}\,\sum_{a}G^{\mu\nu}_{a}\,G^{a}_{\mu\nu}+g_{\mathrm{\gamma}}\,F^{\mu\nu}\,F_{\mu\nu}\right)\,\phi$
(2.1)

where $q$ / $\ell$ / $a$ sum over the different quarks / leptons / gluons,
$\psi_{x}$ are the fermion spinors and $G^{a}_{\mu\nu}$ / $F_{\mu\nu}$ are the
gluon / electromagnetic field tensors. $g_{q}$, $g_{\ell}$ are dimensionless
Yukawa coupling constants. $g_{\mathrm{g}}$ and $g_{\gamma}$ are scalar
couplings with dimension mass-1. We define $\phi$ such that $\phi=0$ gives the
minimum energy in the absence of matter. We assume minimal coupling of $\phi$
to gravity, e.g. no significant $\phi\,R$ term in $\mathcal{L}$. The couplings
are taken as free parameters but clearly must be extremely small for $\phi$ to
have escaped detection. The main goal of this paper is to explore how such
couplings modify stellar structure if such a $\phi$ exists.

If $\phi$ were a constant then one could absorb the $g$’s by suitable field
redefinitions and adjustment of mass matrices and gauge couplings. However we
show that the value of $\phi$ will vary according to the macroscopic
environment (“finite density” effects) resulting in environment dependent SM
particle parameters and properties which differ from the near vacuum values
measured in the laboratory. This is analogous to the Higgs mechanism in SM
electro-weak symmetry breaking in which particles attain mass. Here however
there is no phase transition, the environment dependence is gradual and in
static stars purely spatial. Environment dependence plays a key role in
modified gravity / 5th-force models through a chameleon mechanism which is
analogous but different to that considered here. $\phi$ can contribute to a
long range (5th) force but we will focus on parameter regimes where such
forces do not extend significantly beyond the stellar surface.

Astrophysical phenomena are usually formulated in terms of a few of the
fundamental particles of eq. 2.1 plus other strongly or weakly555meaning the
binding energy much less than the mass of the constituent particles bound
composites of the other elementary constituents. These composite particles
will inherit a coupling to $\phi$ through their constituents. To 1st order in
$\phi$ the form of the effective interaction with $\phi$ for composite
particles is dictated by Lorentz invariance which for unpolarized fermions is

$\displaystyle\mathcal{L}_{\mathrm{eff}}\supset\phi\,\sum_{f}\,g_{f}\,\bar{\psi}_{f}\,\psi_{f}$
(2.2)

where $\psi_{f}$ is the spinor of fermion $f$. In an astrophysical context
particles are rarely significantly polarized, and even if not, a spin coupling
is likely a subdominant effect.666 A non-zero $\bar{\psi}_{f}\,\psi_{f}$ is
analogous to a bulk space-charge density in electromagnetism without
oppositely signed charges to shield long range fields. Bulk space charges
leads to much larger field strengths than spin couplings which, even for fully
polarized materials, contributes only surface and gradient terms as, say, in
magnetic materials. For weakly bound composites such as multi-nucleon nuclei
and atoms the coupling constant will approximately be the sum of the couplings
of the constituents, e.g. for nuclei $g_{(Z,A)}\approx
Z\,g_{\mathrm{p}}+(A-Z)\,g_{\mathrm{n}}$ and neutral atoms add
$+Z\,g_{\mathrm{e}}$. The three couplings, $g_{\mathrm{p}}$, $g_{\mathrm{n}}$,
$g_{\mathrm{e}}$, are the most important for stellar structure.

In principle $g_{\mathrm{p}}$ and $g_{\mathrm{n}}$ are related to the more
fundamental couplings of eq. 2.1. However hadrons are strongly bound
composites of quarks and gluons with no simple relation between $g_{h}$
($h\in$ nucleons, mesons, etc.) and the more fundamental couplings.
Determining this relation would require knowledge of $\phi$-dependent hadronic
structure. This structure has not yet been fully understood even in the SM and
is outside the scope of this paper. The relative size of $g_{h}$ depends on
the relative contribution of quark and gluon couplings which adds further
uncertainty. Nuclear interactions which are important in neutron stars are
mediated by mesons which are mostly pseudo-scalars and cannot have linear
couplings, e.g. $\pi^{0}\,\phi$, as in eq. 2.2 though they can have other
couplings e.g. quadratic $\pi^{0}\,\pi^{0}\,\phi$. We leave discussion of
nuclear interactions to a follow-on paper.

###  2.1 Effective Mean Field Theory

The $\phi$ field as a force carrier whose range is given by the mass of $\phi$
which if small will coherently interact with a very large number, say $N$, of
SM particles. This effectively boosts the interaction energy by a factor of
$N$ over single particle incoherent interactions. Small mass or large $N$ can
lead to interaction energy between $\phi$ and the SM matter exceeding SM
forces even for extremely small couplings. SM forces do not have large
coherent boosts either because of the large mass of the force carriers in the
case of nuclear interactions or charge screening in the case of
electromagnetism. Due to the finite range of the force the boost factor $\sim
N$ increases as the density so manifestations of $\phi$ can dominate at the
large densities of white dwarfs or neutron stars while being very small at
laboratory densities. In vacuum, $\phi$ may fluctuate around zero, but in a
high density environment with large $N$, the lowest energy state will have
very large $\phi$ values, thus in a well cooled and equilibrated star, these
non-zero values, $\tilde{\phi}\equiv\langle\phi\rangle$, will be attained.
Stars are not uniform and $\tilde{\phi}$ will vary within the star. We treat
the SM particles in the same mean field approximation supposing they have
reached some local equilibrium. This mean field approach does not reveal
certain phenomena such as $\phi$ production from inelastic SM particle
scatterings, but these processes are suppressed by $N$ and in most cases of
interest are negligible. Scattering production of $\phi$’s can also be highly
suppressed by Fermi blocking.

The effective action of the system under study is

$\displaystyle S=\int
d^{4}x^{\mu}\,\sqrt{-\mathrm{g}}\,\mathcal{L}_{\mathrm{eff}}\qquad\mathcal{L}_{\mathrm{eff}}=\mathcal{L}_{\mathrm{SM}}+\tilde{g}\,\tilde{\phi}\,\tilde{n}+\mathcal{L}_{\tilde{\phi}}$
(2.3)

where $\mathcal{L}_{\mathrm{SM}}$ is the standard model Lagrangian density,
$\mathcal{L}_{\tilde{\phi}}=\frac{\tilde{\phi}^{;\mu}\,\tilde{\phi}_{;\mu}}{2}-V(\tilde{\phi})$,
and $\mathrm{g}$ is the determinant of the metric tensor (not to be confused
with $\tilde{g}$ which is an unspecified weighted average of the couplings
previously described), $\tilde{n}$ is the corresponding sum of mean field
values, $\langle\bar{\psi}\,\psi\rangle$’s. When a single free fermion
species, say neutrons, dominates, one should choose $\tilde{g}=g_{\mathrm{n}}$
and $\tilde{n}=\langle\bar{\psi}_{\mathrm{n}}\,\psi_{\mathrm{n}}\rangle$.

Note that $\tilde{n}$ has dimensions of particle number density and is related
to the particle number density of nucleons as we shall show. Unlike the frame
dependent volume number density, $\tilde{n}$ is a Lorentz scalar and does not
invoke any reference frame (4-velocity $u^{\mu}$) such as the center of
momentum frame of a fluid. Therefore, there is no possibility of a conserved
current ($(\tilde{n}\,u^{\mu})_{;\mu}=0$) even when particle numbers are
conserved.

##  3 Stellar Structure Equations

The structure of a star is given by the stress-energy $T_{\mu\nu}$. The
Hilbert stress-energy777which may be generalized to include spinor fields,
$T_{\mu\nu}=2\,\frac{\delta\mathcal{L}}{\delta\mathrm{g}^{\mu\nu}}-\mathrm{g}_{\mu\nu}\,\mathcal{L}$,
may be decomposed as

$\displaystyle
T_{\mu\nu}=T^{\mathrm{SM}}_{\mu\nu}+2\,\tilde{g}\,\tilde{\phi}\,\frac{\delta\tilde{n}}{\delta\mathrm{g}^{\mu\nu}}-\tilde{g}\,\tilde{\phi}\,\tilde{n}\,\mathrm{g}_{\mu\nu}+T^{\tilde{\phi}}_{\mu\nu}\
.$ (3.1)

Spherical static stars have time-translation and 3-D rotational isometries
which dictate the form of the $\mathrm{g}_{\mu\nu}$ and $T_{\mu\nu}$, which
may be written in spherical polar coordinates $\\{t,r,\theta,\phi\\}$ as

$\displaystyle\mathrm{g}_{\mu\nu}$ $\displaystyle=$
$\displaystyle\mathrm{diag}\left(e^{2\Phi(r)},\frac{-1}{1-\frac{2\,G\,M(r)}{r}},-r^{2},-r^{2}\,\sin^{2}\theta\right)$
$\displaystyle{T^{\mu}}_{\nu}$ $\displaystyle=$
$\displaystyle\mathrm{diag}\left(\mathcal{E}(r),-P_{\parallel}(r),-P_{\perp}(r),-P_{\perp}(r)\right)$
(3.2)

where $r$ is the circumferential radius, $\mathcal{E}$, $P_{\parallel}$ and
$P_{\perp}$ are respectively the energy density, radial and tangential
pressure in the stellar rest frame. The two metric functions $M(r)$ and
$\Phi(r)$ were chosen to be similar to their Newtonian analogs. Einstein’s
equations may be written as

$\displaystyle M^{\prime}(r)$ $\displaystyle=$ $\displaystyle 4\pi
r^{2}\mathcal{E}(r)$ (3.3) $\displaystyle\Phi^{\prime}(r)$ $\displaystyle=$
$\displaystyle\frac{G\,(M(r)+4\,\pi\,r^{3}P_{\parallel}(r))}{r^{2}\,\left(1-\frac{2\,G\,M(r)}{r}\right)}=\frac{2}{r}\frac{P_{\perp}(r)-P_{\parallel}(r)}{\mathcal{E}(r)+P_{\parallel}(r)}-\frac{P_{\parallel}^{\prime}(r)}{\mathcal{E}(r)+P_{\parallel}(r)}$

where $f^{\prime}(r)\equiv\frac{d}{dr}f(r)$. In the Newtonian limit
($P_{\parallel},P_{\perp}\ll\mathcal{E}$, $G\,M\ll r$), $r$ is the radial
distance from the center, $M(r)$ is the mass within radius $r$ and $\Phi(r)$
is the gravitational potential. As in Newtonian gravity, the zero point of
$\Phi$ is unspecified and can be chosen for convenience. Regular stellar
solutions have $r=0$ at the center of the star, a finite $M_{\star}\equiv
M(r\rightarrow\infty)$ and require $2\,G\,M(r)<r$ so $M(0)=0$. The metric
asymptotes to the Schwarzschild form with mass $M_{\star}$ at large $r$.
$M(r)$ can be interpreted as the gravitational mass inside radius $r$ and
$M_{\star}$ the total gravitational mass of the star.

SM stars have very nearly isotropic pressure i.e.
$|P_{\parallel}-P_{\perp}|\ll
P\equiv\frac{1}{3}\,(P_{\parallel}+2\,P_{\perp})$ in which case the Tolman-
Oppenheimer-Volkoff equation (TOV) equation [3, 4] holds

$\displaystyle
M^{\prime}(r)=4\pi\,r^{2}\,\mathcal{E}(r)\qquad\frac{P^{\prime}(r)}{\mathcal{E}(r)+P(r)}=-\frac{G\,(M(r)+4\,\pi\,r^{3}\,P(r))}{r^{2}\,\left(1-\frac{2\,G\,M(r)}{r}\right)}.$
(3.4)

Cold SM stars have no atmospheres and the surface of the star at radius
$r=R_{\star}$ is defined by $P(R_{\star})=0$. The EoS defines a one parameter
class of stellar solutions. BSM stars may differ from SM stars not only in the
EoS but they may have anisotropic pressure, $P_{\parallel}\neq P_{\perp}$.

###  3.1 Balance of Chemical and Gravitational Forces

For any perfect fluid EoS, $P(\mathcal{E})$, one can define a conserved
“charge” density:
$n\propto\mathrm{exp}\int\frac{d\mathcal{E}}{\mathcal{E}+P(\mathcal{E})}$
which defines as conserved current $(n\,u^{\mu})_{;\mu}=0$ where $u^{\mu}$ is
the 4-velocity of the fluid. This conservation law derives from perfect fluid
energy momentum conservation:
$\left((\mathcal{E}+P)\,u^{\mu}u^{\nu}-P\,g^{\mu\nu}\right)_{;\nu}=0$, which
also gives $\frac{P^{\prime}(r)}{\mathcal{E}(r)+P(r)}=-\Phi^{\prime}(r)$ of
eq. 3.3. If one re-expresses $\mathcal{E}$ and $P$ in terms of $n$, i.e.
$\mathcal{E}(n)$ and $P(n)$, then one can define the chemical potential of
this charge by

$\displaystyle\mu(n)\equiv\mathcal{E}^{\prime}(n)=\frac{\mathcal{E}(n)+P(n)}{n}$
(3.5)

where the last equality, the adiabaticity equation, follows from the
definition of $n$. One sees that $\mu^{\prime}(n)=\frac{P^{\prime}(n)}{n}$ and
$\frac{P^{\prime}(n)}{\mathcal{E}(n)+P(n)}=\frac{\mu^{\prime}(n)}{\mu(n)}$.
Thus
$\Phi^{\prime}(r)=-\frac{P^{\prime}(r)}{\mathcal{E}(r)+P(r)}=-\frac{\mu^{\prime}(r)}{\mu(r)}$
which has solution

$\displaystyle\mu(r)\,e^{\Phi(r)}=\mathrm{constant}\ .$ (3.6)

For hydrostatic equilibrium a spatial decrease in the gravitational potential,
$\Phi$, must be compensated by a spatial increase in chemical potential,
$\mu$. This may be thought of as a balance of chemical and gravitational
forces on the conserved charges. This balance derives from energy momentum
conservation, is true of any hydrostatic configuration whether spherical or
not, and is more general than GR since it doesn’t derive from Einstein’s
equations. We will make use of this balance of forces below where the
conserved charge is baryon number.

The conserved quantity $n$ can be interpreted as the number density of some
“chemical” component of the fluid. The chemical balance of eq. 3.6 derives
from pressure balance of eq.s 3.3 & 3.4 so chemical and pressure balance
encode the same information and are not separate principles. This equivalence
does require adiabaticity, i.e. that there is no energy leakage as one
compresses fluid elements to higher densities, $\mathcal{E}$.

##  4 Ideal Gas + Yukawa (IGY) Equation of State

###  4.1 Free Fermions

In this section, we consider an ideal gas scenario where SM sector contains
free fermions only. For a free fermion species $f$ with (vacuum) mass $m_{f}$
and Yukawa coupling $g_{f}$ to $\phi$, the Dirac Lagrangian is given by

$\displaystyle\mathcal{L}_{\mathrm{eff}}\supset
i\,\bar{\psi}_{f}\,\gamma^{\mu}\,\partial_{\mu}\,\psi_{f}-(m_{f}-g_{f}\,\tilde{\phi})\,\bar{\psi}_{f}\,\psi_{f}\
.$ (4.1)

One sees that in a region with constant $\tilde{\phi}\neq 0$ the fermion mass
term differs from its vacuum value, i.e. $m_{f}\rightarrow
m_{f}-g_{f}\,\tilde{\phi}$. The dispersion relation for fermion $f$ is
$\epsilon_{f}^{2}={\tilde{m}_{f}}^{2}+|\mathbf{p}_{f}|^{2}$ where
$\epsilon_{f}$ and $\mathbf{p}_{f}$ are the particle energy and momentum and

$\displaystyle\tilde{m}_{f}(\tilde{\phi})\equiv|m_{f}-g_{f}\,\tilde{\phi}|$
(4.2)

is the _kinematic mass_ which we have defined to be non-negative. Since the
sign of the mass term in the Dirac Lagrangian is largely conventional the
change of sign of $m_{f}-g_{f}\,\tilde{\phi}$ in eq. 4.1 for
$\tilde{\phi}>m_{f}/g_{f}$ is not worrying. The $\tilde{\phi}$ dependence of
the kinematic mass of the SM fermion has important effects on EoS of the SM
fermions.

As shown in appendix A for an isotropic particle distribution function,

$\displaystyle\langle\bar{\psi}_{f}\psi_{f}\rangle$ $\displaystyle=$
$\displaystyle
n_{f}\,\langle\gamma^{-1}\rangle_{f}=\frac{\mathcal{E}_{f}-3\,P_{f}}{\tilde{m}_{f}(\tilde{\phi})}$
$\displaystyle\mathcal{E}_{f}$ $\displaystyle=$
$\displaystyle\tilde{m}_{f}(\tilde{\phi})\,n_{f}\,\langle\gamma\rangle_{f}$
(4.3) $\displaystyle P_{f}$ $\displaystyle=$
$\displaystyle\tilde{m}_{f}(\tilde{\phi})\,n_{f}\,\langle\frac{\gamma^{2}-1}{3\,\gamma}\rangle_{f}$

where $n_{f}$ is the fermion’s number density, $\gamma$ is the Lorentz factor
and $\langle\cdots\rangle_{f}$ gives the particle weighted average over the
distribution function. $\mathcal{E}_{f}$ and $P_{f}$ are the density and
(isotropic) pressure of these particles including the
$g_{f}\,\langle\bar{\psi}_{f}\psi_{f}\rangle\,\tilde{\phi}$ interaction term.
$n_{f},\langle\gamma^{-1}\rangle_{f}$, $\mathcal{E}_{f}$ and $P_{f}$ are all
frame dependent but their combination giving
$\langle\bar{\psi}_{f}\psi_{f}\rangle$ is not. There is no singularity in
$\langle\bar{\psi}_{f}\psi_{f}\rangle$ when $\tilde{m}_{f}\rightarrow 0$ since
in the massless limit $P_{f}\rightarrow\frac{1}{3}\,\mathcal{E}_{f}$.

Comparing eq. 4.1 with eq. 2.3, the quantity $\tilde{n}$ given rise by the
free fermions reads

$\displaystyle\tilde{n}=\sum_{f}\,\frac{g_{f}}{\tilde{g}}\,n_{f}\langle\gamma^{-1}\rangle_{f}=\sum_{f}\,\frac{g_{f}}{\tilde{g}}\,\frac{\mathcal{E}_{f}-3\,P_{f}}{\tilde{m}_{f}(\tilde{\phi})}\
.$ (4.4)

In a GR context, $\tilde{n}$’s dependence on $\mathcal{E}_{f}-3\,P_{f}$ is
very convenient, as the gravitational field also couples to a fermion through
$\mathcal{E}_{f}$ and $P_{f}$ via their contribution to $T_{\mu\nu}$. Eq. 3.1
now becomes

$\displaystyle
T_{\mu\nu}=T^{\tilde{\phi}}_{\mu\nu}+\sum_{f}\left(\left(\mathcal{E}_{f}+P_{f}\right)\,u^{f}_{\mu}\,u^{f}_{\nu}-P_{f}\,\mathrm{g}_{\mu\nu}+\tilde{P}^{f}_{\mu\nu}\right),$
(4.5)

where $u^{f}_{\mu}$ is the center of momentum 4-velocity and
$\tilde{P}^{f}_{\mu\nu}$ gives the anisotropic pressure tensor which is only
non-zero if the momentum distribution is anisotropic in the $u^{f}_{\mu}$
frame. The $\frac{\delta\tilde{n}}{\delta\mathrm{g}^{\mu\nu}}$ and
$\tilde{g}\,\tilde{\phi}\,\tilde{n}\,\mathrm{g}_{\mu\nu}$ terms in eq. 3.1 is
included in eq. 4.5, when one computes $\mathcal{E}_{f}$ and $P_{f}$ with the
kinematic mass $\tilde{m}_{f}$ instead of $m_{f}$.

In static stars, $u^{f}_{\mu}$ will be the stellar rest frame and we expect
from SM interactions a nearly isotropic momentum distribution in this frame so
$\tilde{P}^{f}_{\mu\nu}\simeq 0$ and

$\displaystyle\mathcal{E}=\mathcal{E}_{\tilde{\phi}}+\sum_{f}\mathcal{E}_{f}\
\qquad P=P_{\tilde{\phi}}+\sum_{f}P_{f}\qquad
P_{\parallel}-P_{\perp}=P^{\tilde{\phi}}_{\parallel}-P^{\tilde{\phi}}_{\perp}\
.$ (4.6)

Only the scalar field contributes to the anisotropic pressure.

###  4.2 Scalar Fields

The $\phi$-sector Lagrangian including the interaction term is give by

$\displaystyle\mathcal{L}_{\mathrm{eff}}\supset\frac{\tilde{\phi}^{;\mu}\,\tilde{\phi}_{;\mu}}{2}-V(\tilde{\phi})+\tilde{g}\,\tilde{\phi}\,\tilde{n}\leavevmode\nobreak\
.$ (4.7)

The classical equation of motion (EoM) and stress-energy is given by

$\displaystyle{\tilde{\phi}^{;\mu}}\phantom{}_{\mu}=-V^{\prime}(\tilde{\phi})+\tilde{g}\,\tilde{n}\qquad
T^{\phi}_{\mu\nu}=\tilde{\phi}_{;\mu}\,\tilde{\phi}_{;\nu}+\mathrm{g}_{\mu\nu}\left(V(\tilde{\phi})-\frac{\tilde{\phi}^{;\xi}\,\tilde{\phi}_{;\xi}}{2}\right)\leavevmode\nobreak\
.$ (4.8)

For $\langle\phi\rangle=0$ to be a stable classical solution in vacuum
($\tilde{n}=0$) and contribute negligibly to the cosmological constant it is
required that $V^{\prime}(0)=0$, $V^{\prime\prime}(0)>0$ and $V(0)\simeq 0$.

In the background of the static spherical star metric of eq. 3.2,
$\tilde{\phi}$ will have the same symmetries as the spacetime, i.e. with no
temporal or azimuthal dependence, or $\tilde{\phi}(r)$. This gives stress-
energy contribution

$\displaystyle\mathcal{E}_{\tilde{\phi}}(r)$ $\displaystyle=$
$\displaystyle+\frac{1-\frac{2\,G\,M(r)}{r}}{2}\,\tilde{\phi}^{\prime}(r)^{2}+V(\tilde{\phi}(r))$
$\displaystyle P_{\tilde{\phi}}(r)$ $\displaystyle=$
$\displaystyle-\frac{1-\frac{2\,G\,M(r)}{r}}{6}\,\tilde{\phi}^{\prime}(r)^{2}-V(\tilde{\phi}(r))$
(4.9) $\displaystyle
P^{\tilde{\phi}}_{\parallel}(r)-P^{\tilde{\phi}}_{\perp}(r)$ $\displaystyle=$
$\displaystyle\left(1-\frac{2\,G\,M(r)}{r}\right)\,\tilde{\phi}^{\prime}(r)^{2}$

including the anisotropic pressure. The wave equation of eq. 4.8 becomes

$\displaystyle\tilde{\phi}^{\prime\prime}(r)=\frac{1}{2}\,\frac{\partial\,\mathrm{ln}\frac{e^{-2\Phi(r)}}{r^{4}\,\left(1-\frac{2\,G\,M(r)}{r}\right)}}{\partial
r}\,\tilde{\phi}^{\prime}(r)+\frac{V^{\prime}(\tilde{\phi}(r))-\tilde{g}\,\tilde{n}(r)}{1-\frac{2\,G\,M(r)}{r}}$
(4.10)

For fixed $M(r)$ and $\tilde{n}(r)$ eq. 4.10 is formally an inhomogeneous ODE
for $\tilde{\phi}(r)$. This is only formal since $\tilde{\phi}(r)$ may also be
reflected in $M(r)$ and $\tilde{n}(r)$. If
$V(\phi)=\frac{1}{2}\,m_{\phi}^{2}\,\phi^{2}$ then
$\tilde{m}_{\phi}(\tilde{\phi})=m_{\phi}$ is constant and since
$V^{\prime}(\tilde{\phi}(r))=m_{\phi}^{2}\,\tilde{\phi}(r)$ this is a linear
equation in $\tilde{\phi}$. More generally eq. 4.10 is formally nonlinear in
$\tilde{\phi}$. For a given local value of $\tilde{\phi}$ we define the
effective local $\phi$ mass by $\tilde{m}_{\phi}(r)\equiv
V^{\prime\prime}(\tilde{\phi}(r))^{1/2}$.

As a convention we take $\phi=0$ to be the vacuum value of $\phi$. We also
impose various conditions on $V(\phi)$:

$\displaystyle V(0)\approx 0\qquad
V^{\prime\prime}(\tilde{\phi})>0\qquad\left|\frac{V^{\prime\prime}(\tilde{\phi})\,V(\tilde{\phi})}{V^{\prime}(\tilde{\phi})^{2}}\right|\mathrel{/\\!\\!\\!\\!\\!\\!\gg}1\qquad\left|\frac{V^{\prime\prime\prime}(\tilde{\phi})\,V^{\prime}(\tilde{\phi})}{V^{\prime\prime}(\tilde{\phi})^{2}}\right|\mathrel{/\\!\\!\\!\\!\\!\\!\gg}1\leavevmode\nobreak\
.$ (4.11)

which need only apply to the values of $\tilde{\phi}$ attained inside the star
under study. The 1st condition is assured by noting the observed cosmological
constant is negligible compared to the density of stars. The 2nd condition
avoids instabilities for large $\tilde{\phi}$ field values. The other
assumptions can be characterized as requiring that $V(\phi)$ not be
“unnaturally flat” and are required from some of the conclusions we obtain
below.

###  4.3 Polarization

For the simplest potential, $V(\phi)=\frac{1}{2}\,m_{\phi}^{2}\,\phi^{2}$, in
Minkowski space the classical $\phi$ EoM is the inhomogeneous Klein-Gordon
equation

$\displaystyle\ddot{\phi}-\nabla^{2}\phi+m_{\phi}^{2}\,\phi=\tilde{g}\,\tilde{n}$
(4.12)

whose classical and quantum phenomenology are well known. For a fixed
(classical) static $\tilde{n}({\bf x})$, the lowest energy (ground state)
solution is a Yukawa potential convolved with $\tilde{n}$:

$\displaystyle\phi({\bf x},t)=\tilde{\phi}({\bf
x})\equiv\frac{\tilde{g}}{4\pi}\int d^{3}{\bf
x^{\prime}}\,\frac{e^{-m_{\phi}\,|{\bf x}-{\bf x^{\prime}}|}}{|{\bf x}-{\bf
x^{\prime}}|}\,\tilde{n}({\bf x^{\prime}})\ .$ (4.13)

We call $\tilde{\phi}$ the “polarization” of $\phi$ by $\tilde{n}$. To this
one can add excitations $\propto e^{i\,({\bf k\cdot x}\pm\sqrt{|{\bf
k}|^{2}+m_{\phi}^{2}}\,t)}$ which can be quantized into $\phi$ particles. Even
with no such particle excitations there is a non-zero polarization:
$\tilde{\phi}\neq 0$. Quantum and classical fluctuations in $\tilde{n}$ should
be small when the number density of particles contributing to the integral is
large, roughly $\tilde{n}\gg m_{\phi}^{3}$ for non- or mildly relativistic
fermions.

A polarization, $\tilde{\phi}$, tied to a static $\tilde{n}$ will exist in any
static curved spacetime though the Green functions are more complicated. If
$\tilde{n}$ varies with time then $\tilde{\phi}$ will readjust itself to the
new ground state on a timescale $\sim m_{\phi}^{-1}$. In a star if the
timescale of change is much slower than $m_{\phi}^{-1}$, then any excess
energy above the ground state will go into coherent $\phi$ and SM matter
oscillations. Such excess energy above a ground state would eventually be
radiated away by SM cooling processes until a static configuration is reached.
If the temporal change has spectral content with angular frequency
$|\omega|>m_{\phi}$ then coherent $\phi$ radiation may be emitted. Since there
is no conserved charge associated with $\tilde{n}$ this can include $\phi$
monopole radiation.

###  4.4 Heavy, Light and Ultralight Scalars

In flat space with $V(\phi)=\frac{1}{2}\,m_{\phi}^{2}\,\phi^{2}$ we see that
in eq. 4.13 if $\tilde{n}$ does not vary significantly over the size of the
kernel, $m_{\phi}^{-1}$, then one can take $\tilde{n}$ out of the integral
obtaining the approximation

$\displaystyle\tilde{\phi}\approx\frac{\tilde{g}\,\tilde{n}}{m_{\phi}^{2}}\ .$
(4.14)

This we refer to as the light scalar approximation or LSA. The LSA corresponds
to the case where the spatial gradients of $\tilde{\phi}$ contribute
negligibly to the scalar wave equation, eq. 4.10, i.e.

$\displaystyle V^{\prime}(\tilde{\phi})$
$\displaystyle\approx\tilde{g}\,\tilde{n}$ (4.15)

which defines the LSA more generally when
$V(\phi)\neq\frac{1}{2}\,m_{\phi}^{2}\,\phi^{2}$. One might need to solve eq.
4.15 for $\tilde{\phi}$ numerically. Spatial gradients add to the energy
density and the ground state polarized field will minimize the gradients as
much as possible. The gradients become negligible when both $\tilde{n}$ and
$\tilde{m}_{\phi}$ not vary significantly over the range of the $\phi$-field,
$\tilde{m}_{\phi}(\tilde{\phi})^{-1}=V_{\phi}(\tilde{\phi})^{-1/2}$. In
particular $|{\bf\nabla}\,\mathrm{ln}\,\tilde{n}|\ll\tilde{m}_{\phi}$ and
$|{\bf\nabla}\,\mathrm{ln}\,\tilde{m}_{\phi}|\ll\tilde{m}_{\phi}$. Using
$\frac{V^{\prime\prime}(\tilde{\phi})}{V^{\prime}(\tilde{\phi})}\nabla\tilde{\phi}\approx\nabla\ln\tilde{n}$
derived from eq. 4.15 these conditions become

$\displaystyle|{\bf\nabla}\,\mathrm{ln}\,\tilde{n}|\ll
V^{\prime\prime}(\tilde{\phi})^{1/2}\qquad\left|\frac{V^{\prime\prime\prime}(\tilde{\phi})\,V^{\prime}(\tilde{\phi})}{V^{\prime\prime}(\tilde{\phi})^{5/2}}\,{\bf\nabla}\,\mathrm{ln}\,\tilde{n}\right|\ll
1\ .$ (4.16)

Given the 1st condition the 2nd will be satisfied so long as the naturalness
assumptions of eq. 4.11 is satisfied valid. Eq. 4.14 also applies to curved
spacetime if the spacetime curvature scale is much smaller than
${\tilde{m}_{\phi}}^{-1}$ (roughly $G\,\mathcal{E}\ll\tilde{m}_{\phi}^{2}$).

No matter the form of $V(\phi)$ the LSA is only valid where the mean field
approximation is valid so the SM matter distribution must be well approximated
as a continuous field, or roughly $\tilde{m}^{3}_{\phi}(\tilde{\phi}(R))\ll
n(R)$, where $n$ is the number density of the relevant fermions. The
constraint gives the upper mass range of the LSA. Scalars more massive than
this we call heavy scalars. Incoherent interactions of $\phi$ with SM
particles will dominate for heavy scalars in contrast to lighter scalars where
the interactions are coherent. We do not consider heavy scalars in this paper.

Ultra-light scalars have masses smaller than the LSA range. In this case the
relation of $\tilde{\phi}(R)$ to $\tilde{n}(R)$ is less trivial and generally
must be determined numerically. If $m(0)\,R_{\star}\lesssim 1$ then the
$\tilde{\phi}\neq 0$ will extend well beyond the star contributing to a finite
range 5th force which may effect orbiting objects. Light scalars and heavy
scalars in contrast to ultra-light scalars have their $\tilde{\phi}$ fields
confined to the star and do not manifest long range forces.

For example in neutron stars the density varies on a scale of $\sim 1\,$km
with nuclear densities in the center and the density of normal liquids or
solids ($n\sim\AA^{-3}$) at the surface. Thus for the LSA to be valid
throughout the star it is required that

$\displaystyle 10^{-10}\,\mathrm{eV}\ll m_{\phi}\ll
10^{3}\,\mathrm{eV}\leavevmode\nobreak\ .$ (4.17)

White dwarfs has similar surface densities but $\sim 10^{3}\times$ larger in
size so the LSA is valid for

$\displaystyle 10^{-13}\,\mathrm{eV}\ll m_{\phi}\ll
10^{3}\,\mathrm{eV}\leavevmode\nobreak\ .$ (4.18)

The partition of scalar masses into heavy, light and ultra-light regimes is
context dependent. The LSA may be valid in the core of a star but not near its
surface.

###  4.5 Light Scalar Equation of State

The $\phi$ stress-energy given by eq. 4.2 depends on two terms
$|\nabla\tilde{\phi}|^{2}$ and $V(\tilde{\phi})$. In the LSA, using eq. 4.15,
the ratio of these two terms is

$\displaystyle\frac{|\nabla\tilde{\phi}|^{2}}{|V(\tilde{\phi})|}=\left|\frac{V^{\prime}(\tilde{\phi})^{2}}{V^{\prime\prime}(\tilde{\phi})\,V(\tilde{\phi})}\,\frac{|\nabla\,\mathrm{ln}\,\tilde{n}|^{2}}{V^{\prime\prime}(\tilde{\phi})}\right|\ll
1\leavevmode\nobreak\ ,$ (4.19)

where we have used the LSA requirement
$|\nabla\,\mathrm{ln}\,\tilde{n}|^{2}\ll|V^{\prime\prime}(\tilde{\phi})|$ and
a naturalness assumption of eq. 4.11. Thus the LSA implies

$\displaystyle\mathcal{E}_{\tilde{\phi}}\approx-P_{\tilde{\phi}}\approx
V(\tilde{\phi})\qquad
P^{\tilde{\phi}}_{\parallel}-P^{\tilde{\phi}}_{\perp}\approx 0$ (4.20)

so even the scalar stress energy is nearly isotropic and one may use the TOV
equation 3.4 to compute stellar structure.

####  4.5.1 IGY Equation of State

To understand the qualitative aspects of how a linearly coupled $\phi$ effects
the EoS of matter consider a cold ideal gas of a single fermion, $f$, Yukawa
coupled to a scalar with the simplest scalar potential:
$V(\phi)=\frac{1}{2}\,m_{\phi}^{2}\,\phi^{2}$. This EoS we denote IGY for
ideal gas + Yukawa. Take $\tilde{g}=g_{f}$ so $\tilde{n}=\tilde{n}_{f}$. Given
eq.s 4.1, 4.4 & 4.14 we find for an arbitrary distribution function

$\displaystyle\mathcal{E}$ $\displaystyle=+V(\tilde{\phi})+\mathcal{E}_{f}=$
$\displaystyle+\frac{1}{2}\,\left(\frac{g_{f}\,n_{f}\,\langle\gamma^{-1}\rangle_{f}}{m_{\phi}}\right)^{2}+\left(m_{f}-\frac{g_{f}^{2}\,n_{f}}{m_{\phi}^{2}}\,\langle\gamma^{-1}\rangle_{f}\right)\,n_{f}\,\langle\gamma\rangle_{f}$
$\displaystyle P$ $\displaystyle=-V(\tilde{\phi})+P_{f}=$
$\displaystyle-\frac{1}{2}\,\left(\frac{g_{f}\,n_{f}\,\langle\gamma^{-1}\rangle_{f}}{m_{\phi}}\right)^{2}+\left(m_{f}-\frac{g_{f}^{2}\,n_{f}}{m_{\phi}^{2}}\,\langle\gamma^{-1}\rangle_{f}\right)\,n_{f}\,\langle\frac{\gamma^{2}-1}{3\,\gamma}\rangle_{f}\
.$ (4.21)

Specializing to an ideal cold degenerate Fermi gas using the dimensionless
functions, $\xi$ and $\varphi$, defined in eq.s C.8 & C.9 we obtain

$\displaystyle\xi(x)=\frac{x^{2}}{3\pi^{2}}\,\langle\gamma\rangle_{f}\qquad\varphi(x)=\frac{x^{2}}{3\pi^{2}}\,\langle\frac{\gamma^{2}-1}{3\,\gamma}\rangle_{f}\qquad\eta(x)=\frac{x^{2}}{3\pi^{2}}\langle\gamma^{-1}\rangle_{f}=\xi(x)-3\,\varphi(x)\
.$ (4.22)

in terms of
$x\equiv\beta_{\mathrm{F}}\gamma_{\mathrm{F}}=p_{\mathrm{F}}/\tilde{m}_{f}$
and $n_{f}=\frac{x^{2}}{3\pi^{2}}\,{\tilde{m}_{f}}^{3}$. Here $p_{\mathrm{F}}$
is the Fermi momentum. It is convenient to parameterize the EoS in terms of
$x$ rather than $n_{f}$.

Combining eq.s 4.2, 4.14 &4.22 gives an expression for the kinematic mass in
terms $x$:

$\displaystyle
r_{\mathrm{m}}(x)\equiv\frac{\tilde{m}_{f}}{m_{f}}=1-\mathfrak{g}_{f}^{2}\,r_{\mathrm{m}}(x)^{3}\,\eta(x)\qquad\mathfrak{g}_{f}\equiv\frac{g_{f}\,m_{f}}{m_{\phi}}\
.$ (4.23)

There is a single positive real solution for $r_{\mathrm{m}}$ to this cubic
equation

$\displaystyle
r_{\mathrm{m}}(x)\equiv\frac{\tilde{m}_{f}}{m_{f}}=\frac{\left(\sqrt{1+\frac{4}{27\,\mathfrak{g}_{f}^{2}\,\eta(x)}}+1\right)^{\frac{1}{3}}-\left(\sqrt{1+\frac{4}{27\,\mathfrak{g}_{f}^{2}\,\eta(x)}}-1\right)^{\frac{1}{3}}}{(2\,\mathfrak{g}_{f}^{2}\,\eta(x))^{\frac{1}{3}}}\in[0,1)\leavevmode\nobreak\
.$ (4.24)

The range of $\tilde{\phi}$ only spans
$0\leq\tilde{\phi}<\tilde{\phi}_{\mathrm{max}}^{f}\equiv m_{f}/g_{f}$ so
$r_{\mathrm{m}}(x)\geq 0$. In this model $g_{f}$ only enters the EoS through
$\mathfrak{g}_{f}$, the effective coupling.

Combining eq.s 4.5.1, 4.22 & 4.24 the IGY EoS is given by

$\displaystyle\frac{\mathcal{E}_{\mathrm{IGY}}(x)}{m_{f}^{4}}$
$\displaystyle=r_{\mathrm{m}}(x)^{4}\,\xi(x)+\mathcal{V}_{\mathrm{IGY}}(x)$
$\displaystyle\frac{n_{f}(x)}{m_{f}^{3}}$
$\displaystyle=r_{\mathrm{m}}(x)^{3}\,\frac{x^{3}}{3\,\pi^{2}}$
$\displaystyle\mathcal{V}_{\mathrm{IGY}}(x)$
$\displaystyle=\frac{{\mathfrak{g}_{f}}^{2}}{2}\,r_{\mathrm{m}}(x)^{6}\,\eta(x)^{2}$
$\displaystyle\frac{P_{\mathrm{IGY}}(x)}{m_{f}^{4}}$
$\displaystyle=r_{\mathrm{m}}(x)^{4}\,\varphi(x)-\mathcal{V}_{\mathrm{IGY}}(x)$
$\displaystyle\frac{\tilde{n}_{f}(x)}{m_{f}^{3}}$
$\displaystyle=r_{\mathrm{m}}(x)^{3}\,\eta(x)$
$\displaystyle\frac{E_{\mathrm{F}}(x)}{m_{f}}$
$\displaystyle=r_{\mathrm{m}}(x)\,\sqrt{1+x^{2}}\ .$ (4.25)

This is simplified further in appendix D. Here $\mathcal{V}_{\mathrm{IGY}}$
contains the $\tilde{\phi}$ contribution to the energy density $\mathcal{E}$
and pressure $P$. EoSs are usually expressed as $P(\mathcal{E})$ but there is
no analytic expression for either of these in IGY and the IGY EoS is only
defined implicitly by eq.s 4.5.1. Also given is the fermion Fermi energy,
$E_{\mathrm{F}}\equiv\sqrt{{\tilde{m}_{f}}^{2}+p_{\mathrm{F}}^{2}}$. One can
verify adiabaticity

$\displaystyle\mu_{f}(x)\equiv\frac{\partial\mathcal{E}_{\mathrm{IGY}}}{\partial
n_{f}}=\frac{{\mathcal{E}_{\mathrm{IGY}}}^{\prime}(x)}{{\tilde{n}_{f}}^{\prime}(x)}=\frac{\mathcal{E}_{\mathrm{IGY}}(x)+P_{\mathrm{IGY}}(x)}{\tilde{n}_{f}(x)}=E_{\mathrm{F}}(x)$
(4.26)

which is the 1st law of thermodynamics (energy conservation)
$dU=-P\,dV+T\,dS+\mu_{f}\,dN_{f}$ for adiabatic compression ($dS=0$) when the
number of fermions is conserved ($dN_{f}=0)$. Here $U$, $V$ and $N_{f}$ are
the energy, volume and number of fermions for a fluid element and $\mu_{f}$ is
the fermion chemical potential. For IGY $\mu_{f}=E_{\mathrm{F}}$ since
fermions are added to or subtracted from the Fermi surface.

The $\phi$ contribution to the the EoS in eq.s 4.20 & 4.5.1 has the form of a
cosmological constant which goes to zero at low densities, $n_{f}\rightarrow
0$, and saturates,
$\mathcal{V}_{\mathrm{IGY}}\rightarrow\frac{1}{2\,\mathfrak{g}_{f}^{2}}$, at
high densities, $n_{f}\rightarrow\infty$. Since the fermion density and
pressure grow unbounded as $n_{f}\rightarrow\infty$ the $\phi$ contribution
becomes negligible at both high and low densities. It is only at intermediate
densities, as the fermion transitions from non-relativistic to ultra-
relativistic behavior, that the Yukawa coupling will effect the EoS. The large
decrease in the fermion kinematic mass at very high density has no effect
since the fermion mass is insignificant for ultra-relativistic particles. The
fermions do become relativistic at lower densities the larger the value of
$\mathfrak{g}_{f}$ is due to the more rapid decrease in $\tilde{m}_{f}$ with
density. This combined with the negative scalar field pressure which comes to
dominate at intermediate densities results in a large $\mathfrak{g}_{f}$
dependence of the EoS for the range of intermediate densities. The asymptotic
behavior of the IGY EoS is independent of $\mathfrak{g}_{f}$:

$\displaystyle\mathcal{E}_{\mathrm{IGY}}\rightarrow\begin{cases}m_{f}\,n_{f}&\tilde{n}_{f}\rightarrow
0\\\
\left(\frac{9\pi}{8}\right)^{\frac{2}{3}}\,{n_{f}}^{\frac{4}{3}}&\tilde{n}_{f}\rightarrow\infty\end{cases}\qquad
P_{\mathrm{IGY}}\rightarrow\begin{cases}\frac{1}{5}\,(3\,\pi^{2})^{\frac{2}{3}}\,\frac{{n_{f}}^{\frac{5}{3}}}{m_{f}}&\tilde{n}_{f}\rightarrow
0\\\
\frac{1}{3}\,\left(\frac{9\pi}{8}\right)^{\frac{2}{3}}\,{n_{f}}^{\frac{4}{3}}&\tilde{n}_{f}\rightarrow\infty\end{cases}\
,$ (4.27)

the same as a degenerate Fermi gas with zero Yukawa coupling. At high density
the fermion kinematic mass becomes small but does not go negative since
$\tilde{\phi}>\frac{m_{f}}{g_{f}}$.

Figure 1: Plotted are the ideal gas + Yukawa (IGY) equation-of-state (EoS),
$P(\mathcal{E})$, from eq. 4.5.1 for various values of the coupling
$\mathfrak{g}_{f}\equiv\frac{g_{f}\,m_{f}}{m_{\phi}}$ as indicated in the
legend. Negative pressures are plotted as dashed curves. The partially
obscured black curve is the idea gas EoS with no Yukawa coupling.

The IGY EoS, $P(\mathcal{E})$, is plotted for various values of
$\mathfrak{g}_{f}$ in fig. 1. A significant suppression of pressure relative
to the $\mathfrak{g}_{f}=0$ EoS occurs for $\mathfrak{g}_{f}\gtrsim 3$ and
only for a range of $\mathcal{E}$ as illustrated. Pressure suppression is more
pronounced and grows more rapidly for larger $\mathfrak{g}_{f}$ leading to
larger suppression at lower densities and negative pressures for
$\mathfrak{g}_{f}\gtrsim 5$ . Pressure suppression is the main consequence of
a Yukawa coupling and is caused by the combination of decreasing
$\tilde{m}_{f}$ and the negative contribution of the scalar field polarization
$\tilde{\phi}$.

##  5 Equation of State Anomalies

The relativistic fluid described by the IGY EoS of eq. 4.5.1 has many
conventional properties

  * •

energy density increases with particle density:
$\frac{\partial\mathcal{E}}{\partial n_{f}}>0$,

  * •

no anomalously large sound speeds: $\frac{\partial
P}{\partial\mathcal{E}}<\frac{1}{3}$,

  * •

weak energy condition: $\mathcal{E}\geq 0$,

  * •

null energy condition: $\mathcal{E}+P\geq 0$,

  * •

dominant energy condition: $\mathcal{E}\geq|P|$ and

  * •

adiabaticity: $\frac{\partial\mathcal{E}}{\partial
n_{f}}=\frac{\mathcal{E}+P}{n_{f}}$.

Furthermore as discussed previously the small and large density behavior are
precisely those of an ideal cold Fermi gas of a single fermion species as
indicated above.

In addition to conventional behavior the IGY EoS exhibits anomalous behavior:

  1. 1.

$\frac{\partial P}{\partial\mathcal{E}}<0$ for
$\mathfrak{g}_{f}>\mathfrak{g}_{\mathrm{crit}}=4.987\ldots$ and
$\mathcal{E}^{-}_{\sim}(\mathfrak{g}_{f})\leq\mathcal{E}\leq\mathcal{E}^{+}_{\sim}(\mathfrak{g}_{f})$,

  2. 2.

$P<0$ for $\mathfrak{g}_{f}>\mathfrak{g}_{0}=5.124\ldots$ and
$\mathcal{E}^{-}_{0}(\mathfrak{g}_{f})\leq\mathcal{E}\leq\mathcal{E}^{-}_{0}(\mathfrak{g}_{f})$,

  3. 3.

$\mathcal{E}+3P<0$ for $\mathfrak{g}_{f}>28.094\ldots$ gravitational repulsion
- a violation of the strong energy condition.

each of which persist only over a finite range of densities for a given
$\mathfrak{g}_{f}$. While these anomalous states of matter do exist for IGY
(some hydrodynamically unstable for wavenumbers $k\lesssim\,m_{\phi}$) it is
argued below that they do not actually occur in a stellar environment. While
this pressure behavior is anomalous for stellar matter it is not uncommon in
condensed matter physics. As with condensed matter qualitatively new phenomena
can arise because of these anomalies and can be understood in terms of
different phases of matter as is now described.

###  5.1 Forbidden, Soft and Hard Matter

As one increases $\mathfrak{g}_{f}$ first anomaly 1 occurs first, then 1 & 2
and finally 1 & 2 & 3\. The onset of these pressure anomalies occur at
critical point: $\mathfrak{g}_{f}=\mathfrak{g}_{\mathrm{crit}}$ and
$\mathcal{E}_{\mathrm{crit}}=\mathcal{E}^{-}_{\sim}(\mathfrak{g}_{\mathrm{crit}})=\mathcal{E}^{+}_{\sim}(\mathfrak{g}_{\mathrm{crit}})=0.006197...$
when the fermions are non-relativistic and when the $\phi$ contributes only a
small fraction of the total: $\mathcal{E}_{\phi}/\mathcal{E}\approx 0.06$.
This is possible because the fermion pressure is so small,
$P_{f}/\mathcal{E}_{f}\approx 0.07$ whereas
$P_{\tilde{\phi}}/\mathcal{E}_{\tilde{\phi}}=-1$ in the LSA. The rapid onset
of anomalous behavior with increasing $\mathfrak{g}_{f}$ means that stellar
structure can be very sensitive to the value of this parameter.

For $\mathfrak{g}_{f}<\mathfrak{g}_{\mathrm{crit}}$ matter at all densities
behave normally (no anomalies) which we refer to as normal matter. For
$\mathfrak{g}_{f}>\mathfrak{g}_{\mathrm{crit}}$ some points in the
$\mathfrak{g}_{f}$-$\mathcal{E}$phase plane cannot exist inside a star. Such
phases of matter we refer to as forbidden and are indicated in fig. 2. In
order to be hydrodynamically stable one requires $\frac{\partial
P}{\partial\mathcal{E}}\geq 0$ which is not satisfied for
$\mathfrak{g}_{f}>\mathfrak{g}_{\mathrm{crit}}$ in the density range
$\mathcal{E}^{-}_{\sim}(\mathfrak{g}_{f})\leq\mathcal{E}\leq\mathcal{E}^{+}_{\sim}(\mathfrak{g}_{f})$.
This instability region is forbidden. Since a star has zero pressure at its
surface, as indicated by the TOV eq. 3.4, the gravitational overburden will
cause the pressure to rise as one goes from the surface toward the center and
it follows that $P\geq 0$ in the star. This is not satisfied for
$\mathfrak{g}_{f}>\mathfrak{g}_{0}$ in the density range
$\mathcal{E}^{-}_{0}(\mathfrak{g}_{f})\leq\mathcal{E}\leq\mathcal{E}^{-}_{0}(\mathfrak{g}_{f})$.
The negative pressure region is also forbidden. Since
$\mathcal{E}^{-}_{\sim}(\mathfrak{g}_{f})<\mathcal{E}^{-}_{0}(\mathfrak{g}_{f})<\mathcal{E}^{+}_{\sim}(\mathfrak{g}_{f})<\mathcal{E}^{+}_{0}(\mathfrak{g}_{f})$
the instability region and the negative pressure region overlap forming the
entire forbidden region. These two regions in the
$\mathfrak{g}_{f}$-$\mathcal{E}$ plane are shown in fig. 2.

For $\mathfrak{g}_{f}>\mathfrak{g}_{\mathrm{crit}}$ the forbidden region
separates a low density phase of matter that we call soft from a higher
density phase of matter that we call hard. This terminology derives from the
quantity $\frac{\partial\,\mathrm{ln}P}{\partial\,\mathrm{ln}\mathcal{E}}$
which is a measure of hardness of an EoS. The hard phase can have large values
$\frac{\partial\,\mathrm{ln}P}{\partial\,\mathrm{ln}\mathcal{E}}$ which
becomes infinite at $\mathcal{E}=\mathcal{E}_{0}^{+}$ whereas soft matter has
$0\leq\frac{\partial\,\mathrm{ln}P}{\partial\,\mathrm{ln}\mathcal{E}}\leq\frac{5}{3}$.

Figure 2:  Shown is a “phase diagram” for the IGY equation of state,
$P(\mathcal{E})$, as a function of the energy density $\mathcal{E}$ in units
of ${m_{f}}^{4}$ and the effective coupling
$\mathfrak{g}_{f}\equiv\frac{g_{f}\,m_{f}}{m_{\phi}}$. Since IGY is restricted
to cold ($T=0$) matter the phase is a function of one parameter,
$\mathcal{E}$, but this dependence differs for different couplings
$\mathfrak{g}_{f}$. Shown are normal, forbidden, soft and hard phases as
described in the text. The forbidden phase is not to be found in stars either
because of hydrodynamic instability, $\frac{\partial
P}{\partial\mathcal{E}}<0$ (green shaded) or negative pressure (gray shaded).
The onset of instability occurs at
$(\mathcal{E}_{\mathrm{crit}},\mathfrak{g}_{\mathrm{crit}})$ and negative
pressure at $(\mathcal{E}_{0},\mathfrak{g}_{0})$; these two anomalies are
restricted to densities
$\mathcal{E}_{\sim}^{-}<\mathcal{E}<\mathcal{E}_{\sim}^{+}$ and
$\mathcal{E}_{0}^{-}<\mathcal{E}<\mathcal{E}_{0}^{+}$, respectively; with
lower allowed densities referred to as soft and higher allowed densities as
hard. If
$\mathfrak{g}_{\mathrm{crit}}<\mathfrak{g}_{f}<\mathfrak{g}_{\mathrm{vac}}$
hard and soft matter can form a stable boundary in pressure and chemical
equilibrium at soft density $\mathcal{E}^{-}_{\mathrm{pb}}$ and hard density
$\mathcal{E}^{+}_{\mathrm{pb}}$ (pink dashed curve). If
$\mathfrak{g}_{f}\geq\mathfrak{g}_{\mathrm{vac}}$ hard matter with density
$\mathcal{E}_{0}^{+}$ (orange dashed curve) can form a stable boundary in
pressure equilibrium with vacuum. Soft matter in the shaded yellow region is
metastable to a transition to hard matter. If
$\mathfrak{g}_{f}<\mathfrak{g}_{\mathrm{crit}}$ there is no anomalous behavior
which we refer to as normal.

If $\mathfrak{g}_{f}\leq\mathfrak{g}_{\mathrm{crit}}$ then only a normal phase
exists, all $\mathcal{E}$ are possible since $\frac{\partial
P}{\partial\mathcal{E}}\geq 0$ and $P>0$. In this $\mathfrak{g}_{f}$-range a
star can contain any $\mathcal{E}$ allowed by gravitational stability. This is
not true if $\mathfrak{g}_{f}>\mathfrak{g}_{\mathrm{crit}}$ where a star can
consist of either a soft phase or a hard phase or both but there is a
forbidden range of $\mathcal{E}$. Soft, forbidden and hard phases occur
whenever $P(\mathcal{E})$ increases from zero at small $\mathcal{E}$,
decreases for intermediate $\mathcal{E}$ and increases again for large
$\mathcal{E}$. Below we describe how a hard and soft phase can coexist in
hydrostatic equilibrium and how a hard phase can interface with vacuum in
hydrostatic equilibrium with no intermediate soft phase.

The phases of Yukawa coupled cold matter is somewhat analogous to the behavior
of liquid and gas phases for ordinary materials such as water: the hard phase
is analogous to a liquid phase and the soft phase analogous to a gas phase.
These two phases only exist in the coupling range
$\mathfrak{g}_{f}>\mathfrak{g}_{\mathrm{crit}}$. One can consider a hard phase
as having certain (equilibrium) “vapor pressure” for the soft phase with which
it is in equilibrium. For a small range of couplings
$\mathfrak{g}_{\mathrm{crit}}<\mathfrak{g}_{f}<\mathfrak{g}_{\mathrm{vac}}$
this vapor pressure is non-zero allowing hard and soft phases to coexist
($\mathfrak{g}_{\mathrm{vac}}$ is the value of $\mathfrak{g}_{\mathrm{n}}$
where $\mu_{f}(\mathcal{E}_{0}^{+})=\epsilon_{f}(\mathcal{E}_{0}^{+})=m_{f})$.
It is possible that the soft phase can exist in a metastable state with
pressure larger than the equilibrium vapor pressure. This is analogous to
supersaturation of water vapor ($>100\%$ humidity). The equilibrium vapor
pressure is 0 if $\mathfrak{g}_{f}>\mathfrak{g}_{\mathrm{vac}}$ and all low
density soft matter is metastable in this case. There is a large barrier for
the metastable soft phase to transition to lower energy stable hard phase. One
way to seed such a phase transition is to compress the soft phase enough to
reach the instability region $\mathcal{E}_{f}>\mathcal{E}^{-}_{\sim}$ over a
volume greater than $m_{\phi}^{-3}$. The instability will grow to produce some
hard matter which can seed the transition of all the connected metastable soft
matter to hard matter.

###  5.2 Phase Boundaries

Phase boundaries allow matter of very different densities to exist alongside
each other, allowing stars or other objects in hydrostatic equilibrium to
avoid the forbidden unstable region. One can understand the phase boundaries
on macroscopic scales in terms of jump conditions between two allowed phases
arrived at from general principles or one can derive these jump conditions
from a microscopic analysis. We begin with a microscopic description on the
length scales of the force carrier, ${m_{\phi}}^{-1}$ abandoning the LSA which
is only valid on scales $\gg{m_{\phi}}^{-1}$.

A phase boundary is a hydrostatic configuration with planar symmetry
separating matter with two asymptotic bulk properties. The transition layer is
assumed to be narrow enough so that gravitational tidal forces from one side
to the other are negligible. In flat space stress-energy conservation,
${T^{\mu\nu}}_{;\nu}=0$, for hydrostatic solutions have constant pressure,
${\bf\nabla}P=0$. In terms of the scalar and fermionic pressure from eq.s 4.1,
4.8 & 4.5.1 the pressure profile is given by

$\displaystyle
P=P_{f}(z)+P^{\parallel}_{\phi}(z)=(m_{f}-g_{f}\,\tilde{\phi}(z))^{4}\,\varphi(x(z))+\frac{1}{2}\,\left(\tilde{\phi}^{\prime}(z)^{2}-{m_{\phi}}^{2}\,\tilde{\phi}(z)^{2}\right)$
(5.1)

where $z$ is the Euclidean coordinate orthogonal to the transition plane and
${}^{\prime}\equiv\frac{d}{dz}$. Using eq. 4.12 the equation of hydrostatic
equilibrium becomes
$P_{f}^{\prime}(z)=\mathfrak{g}_{f}\,\tilde{n}_{f}(z)\,\tilde{\phi}^{\prime}(z)$
which balances the scalar field gradient force by a pressure gradient in the
fermion fluid. From this one can show (for an ideal cold fermion gas) that the
the Fermi energy,
$E_{\mathrm{F}}=(m_{f}-g_{f}\,\tilde{\phi}(z))\,\sqrt{1+x(z)^{2}}$, is
constant, relating $x(z)$ to $\tilde{\phi}(z)$. Choosing the $+z$ to be the
direction in which $\tilde{\phi}$ increases one obtains the ODE describing
planar hydrostatic equilibrium without gravity:

$\displaystyle\tilde{\phi}^{\prime}(z)=\sqrt{2\,P-2\,(m_{f}-g_{f}\,\tilde{\phi}(z))^{4}\,\varphi\left[\sqrt{\left(\frac{E_{\mathrm{F}}}{m_{f}-g_{f}\,\tilde{\phi}(z)}\right)^{2}-1}\right]+{m_{\phi}}^{2}\,\tilde{\phi}(z)^{2}}$
(5.2)

whose solution is characterized by the two physical constants $P$ and
$E_{\mathrm{F}}$ plus a third which translates the configuration in $z$. Only
when $P$, $E_{\mathrm{F}}$ and $\tilde{\phi}$ take on values of any IGY
solution (eq.s 4.5.1 & 4.14) will $\tilde{\phi}^{\prime}(z)=0$ so only IGY
EoS’s can describe any asymptotic “bulk” properties of cold matter.

Figure 3:  This graph illustrates the phase boundaries which are allowed by
the IGY EoS. Superposed on an expanded view of the phase diagram of fig. 2 are
iso-pressure ($P$ \- black) and iso-chemical potential ($\mu_{f}$ \- dashed
purple) curves at values of these quantities which form stable interfaces
between soft and hard matter at four values of the effective coupling
$\mathfrak{g}_{f}=g_{f}m_{f}/m_{\phi}$. The values of $\mathfrak{g}_{f}$, $P$
and $\mu_{f}$ are given in the inset. The condition for a phase boundary is
that the pressure and chemical potential do not change across the boundary for
a fixed value of $\mathfrak{g}_{f}$. One can see that this condition is
satisfied for these values. These phase boundaries allow the matter in a star
to increase from low densities along the blue horizontal lines and then jump
from the pink dot on the left to the pink dot on the right at much higher
density from which the density can increase further. This jump is from soft
matter to hard matter circumventing the forbidden region. Soft/hard phase
boundaries do not exist when $\mathfrak{g}_{f}>\mathfrak{g}_{\mathrm{vac}}$
but in this case hard matter can form a stable interface with vacuum. This
vacuum/hard matter interface is indicated by the orange squares where the
density jumps from zero to a larger (hard) value and further increases inside
the star along the magenta lines. The pink dots lie along the pink dashed
$\mathcal{E}_{\mathrm{pb}}^{\pm}$ curves in fig. 2 while the orange squares
lie along the orange dashed $\mathcal{E}_{0}^{+}$ curve.

####  5.2.1 Connecting Hard and Soft Matter

A trivial solution to eq. 5.2 is any uniform IGY EoS where $P$,
$E_{\mathrm{F}}$ and $\tilde{\phi}$ are given by the same IGY phase but this
is not a boundary between two different phases. Values of $P$ and
$E_{\mathrm{F}}$ corresponding to phase boundaries must have different
asymptotic bulk properties corresponding to two different phases, i.e. two
different $\tilde{\phi}$ values with the same $P$ and $E_{\mathrm{F}}$. For an
adiabatic EoS the chemical potential have extrema when
$0=\frac{\partial\mu_{f}}{\partial n_{f}}=\frac{1}{n_{f}}\,\frac{\partial
P}{\partial n_{f}}$ or $P^{\prime}(\mathcal{E})=0$ (since
$\frac{\partial\mathcal{E}}{\partial n_{f}}>0$). For normal matter
($\mathfrak{g}_{f}<\mathfrak{g}_{\mathrm{crit}}$) $P^{\prime}(\mathcal{E})>0$
and thus $\frac{\partial\mu_{f}}{\partial n_{f}}>0$ so no two phases have the
same $P$ or the same $\mu_{f}$ and a phase boundary is not possible. However
when $\mathfrak{g}_{f}>\mathfrak{g}_{\mathrm{crit}}$ both $P$ and $\mu_{f}$
have a maximum and a minimum at $\mathcal{E}=\mathcal{E}_{\sim}^{-}$ and
$\mathcal{E}_{\sim}^{+}$, respectively so different phases can have the same
$P$ and/or the same $E_{\mathrm{F}}$ and phase boundaries are possible.
Mathematically, from eq. 4.5.1, phase boundaries are described by solutions of

$\displaystyle-\frac{1}{2}\,\mathfrak{g}_{f}^{2}\,r_{\mathrm{m}}(x_{-})^{6}\,\eta(x_{-})^{2}+r_{\mathrm{m}}(x_{-})^{4}\,\varphi(x_{-})$
$\displaystyle=$
$\displaystyle-\frac{1}{2}\,\mathfrak{g}_{f}^{2}\,r_{\mathrm{m}}(x_{+})^{6}\,\eta(x_{+})^{2}+r_{\mathrm{m}}(x_{+})^{4}\,\varphi(x_{+})$
$\displaystyle r_{\mathrm{m}}(x_{-})\,\sqrt{1+{x_{-}}^{2}}$ $\displaystyle=$
$\displaystyle r_{\mathrm{m}}(x_{+})\,\sqrt{1+{x_{+}}^{2}}$ (5.3)

where $x_{\pm}\equiv x(\pm\infty)$ and $x_{-}<x_{+}$. The asymptotic energy
densities of these solutions are

$\displaystyle\mathcal{E}^{\pm}_{\mathrm{bp}}(\mathfrak{g}_{f})\equiv\mathcal{E}(\pm\infty)\
=m_{f}^{4}\,\left(+\frac{1}{2}\,\mathfrak{g}_{f}^{2}\,r_{\mathrm{m}}(x_{\pm})^{6}\,\eta(x_{\pm})^{2}+r_{\mathrm{m}}(x_{\pm})^{4}\,\xi(x_{\pm})\right)$
(5.4)

which can differ by orders of magnitude. $\mathcal{E}^{\pm}_{\mathrm{bp}}$ lie
along the pink dashed curves in fig. 2 which meet at the critical point. There
is a single solution of eq.s 5.3 only for
$\mathfrak{g}_{\mathrm{crit}}<\mathfrak{g}_{f}\leq\mathfrak{g}_{\mathrm{vac}}$.
Fig. 3 shows four examples of such solutions for different values of
$\mathfrak{g}_{f}$. The width of the transition layer is several
${m_{\phi}}^{-1}$ which by assumption is much smaller than the length scale of
the matter configuration under the LSA. The asymptotic $x_{-}$ IGY phase is
always soft matter and the asymptotic $x_{+}$ is always hard matter so these
phase boundaries allow the fluid to jump over the forbidden region. The phase
boundaries occur at particular values of soft and hard matter densities which
are $\mathfrak{g}_{f}$ dependent. The pressure at these phase boundaries are
always positive.

One could have foregone this microscopic analysis since from general
principles a stable phase boundary should obey the Gibbs jump conditions 1)
pressure balance: the pressure of the two phases must be equal, for
hydrostatic equilibrium and 2) chemical equilibrium: the chemical potential of
the fermions (the Fermi energy in this case) must be equal otherwise fermions
would flow into the phase with lower $\mu_{f}$ causing the boundary to move.
This macroscopic analysis at the level of LSA only requires solutions of the
jump conditions of eq. 5.3 and would treat the phase boundary as being
discontinuous.

Figure 4: Plotted is the profile of a 1-D phase boundary between soft with
hard matter for $\mathfrak{g}_{f}=5.2$ as described in §5.2.1. A bulk IGY soft
state is reached as $z\rightarrow-\infty$ and bulk IGY hard matter state as
$z\rightarrow+\infty$. The red, blue and brown curves give the profiles of the
fermion energy density, the scalar field energy density and the scalar field;
all as a function of the distance, $z$, from a nominal center of the phase
boundary. The two horizontal lines give the asymptotic soft and hard matter
total energy density.

####  5.2.2 Connecting Hard Matter with Vacuum

The largest coupling, $\mathfrak{g}_{f}$, for a soft/hard phase boundary is
$\mathfrak{g}_{\mathrm{vac}}$ at which one finds $P=0$,
$E_{\mathrm{F}}=m_{f}$,
$\mathcal{E}^{-}_{\mathrm{bp}}(\mathfrak{g}_{\mathrm{vac}})=0$ and
$\mathcal{E}^{+}_{\mathrm{bp}}(\mathfrak{g}_{\mathrm{vac}})=\mathcal{E}^{+}_{0}(\mathfrak{g}_{\mathrm{vac}})$;
so this is essentially a phase boundary of hard matter with vacuum.
Hard/vacuum interfaces do not exist for
$\mathfrak{g}_{f}<\mathfrak{g}_{\mathrm{vac}}$ because $E_{\mathrm{F}}>m_{f}$
and hence fermions would flow into the vacuum. Hard/vacuum interfaces do exist
for $\mathfrak{g}_{f}>\mathfrak{g}_{\mathrm{vac}}$ because
$E_{\mathrm{F}}<m_{f}$ and while fermions would flow from the vacuum there are
no fermions in vacuum to flow, i.e. one does not have to have chemical
equilibrium with vacuum. Pressure equilibrium is still required which for
vacuum interface is $P=0$. For IGY these two conditions are only satisfied on
the curve $\mathcal{E}=\mathcal{E}^{+}_{0}(\mathfrak{g}_{\mathrm{vac}})$ for
$\mathfrak{g}_{f}>\mathfrak{g}_{\mathrm{vac}}$ which can form a stable
interface with vacuum such as on the surface of a star. These vacuum interface
states are shown as the orange dashed curve in fig. 2 and orange squares in
fig. 3.

A microscopic description of the hard/vacuum interface is a 1D hydrostatic
configuration described by eq. 5.2 with $P=0$ and $E_{\mathrm{F}}\leq m_{f}$.
In this case the fermion energy density goes to zero (the argument of
$\varphi$ goes to zero) when $\tilde{\phi}=(m_{f}-E_{\mathrm{F}})/g_{f}$. One
can define $z=0$ at this position. Retaining the $\tilde{\phi}^{\prime}(z)>0$
convention: for $z>0$ there are fermions while for $z<0$ there are no fermions
but $\tilde{\phi}\neq 0$. For $z<0$ eq. 5.2 reduces to
$\tilde{\phi}^{\prime}(z)=m_{\phi}\,\tilde{\phi}(z)$ and matching the $z<0$
and $z>0$ solutions for $\tilde{\phi}$ one finds
$\tilde{\phi}(z)=\frac{m_{f}-E_{\mathrm{F}}}{g_{f}}\,e^{m_{\phi}\,z}$ for
$z\leq 0$. In order for the $z>0$ solution to asymptote to a uniform bulk
state, $\tilde{\phi}^{\prime}(+\infty)=0$, one requires that the Fermi energy
matches that of the IGY hard matter zero pressure state which is
$E_{\mathrm{F}}=\mathcal{E}_{0}^{+}/n_{0}^{+}$. Thus such an interface has an
exponential “$\phi$ atmosphere”,

$\displaystyle\tilde{\phi}(z<0)=\frac{1}{g_{f}}\,\left(m_{f}-\frac{\mathcal{E}_{0}^{+}(\mathfrak{g}_{f})}{n_{0}^{+}(\mathfrak{g}_{f})}\right)\,e^{m_{\phi}\,z}\qquad
n_{f}(z<0)=0$ (5.5)

where $n_{0}^{+}$ is the fermion number density when
$\mathcal{E}=\mathcal{E}_{0}^{+}$. There is no simple expression for the $z>0$
solution but an example numerical solution is shown in fig. 5. Again we we see
that the transition layer is only a few ${m_{\phi}}^{-1}$ thick and in the LSA
this phase boundary can be treated as discontinuous.

Figure 5: Plotted is the profile of a 1-D phase boundary of vacuum with hard
matter for $\mathfrak{g}_{f}=6$ as described in §5.2.2. The blue, red and
brown curves give the profiles of the fermion energy density, the scalar field
energy density and the scalar polarization; all as a function of the distance,
$z$, from the fermion surface. For $z\leq 0$ the fermion energy density,
$\mathcal{E}_{f}$, is zero while the scalar field falls off $\propto
e^{m_{\phi}\,z}$ and the scalar field energy density falls off $\propto
e^{2\,m_{\phi}\,z}$. Pure vacuum is reached as $z\rightarrow-\infty$ and bulk
hard matter IGY state is reached as $z\rightarrow+\infty$. The horizontal line
gives the total energy density of the asymptotic bulk hard matter.

####  5.2.3 Two Phase Instability

These two types of phase boundaries clarify what happens to hydrodynamically
unstable cold matter which occurs when
$\mathcal{E}_{\sim}^{-}<\mathcal{E}<\mathcal{E}_{\sim}^{+}$ for
$\mathfrak{g}_{f}>\mathfrak{g}_{\mathrm{crit}}$. Consider a thought experiment
where one sets up cold Yukawa coupled fermions confined to a finite volume in
a uniform unstable state. Since it is unstable density inhomogeneities will
grow releasing free energy. If one then allows these fermions to cool back to
zero temperature, assuming the volume is $\gg{m_{\phi}}^{-3}$ but not large
enough for gravity to play a role, the matter will reach an inhomogeneous two
phase equilibrium separated by the phase boundaries described above. If
$\mathfrak{g}_{\mathrm{crit}}<\mathfrak{g}_{f}<\mathfrak{g}_{\mathrm{vac}}$
then the two phases are a soft phase with energy density
$\mathcal{E}_{\mathrm{pb}}^{-}$ and a hard matter phase with much higher
energy density $\mathcal{E}_{\mathrm{pb}}^{+}$. This is analogous to a liquid
with non-zero equilibrium vapor pressure. If
$\mathfrak{g}_{f}>\mathfrak{g}_{\mathrm{vac}}$ the two phases will be vacuum
and hard matter with energy density $\mathcal{E}_{0}^{+}$. This is analogous
to a liquid with zero equilibrium vapor pressure. One could imagine such a two
phase instability to occur during the collapse of a cold star. In this paper
we have focused on cold matter but have shown that the Yukawa force decreases
as the energy of the particles increases. Therefore one could also imagine
that the onset of this instability could occur as a star cools.

####  5.2.4 Metastable Matter

States with the lowest energy per conserved fermion are energetically
preferred. This is determined by the chemical energy:
$\epsilon_{f}\equiv\frac{\mathcal{E}}{n_{f}}$ which does not include
gravitational binding energy. For an adiabatic EoS the chemical energy extrema
are given by $0=\frac{\partial\epsilon_{f}}{\partial n_{f}}=P/n_{f}^{2}$, i.e.
at $P=0$ points. For IGY these densities are $\mathcal{E}=0$ and for
$\mathfrak{g}_{f}\geq\mathfrak{g}_{0}$ also
$\mathcal{E}=\mathcal{E}_{0}^{\pm}$. The chemical energy is minimized,
maximized, and minimized at $\mathcal{E}=0,\,\mathcal{E}_{0}^{-}$ and
$\mathcal{E}_{0}^{+}$, respectively. At these extrema the chemical energy is
equal to the chemical potential: $\epsilon_{f}=\mu_{f}$. The low density
minimum is $\epsilon_{f}=m_{f}$ while the high density minimum is
$\epsilon_{f}=\mu_{f}=\mathcal{E}_{0}^{+}/n_{0}^{+}$. As indicated above
$\mathcal{E}_{0}^{+}<m_{f}\,n_{0}^{+}$ for
$\mathfrak{g}_{f}\geq\mathfrak{g}_{\mathrm{vac}}$ in which case
$\mathcal{E}_{0}^{+}$ is energetically preferred over a low density state.
Since for soft matter $\epsilon_{f}$ increases with density
$\mathcal{E}_{0}^{+}$ is energetically preferred to all soft matter for
$\mathfrak{g}_{f}\geq\mathfrak{g}_{\mathrm{vac}}$. This energetic preference
means that energy would be released by transforming soft matter into hard
matter. However the soft matter does not spontaneously transform to this lower
energy state since it is hydrodynamically stable. We say that soft matter is
metastable if $\mathfrak{g}_{f}\geq\mathfrak{g}_{\mathrm{vac}}$. Metastability
is indicated by yellow in fig. 2. Note that this implies that SM matter would
be metastable at low densities if it is coupled to a scalar field with a
sufficiently strong Yukawa coupling!

Metastability also occurs when
$\mathfrak{g}_{\mathrm{crit}}<\mathfrak{g}_{f}<\mathfrak{g}_{\mathrm{vac}}$
where soft/hard phase boundaries are possible. At the phase boundary
$\mathcal{E}$ jumps from $\mathcal{E}_{\mathrm{bp}}^{-}$ to
$\mathcal{E}_{\mathrm{bp}}^{+}$, jumping over soft, forbidden and hard matter.
All of the denser soft matter in the range
$\mathcal{E}_{\mathrm{pb}}^{-}<\mathcal{E}<\mathcal{E}_{\sim}^{-}$ have larger
chemical potentials than those of hard matter with
$\mathcal{E}\geq\mathcal{E}_{\mathrm{pb}}^{+}$ so the fermions would flow from
the soft to hard matter if the former came into contact with the latter.
However this denser soft matter is hydrodynamically stable and would not
spontaneously transition into hard matter. In this sense the denser soft
matter is also metastable as indicated by yellow in fig. 2.

The features of the IGY EoS above the critical point: soft/forbidden/hard
matter, phase boundaries, two phase instabilities and metastable soft phases
are generic features of EoSs with an anomalous $P(\mathcal{E})$ which rises,
falls and then rises again.

##  6 IGY Prototype Yukawa Coupled Neutron Star

Here we compute the structure of a spherical hydrostatic object consisting of
matter described by the IGY EoS. The IGY EoS has relevance to the structure of
neutron stars where density and pressure are both dominated by a single
fermion species, the neutron, at least over some range of densities. For this
reason we take $f\rightarrow\mathrm{n}$ in this section and quote physical
values for quantities based on the neutron mass. IGY treats baryonic matter as
an ideal cold gas of neutron just as was done in the foundational work of
Oppenheimer and Volkoff [4] (OV). The OV model corresponds to an IGY neutron
star with $\mathfrak{g}_{\mathrm{n}}=0$. IGY and OV have numerous shortcomings
in describing neutron stars. At low densities they do not incorporate the
pressure of muons and electrons or the changing abundance of protons and
nuclei with pressure. At higher densities they do not include hadronic
interactions which stiffen the EoS; and will eventually lead to the
disassociation of neutrons into quarks. Thus IGY like OV provides only a rough
approximation to a neutron star structure over a limited range of densities.
The purpose of this exercise is to illustrate the phenomenology of stellar
structure which can be caused by a Yukawa coupling. The expectation is that
qualitative features of IGY prototype neutron stars will also manifest in
stars with general and more realistic EoSs with Yukawa coupling to a light
scalar.

The different phases of IGY matter identified above and the ability to form
stable phase boundaries have a profound effect on stellar structure. This
depends on the strength of the Yukawa coupling. One may define three ranges of
$\mathfrak{g}_{f}$ values according to the allowed phases of matter which
hydrostatic configurations may contain:

  * •

weak coupling: $\mathfrak{g}_{f}\leq\mathfrak{g}_{\mathrm{crit}}=4.987\ldots$

    * -

normal stars: all normal matter

  * •

moderate coupling:
$\mathfrak{g}_{\mathrm{crit}}<\mathfrak{g}_{f}<\mathfrak{g}_{\mathrm{vac}}$

    * -

soft stars: all $\mathcal{E}<\mathcal{E}_{\sim}^{-}$ soft matter

    * -

hybrid stars: outer $0<\mathcal{E}<\mathcal{E}_{\mathrm{pb}}^{-}$ soft shell
with inner $\mathcal{E}>\mathcal{E}_{\mathrm{pb}}^{+}$ hard core

  * •

strong coupling: $\mathfrak{g}_{f}\geq\mathfrak{g}_{\mathrm{vac}}=5.681\ldots$

    * -

soft stars: all $\mathcal{E}<\mathcal{E}_{\sim}^{-}$ soft matter

    * -

hard stars: all $\mathcal{E}>\mathcal{E}_{0}^{+}$ hard matter

The densities referred to here are in the macroscopic light scalar
approximation, LSA, not the microscopic description illustrated in fig.s 4 &
5. IGY hybrid stars are analogous to hybrid quark stars [5] which have have
been extensively studied.

Since the IGY has isotropic pressure, $P(\mathcal{E})$, the spherical stellar
structure, $\mathcal{E}(r)$, is given by solutions to eq. 3.4 with boundary
conditions $\mathcal{E}^{\prime}(0)=0$ and $P(\mathcal{E}(R_{\star}))=0$. Zero
surface pressure requires that $\mathcal{E}(R_{\star})=0$ for all solutions
except hard stars which have $\mathcal{E}(R_{\star})=\mathcal{E}_{0}^{+}$.
Solutions can be parameterized by the effective Yukawa coupling,
$\mathfrak{g}_{\mathrm{n}}$, and the energy density at the center,
$\mathcal{E}(0)$.

We have determined the stellar structure, $\mathcal{E}(r)$, $P(r)$,
$M_{\star}$, $R_{\star}$, etc. for thousands of IGY stellar models
parameterized by $\mathfrak{g}_{\mathrm{n}}$ and $\mathcal{E}(0)$. Below we
will only present a sample of these results which we believe is large enough
for the reader to get a sense of the phenomenology to be able to interpolate
in this parameter space. We also briefly explain why the solutions behave as
they do. We begin with some gross features of the solutions and then go into
more details.

Not all hydrostatic solutions will be stable which further limits viable
hydrostatic models of stars. The hydrodynamical instability identified in 5.1
is for a uniform medium neglecting gravity. Stars are not uniform and may have
other hydrodynamic or gravitational instabilities. Below we treat such global
instabilities in only a perfunctory way by determining its nominal stability
according to the sign of $\frac{dM_{\star}}{d\mathcal{E}(0)}$ [6]. Negative
values are indicative of instability to gravitational core collapse. A
definitive determination of stability by normal mode analysis is considered
too cumbersome to apply to the large parameter space of stellar models under
consideration.

###  6.1 Mass and Radius

The most easily determined (measured) quantity for observed neutron stars is
their rotation rates. However this does not bear much relation to the star’s
EoS so it is not most relevant for this paper. A more relevant quantity is the
gravitational mass of neutron stars, $M_{\star}$, which has been determined
with some accuracy for a large number of stars. $M_{\star}$ combined with
direct measurement of neutron star radii, $R_{\star}$ (currently few in
number) will provide strong constraints on the neutron stars EoS.
Gravitational wave observations of merging neutron stars can also measure
combined tidal deformablilty which depends on the density profile within the
star with a strong dependence on $R_{\star}$. A complete
$M_{\star}(R_{\star})$ relation (accounting for rotation) would determine the
EoS $\mathcal{E}(P)$ if, as expected, the nuclear matter EoS is universal. A
theoretical $M_{\star}(R_{\star})$ curve (usually) gives a maximal mass,
$M_{\star}^{\mathrm{max}}$, and minimum radius, $R_{\star}^{\mathrm{min}}$,
which are generally achieved simultaneously in the maximal star. Prediction of
these quantities are of particular interest as it allows one to exclude
candidate EoSs using single star measurements of either $M_{\star}$ or
$R_{\star}$.

With this motivation we start the exposition of solutions to the TOV equation
for IGY stars by plotting $M_{\star}(R_{\star})$ for various values of
$\mathfrak{g}_{\mathrm{n}}$ while taking the fermion to be a neutron. For
clarity we plot $M_{\star}(R_{\star})$ for weak, moderate and strong coupling
separately; in fig.s 6, 7 & 8, respectively. The black
$\mathfrak{g}_{\mathrm{n}}=0$ curve in fig. 6 is the OV model which is the
“standard model” curve with which to compare the other curves. Irrespective of
$\mathfrak{g}_{\mathrm{n}}$ all IGY models at small masses only sample the
$\Gamma\equiv\frac{d\mathrm{ln}P}{d\mathrm{ln}\mathcal{E}}\approx\frac{5}{3}$
part of the EoS in fig. 1. Thus IGY $M_{\star}(R_{\star})$ asymptotes to a
$\Gamma=\frac{5}{3}$ polytrope (polytropic index $\frac{3}{2}$) for small
masses / low central densities where $M_{\star}\propto R_{\star}^{-3}$. These
are the stellar models corresponding to the points on the curves extending to
the lower right in fig.s 6, 7 & 8. It is $M_{\star}(R_{\star})$ where the EoS
deviates from this polytropic limit which depends on
$\mathfrak{g}_{\mathrm{n}}$.

Most EoSs, including IGY, asymptote to ultra-relativistic matter with
$\Gamma\rightarrow 1$, i.e. become softer, at large enough density. Such a
small $\Gamma$ at the center of a star indicates instability to core collapse.
Thus there is a maximum central density, $\mathcal{E}_{\mathrm{max}}$, for
stable stellar solutions and the star central density
$\mathcal{E}_{\mathrm{max}}$ yields the maximal star with maximum mass,
$M_{\star}^{\mathrm{max}}$. As can be seen from fig. 1 for
$\mathfrak{g}_{\mathrm{n}}\gtrsim 4$ the EoS can also become soft at
intermediate densities. This soft EoS at intermediate densities leads to
nominally unstable solutions (plotted in gray) spanning a range of $R_{\star}$
for $4\lesssim\mathfrak{g}_{\mathrm{n}}<\mathfrak{g}_{\mathrm{vac}}$ in fig.s
6 & 7.

Figure 6: $M_{\star}$ vs. $R_{\star}$ for spherical hydrostatic configurations
of matter with the IGY equation of state for the values of the effective
coupling,
$\mathfrak{g}_{\mathrm{n}}\equiv\frac{g_{\mathrm{n}}\,m_{\mathrm{n}}}{m_{\phi}}$,
as indicated in the legend. The couplings in this graph are limited to the
range of weak couplings
$0\leq\mathfrak{g}_{\mathrm{n}}\leq\mathfrak{g}_{\mathrm{crit}}$. $M_{\star}$
is the total gravitational mass and $R_{\star}$ the circumferential radius.
Solutions which are nominally gravitationally unstable are plotted in light
gray. The dot gives the critical star solution where where
$\mathfrak{g}_{\mathrm{n}}=\mathfrak{g}_{\mathrm{crit}}$ and the central
density is $\mathcal{E}(0)=\mathcal{E}_{\mathrm{crit}}$ where
$P^{\prime}(\mathcal{E})=P^{\prime\prime}(\mathcal{E})=0$. This unstable
solution is the separatrix between normal stars which are illustrated here and
the soft and hybrid stars which occur for larger coupling.

For weak couplings, $\mathfrak{g}_{\mathrm{n}}<\mathfrak{g}_{\mathrm{crit}}$,
since $\Gamma>0$ one can obtain the pressure gradient, $P^{\prime}(r)$,
required by eq. 3.4 to support the star with a sufficiently large density
gradient, $\mathcal{E}^{\prime}(r)$. For moderate couplings this is not always
possible since at intermediate densities $\Gamma<0$ for
$\mathfrak{g}_{\mathrm{n}}>\mathfrak{g}_{\mathrm{crit}}$. As one transitions
from the weak coupling stellar solutions in fig. 6 to the moderate coupling
solutions in fig. 7 these large density gradients are replaced by the
discontinuous soft/hard phase interface described in §5.2.1. For moderate
couplings there is not only an instability gap between large and small radius
solutions but there are two separate branches of solutions to the the TOV
equation. In fig. 7 these two branches bifurcate at a stellar solution
indicated by a dot which are stars with central density
$\mathcal{E}(0)=\mathcal{E}_{\mathrm{pb}}^{-}(\mathfrak{g}_{\mathrm{n}})$. The
larger $R_{\star}$ branch are soft stars consisting only of soft matter. The
soft matter branch extends beyond the bifurcation point to larger mass and
smaller radius. In this extension the center consists of chemically metastable
soft matter described in §5.2.4. The other branch extending from the
bifurcation point, tracing the outline of a swan, are hybrid stars consisting
of hard matter in the center and soft matter on the outside separated by a
soft/hard interface. This interface is at the center of the star at the tip of
the swan’s tail and moves outward along this curve but never reaches the
surface. The swan’s tail region is nominally gravitationally unstable but
becomes stable in the swan’s back at lower mass and radius proceeding up the
neck of the swan to the top of it’s head where the solution again becomes
gravitationally unstable at the maximal star.

Figure 7: $M_{\star}$ vs. $R_{\star}$ with the IGY equation of state in the
range of moderate couplings
$\mathfrak{g}_{\mathrm{crit}}\leq\mathfrak{g}_{\mathrm{n}}\leq\mathfrak{g}_{\mathrm{vac}}$.
There are two branches of the $M_{\star}(R_{\star})$ solutions which coincide
at bifurcation solutions denoted by the dots . The larger radius
$M_{\star}(R_{\star})$ curve consists only of stars made of soft matter which
we call soft stars. The other curve emanating from the dot (which resembles
the outline of a swan) are hybrid stars which consist of hard matter in the
center and soft matter on the outside joined together by the soft/hard
interface described in §5.2.1. The gray parts of the curve are solutions which
are nominally gravitationally unstable.

As $\mathfrak{g}_{\mathrm{n}}\rightarrow\mathfrak{g}_{\mathrm{vac}}$ from
below the bifurcation point moves to $M_{\star}\rightarrow 0$ and
$R_{\star}\rightarrow\infty$. The swan’s back and tail also move to
$M_{\star}\rightarrow 0$ and disappear leaving only the head and the neck. For
strong coupling, $\mathfrak{g}_{\mathrm{n}}>\mathfrak{g}_{\mathrm{vac}}$,
there is no stable hard/soft interface, so hybrid stars cannot exist, and
there is no bifurcation point. Instead the two solution branches separate into
two disconnected solutions as shown in fig. 8. A gap exists between these two
solutions but not because of gravitational instability; rather there are no
hydrostatic solutions, stable or unstable, in this gap. The large $R_{\star}$
solutions are soft stars consisting solely of soft matter. Soft star solutions
as seen for moderate couplings continue into the hard coupling regime. The
maximum soft star mass shrinks with $\mathfrak{g}_{\mathrm{n}}$, increasing
the gap between the two branches, since the energy density becomes softer at
lower density for larger $\mathfrak{g}_{\mathrm{n}}$ as can be seen from fig.
1. The small $R_{\star}$ solutions are hard stars consisting solely of hard
matter with a vacuum interface on the surface as described in §5.2.2. The hard
star $M_{\star}(R_{\star})$ curves are cane shaped: the swan’s neck morphed
into the length of the cane and the swan’s head into the cane’s handle. The
top of the handle is the maximal star while the length of the cane below the
handle describe constant density spheres with
$M_{\star}\propto{R_{\star}}^{3}$ extending to very low mass and radius. We
call solutions in the constant density regime hard drops and discuss them
further in §6.2. Hard stars fill the $R_{\star}$-$M_{\star}$ plane
approximately up to a maximum density giving $M_{\star}\lesssim
0.0064\,\left(\frac{R_{\star}}{\mathrm{km}}\right)^{3}\,M_{\odot}$ and a
maximum surface gravitational potential giving $M_{\star}\lesssim
0.2\,\frac{R_{\star}}{\mathrm{km}}\,M_{\odot}$ (see §6.2 & §6.4). For any
normal, soft or hybrid star there will be a hard star with larger
$\mathfrak{g}_{\mathrm{n}}$ but the same $M_{\star}$ and $R_{\star}$.

Figure 8: $M_{\star}$ vs. $R_{\star}$ with the IGY equation of state in the
range of strong couplings
$\mathfrak{g}_{\mathrm{n}}\geq\mathfrak{g}_{\mathrm{vac}}$. There are two
separated branches of the $M_{\star}(R_{\star})$ solutions. The curves on the
right soft star solutions consisting only of soft matter. The curves on the
left (which resemble canes) are hard stars which consist only of hard matter
and are joined on the surface to vacuum by the vacuum/hard interface described
in §5.2.2. The gray parts of the curve are solutions which are nominally
gravitationally unstable.

###  6.2 Hard Drops

Hard drops form the low mass end of the hard star solutions for strong
coupling, $\mathfrak{g}_{\mathrm{n}}\geq\mathfrak{g}_{\mathrm{vac}}$. Hard
drops are not primarily bound by gravity but rather by the attractive $\phi$
mediated force between neutrons. The structure of hard drops is simple: on the
surface there is a vacuum/hard interface as described in §5.2.2 while in the
interior there is uniform hard matter. Since gravity is unimportant the hard
matter is in hydrostatic equilibrium with zero pressure vacuum and has
nearly888Similarly to a water drop surface tension and finite size leads to an
interior pressure slightly in excess of the ambient pressure. zero pressure
itself, so the energy density takes the value $\mathcal{E}_{0}^{+}$. Thus each
of these solutions have roughly constant energy density throughout the star:
$\mathcal{E}(R)\approx\mathcal{E}_{0}^{+}$, so
$M_{\star}\approx\frac{4\pi}{3}\,\mathcal{E}_{0}^{+}\,{R_{\star}}^{3}$. Hard
drop sizes extend down to the range of the $\phi$-force,
$R_{\star}\gtrsim{m_{\phi}}^{-1}$, which is the range of validity of the LSA
(light scalar approximation). Depending on $m_{\phi}$ the minimum hard drop
size can be very much smaller than the $\sim 10\,$km size of neutrons stars.
To $<15\%$ accuracy we find for IGY hard drops the density and chemical energy
are $\mathcal{E}_{0}^{+}\approx
2\,\frac{\mathfrak{g}_{\mathrm{n}}-4}{{\mathfrak{g}_{\mathrm{n}}}^{3}}\,{m_{\mathrm{n}}}^{4}$
and
$\epsilon_{\mathrm{n}}\equiv\frac{\mathcal{E}_{0}^{+}}{n_{0}^{+}}\approx\sqrt{\frac{\mathfrak{g}_{\mathrm{vac}}}{\mathfrak{g}_{\mathrm{n}}}}\,m_{\mathrm{n}}$.
Thus the mass-radius relation is
$M_{\star}\approx\frac{8\pi\,(\mathfrak{g}_{\mathrm{n}}-4)}{3\,{\mathfrak{g}_{\mathrm{n}}}^{3}}\,{m_{\mathrm{n}}}^{4}\,{R_{\star}}^{3}$
up to the maximum density
$\frac{3\,M_{\star}}{{4\pi\,R_{\star}}^{3}}<\frac{\mathcal{E}_{\mathrm{vac}}}{n_{\mathrm{vac}}}\approx\frac{{m_{\mathrm{n}}}^{4}}{60}\approx
0.0015\,\,M_{\odot}/\mathrm{km}^{3}$ when
$\mathfrak{g}_{\mathrm{n}}=\mathfrak{g}_{\mathrm{vac}}$. The binding energy
per neutron is $m_{\mathrm{n}}-\epsilon_{\mathrm{n}}\approx
m_{\mathrm{n}}\,\left(1-\sqrt{\frac{\mathfrak{g}_{\mathrm{vac}}}{\mathfrak{g}_{\mathrm{n}}}}\right)$
which can be a large fraction $m_{\mathrm{n}}$.

The hard drop name derives from the similarity to drops of a normal liquid
which consist of atoms or molecules bound together by local forces in
hydrostatic equilibrium with the surrounding medium. Hard drops like normal
drops will have a surface tension which adjusts the surface to minimize its
curvature but as this is only a surface effect any configuration with the same
volume, whether in different shapes or in multiple drops, have nearly
identical energy. The gravitational force is long range and although it plays
a subdominant role in binding a hard drop it may nevertheless be the dominant
force in keeping the drop spherical and in one piece; or making it aspherical
in the presence of tidal forces.

Small hard drops are not sufficiently massive to be candidates for observed
neutron stars but one might imagine that the debris field produced during
neutron star formation or neutron star mergers might contain such drops. A
hard drop thus produced would form a seed for conversion to hard matter of any
low density matter it comes into contact, since for strong couplings the soft
matter is always metastable to conversion to hard matter. Since the binding
energy released per baryon is very large one would expect such a conversion to
produce observable phenomena. Accretion of the interstellar medium onto hard
drops would cause them to glow. Interaction with stars or planets could
release an enormous amount of energy, however since the radiation emitted
would be Eddington limited, $L_{\gamma}<3\times
10^{4}\,\frac{M}{M_{\odot}}\,L_{\odot}$, they would not be exceedingly
luminous. Super-Eddington luminosity would likely blow away much of the
stellar matter before it could convert.

The fact that the dark matter density is only a factor of $\sim 5\times$
larger than the baryon density is suggestive that these two cosmological
remnants are closely related. One should consider hard drops (possibly of
quarks rather than nucleons) as a candidate for dark matter since then the
dark and baryonic matter would have a common origin in a cosmological baryon
asymmetry. Hard drops would carry most of this asymmetry but in a sequestered
form. How hard drops could form in the early universe is a matter of
speculation. Hard drops and stars share a number of similarities to strange
quark nuggets [7] and stars [8], preon stars [9] and neutrino nuggets [1].

###  6.3 Radius Gap

For $\mathfrak{g}_{\mathrm{n}}\gtrsim 4$ IGY stars have a gap in the allowed
radii. In this radius range solutions to the TOV equations are either unstable
or do not exist. A gap can exist for weak coupling normal stars (fig. 6) and
always exists for moderate coupling hybrid stars (fig. 7). The gap also always
exists for strong coupling between soft and hard stars (fig. 8). As can be
seen from these figures radius gaps are accompanied by stable solutions with
double valued $R_{\star}(M_{\star})$. Fig. 9 shows the gap as a function of
$\mathfrak{g}_{\mathrm{n}}$. The radius gap appears at
$\mathfrak{g}_{\mathrm{n}}\gtrsim 4$ and spans nearly an order of magnitude in
radius for $\mathfrak{g}_{\mathrm{n}}>\mathfrak{g}_{\mathrm{vac}}$.

Figure 9:  Plotted as a function of the coupling
$\mathfrak{g}_{\mathrm{n}}\equiv\frac{g_{\mathrm{n}}\,m_{\mathrm{n}}}{m_{\phi}}$
is the range of stellar radii, $R_{\star}$, of nominally stable solutions to
the TOV equation given the IGY EOS. The IGY EoS is an idealization and the
excluded regions plotted here (in black) will vary for a more realistic EoS
though we expect the qualitative geometry of the excluded regions to be
insensitive to improvements in the EoS. The condition for (nominal) stability
used here is $\frac{dM_{\star}}{d\mathcal{E}(0)}>0$. The labels give the
different classes of stars as described in the text. The arrows point in the
“general” direction of increasing $M_{\star}$. This directionality cannot be
made precise because multiple solutions with the same $R_{\star}$ may have
different $M_{\star}$ and different signs of $\frac{dM_{\star}}{dR_{\star}}$.

###  6.4 Maximal Stars

The maximum measured neutron star mass has provided an important constraint on
nuclear EoSs. The theoretical maximal star, the spherical hydrostatic
configuration with the maximum mass, is a local maximum of
$M_{\star}(R_{\star})$, defining a boundary between the nominally stable and
unstable solutions. For IGY stars one can identify these points on the curves
in fig.s 6, 7, & 8: the maximal star is at the top of the swan’s head for
moderate coupling and the top of the cane’s handle for strong coupling. Where
there is a radius gap the maximal star lies below the gap. Denote the maximal
star’s mass and radius by $M_{\star}^{\mathrm{mx}}$ and
$R_{\star}^{\mathrm{mx}}$. Adjacent solutions with
$R_{\star}<R_{\star}^{\mathrm{mx}}$ are nominally unstable and solutions with
$R_{\star}>R_{\star}^{\mathrm{mx}}$ are nominally stable.

Fig. 10 plots $M_{\star}^{\mathrm{mx}}$ and $R_{\star}^{\mathrm{mx}}$ for
selected values of $\mathfrak{g}_{\mathrm{n}}$ from $0$ to $27$. Signification
deviations from the maximal OV star occurs only for
$\mathfrak{g}_{\mathrm{n}}>1$; for $\mathfrak{g}_{\mathrm{n}}\lesssim 5$ the
Yukawa coupling results in a decrease in $M_{\star}^{\mathrm{mx}}$ of up to
$\sim 40\%$ and a decrease in $R_{\star}^{\mathrm{mx}}$ by up to $\sim 70\%$.
This trend is rapidly reversed at $\mathfrak{g}_{\mathrm{n}}\sim 5$ where the
maximal star’s mass and radius begins to rise. For
$\mathfrak{g}_{\mathrm{n}}>6$ the maximal star has
$\frac{M_{\star}^{\mathrm{mx}}}{M_{\odot}}\approx\frac{R_{\star}^{\mathrm{mx}}}{5\,\mathrm{km}}\approx\frac{\mathfrak{g}_{\mathrm{n}}}{12}$
which can be compared with the OV maximal star with
$M_{\star}^{\mathrm{mx}}\approx 0.7\,M_{\odot}$ and
$R_{\star}^{\mathrm{mx}}\approx 9\,$km. Maximal hard stars exceeds the OV mass
for $\mathfrak{g}_{\mathrm{n}}\gtrsim 10$ while the maximal hard star radius
exceeds the OV value only for $\mathfrak{g}_{\mathrm{n}}\gtrsim 23$. Depending
on $\mathfrak{g}_{\mathrm{n}}$ maximal IGY stars may be more or less massive
than the maximal OV star and larger or smaller than the maximal OV star. No
maximal IGY star is both less massive and larger in size than the maximal OV
star. $\frac{G\,M_{\star}^{\mathrm{mx}}}{R_{\star}^{\mathrm{mx}}}\approx 0.26$
for maximal hard stars in comparison to 0.114 for the maximal OV star and
$\frac{4}{9}$ which is Buchdahl’s bound on this quantity [10].

Figure 10:  Plotted for various $\mathfrak{g}_{\mathrm{n}}$ is the radius and
mass of the maximal IGY star, the spherical hydrostatic solution configuration
with maximal mass. Each points is labelled by
$\mathfrak{g}_{\mathrm{n}}=\frac{g_{\mathrm{n}}\,m_{\mathrm{n}}}{m_{\phi}}$.
Note $\mathfrak{g}_{\mathrm{crit}}\approx 4.987$ and
$\mathfrak{g}_{\mathrm{vac}}\approx 5.681$. To emphasis the change in stellar
properties caused by non-zero Yukawa couplings the left and bottom axes give
mass and radius in relative to that of the $\mathfrak{g}_{\mathrm{n}}=0$ OV
star represented by the large gray dot. The right and top axes give mass in
radius in physical units using $m_{f}=m_{\mathrm{n}}$ appropriate for a
neutron star. None of these models account for nuclear interactions which
would increase both the mass and the radius of a more realistic neutron star.
The blue and black shaded regions are disallowed:
$\frac{G\,M_{\star}}{R_{\star}}\geq\frac{4}{9}$ and
$\frac{G\,M_{\star}}{R_{\star}}\geq\frac{1}{2}$. The former is an upper limit
derived for an infinitely stiff EoS [10] and the latter are black holes.

###  6.5 Density Profiles

To illustrate the different types of density profiles IGY predicts in fig. 11
is plotted the radial density profile for stars all of the same mass
$M_{\star}=0.3\,M_{\odot}$ but for different $\mathfrak{g}_{\mathrm{n}}$. For
this $M_{\star}$ there is only a single stable configuration for all
$\mathfrak{g}_{\mathrm{n}}$ which would not be true of smaller $M_{\star}$.
Shown are four normal stars with weak coupling, two hybrid stars with moderate
coupling and seven hard stars with strong coupling. The central density,
$\mathcal{E}(0)$, increases with $\mathfrak{g}_{\mathrm{n}}$ by $\sim
100\times$ up to moderate couplings of $\mathfrak{g}_{\mathrm{n}}\sim 5$, and
then decreases. The stellar radius, $R_{\star}$, decreases as $\mathcal{E}(0)$
increases as one would expect. In the IGY model $\mathcal{E}(0)\rightarrow 0$
and $R_{\star}\rightarrow\infty$ as
$\mathfrak{g}_{\mathrm{n}}\rightarrow\infty$. As one increases
$\mathfrak{g}_{\mathrm{n}}$ the profile approaches a step function:
$\mathcal{E}(r)=\mathcal{E}(0)\,\Theta(R_{\star}-r)$ as one obtains for hard
drops. For the hybrid stars illustrated the density decrease at the soft/hard
interface is $\sim 2\times$ and $\sim 100\times$ for
$\mathfrak{g}_{\mathrm{n}}=5$ and $\mathfrak{g}_{\mathrm{n}}=5.5$,
respectively. The hard matter core extends out to 83.7% and 95.6% of the
stellar radius and contains 98.4% and 99.97% of the stellar mass for these two
couplings, The soft matter envelope not only has small densities but steeply
falling densities. The density inside hard stars does not vary more than a
factor of 2 from it’s surface value
$\mathcal{E}(R_{\star})=\mathcal{E}_{0}^{+}$ and is very nearly constant for
$\mathfrak{g}_{\mathrm{n}}\gg\mathfrak{g}_{\mathrm{vac}}$ which is the hard
drop limit.

Figure 11:  The radial density profile, $\mathcal{E}(r)$, for IGY stars of
gravitational mass $M_{\star}=0.3\,M_{\odot}$ with values of the
$\mathfrak{g}_{\mathrm{n}}=\frac{g_{\mathrm{n}}\,m_{\mathrm{n}}}{m_{\phi}}$ as
indicated in the legend. For $\mathfrak{g}_{\mathrm{n}}=0,\,\ldots,\,4$ these
are normal stars. For $\mathfrak{g}_{\mathrm{n}}=5,\,5.5$ these are hybrid
stars with the density discontinuity inside the star illustrated by a gap in
the $\mathcal{E}(r)$ curve. For $\mathfrak{g}_{\mathrm{n}}=6,\,\ldots,\,30$
these are hard stars with a discontinuous density drop to zero (vacuum or
$\mathcal{E}=0$) at the end of the curves. Density discontinuity is a
macroscopic description; see §5.2 for microscope description. For
$\mathfrak{g}_{\mathrm{n}}>30$ the uniform density hard drop description of
§6.2 is accurate.

###  6.6 Binding Energy

Another important parameter of stellar remnants is the number of baryons they
contain which we denote by $N_{\star}$. Taking $m_{\mathrm{n}}$ as a
conventional “rest mass” per baryon we can define the “rest mass” of a star by
$m_{\mathrm{n}}\,N_{\star}$. This is plotted as a function of
$\mathfrak{g}_{\mathrm{n}}$ in fig. 12a for IGY neutron stars with
gravitational mass $0.3\,M_{\odot}$ and for maximal stars. One sees that the
baryon content of a fixed gravitational mass stellar remnant increases for
large $\mathfrak{g}_{\mathrm{n}}$ due to the increased binding energy due of
the attractive $\phi$ force. For strong coupling the maximal stars are hard
stars which are predominantly bound by $\phi$ forces and can attain both large
gravitational mass and large baryon content.

One can define the average binding energy per baryon by $b_{\star}\equiv
m_{\mathrm{n}}-M_{\star}/N_{\star}$ (using the conventional rest
mass,$m_{\mathrm{n}}$ ). This is plotted in fig. 12b for $0.3\,M_{\odot}$ and
maximal stars. For large $\mathfrak{g}_{\mathrm{n}}$ the $0.3\,M_{\odot}$ (or
any fixed $M_{\star}$) stars approaches the hard drop limit: $b_{\star}\approx
m_{\mathrm{n}}\,\left(1-\sqrt{\frac{\mathfrak{g}_{\mathrm{vac}}}{\mathfrak{g}_{\mathrm{n}}}}\right)$
where the gravitational binding energy becomes negligible (see §6.2). Maximal
stars have only a slightly larger $b_{\star}$ due to non-negligible
gravitational binding energy.

The large binding energy of hard stars is reflected in the energy liberated
during their formation. Denote the gravitational mass and baryon number of the
progenitor star(s) by $M_{0}$ and $N_{0}$, respectively, so the binding energy
per baryon of the progenitors is $b_{0}\equiv m_{\mathrm{n}}-M_{0}/N_{0}$. One
can decompose the total energy released when the progenitor transitions to the
remnant into the “rest mass” ejected, $M_{\mathrm{ej}}$, and the “radiated”
energy, $E_{\mathrm{rad}}$. These two components are given by999These
equations slightly overestimate the actual rest mass ejected and underestimate
the actual radiated energy because the ejecta is largely metals with rest mass
per baryon a few MeV less than $m_{\mathrm{n}}$.

$\displaystyle M_{\mathrm{ej}}$ $\displaystyle=$ $\displaystyle
m_{\mathrm{n}}\,(N_{0}-N_{\star})$ (6.1) $\displaystyle E_{\mathrm{rad}}$
$\displaystyle=$ $\displaystyle
M_{0}-M_{\star}-M_{\mathrm{ej}}=b_{\star}\,N_{\star}-b_{0}\,N_{0}=\frac{b_{\star}}{m_{\mathrm{n}}-b_{\star}}\,M_{\star}-\frac{b_{0}}{m_{\mathrm{n}}-b_{0}}\,M_{0}\
.$

Both $M_{\mathrm{ej}}$ and $E_{\mathrm{rad}}$ are potentially observable. The
former, say, by observations of the SNe produced, the latter by measurements
of neutrinos, gravitational radiation and photons emitted. Neutrino emission
constitutes the largest fraction of $E_{\rm rad}$. A small outflow kinetic
energy is also included in $E_{\rm rad}$ because much of it will be
transformed to heating gas which will cool radiatively by photons. In most
case the progenitor is a white dwarf or a nuclear burning star at the end of
its lifetime in which case $b_{0}\,N_{0}\ll b_{\star}\,N_{\star}$ so
$E_{\mathrm{rad}}\approx\frac{b_{\star}}{m_{\mathrm{n}}-b_{\star}}\,M_{\star}$.
This approximation for $E_{\mathrm{rad}}$ is plotted as a function of
$\mathfrak{g}_{\mathrm{n}}$ in fig. 12c.

From fig. 12c figure one sees that $E_{\rm rad}$ increases rapidly with
increasing $\mathfrak{g}_{\mathrm{n}}$. For $\mathfrak{g}_{\mathrm{n}}$
significantly greater than $\mathfrak{g}_{\mathrm{vac}}$ when the stellar
structure is close to that of hard drops so $E_{\rm
rad}\approx\left(\sqrt{\frac{\mathfrak{g}_{\mathrm{n}}}{\mathfrak{g}_{\mathrm{vac}}}}-1\right)\,M_{\star}$
and the maximum energy radiated by the maximal star is

$\displaystyle
E_{\mathrm{rad}}^{\mathrm{mx}}\approx\frac{\mathfrak{g}_{\mathrm{n}}}{12}\,\left(\sqrt{\frac{\mathfrak{g}_{\mathrm{n}}}{\mathfrak{g}_{\mathrm{vac}}}}-1\right)\,M_{\odot}=8.5\times
10^{53}\,\mathrm{erg}\,\times\,\frac{\mathfrak{g}_{\mathrm{n}}}{\mathfrak{g}_{\mathrm{vac}}}\,\left(\sqrt{\frac{\mathfrak{g}_{\mathrm{n}}}{\mathfrak{g}_{\mathrm{vac}}}}-1\right)$
(6.2)

for $\mathfrak{g}_{\mathrm{n}}>\mathfrak{g}_{\mathrm{vac}}$. The source of
this extra energy is the large binding energy associated with the scalar
mediated attraction between nucleons. For these larger couplings the radiated
energy greatly exceeds that produced by formation of a SM OV star. Thus we
find that Yukawa coupled neutron stars produce anomalously large amounts of
radiated energy, mostly in the form of a neutrino burst, during their
formation. Measuring neutrino bursts from a star collapsing into a neutron
star probes the scalar coupling of nucleons without direct observation of
either the progenitor or remnant star.

(a) Rest mass ($\equiv m_{\mathrm{n}}\,N_{\mathrm{b}}$). (b) Average binding
energy per baryon
($b_{\star}\equiv\frac{m_{\mathrm{n}}\,N_{\mathrm{b}}-M_{\star}}{N_{\mathrm{b}}}$).

(c) Energy released during formation ($E_{\mathrm{rad}}\approx
b_{\star}\,N_{\star}$).

Figure 12: Plotted as a function of
$g_{\mathrm{n}}\equiv\frac{g_{\mathrm{n}}\,m_{\mathrm{n}}}{m_{\phi}}$ for IGY
neutron stars are the three quantities stated in the sub-captions. The lower
curves are for stars with gravitational mass $M_{\star}=0.3\,M_{\odot}$ and
the upper curves are for stars with the maximum mass (see §6.4)

###  6.7 IGY Star Anomalies

In §6.1-§6.6 we describe properties of stars with the IGY EoS . These stellar
solutions have a number of qualitative and quantitative features which differ
from what one finds in most if not all SM EoSs. These are all related to the
EoS anomalies of §5 though some appear at even weaker coupling than the EoS
anomalies. All of these anomalies involve observable quantities which allow
one to limit or detect scalar coupling of nucleons by astronomical observation
of neutron stars.

Here we list several mass/radius anomalies:

  * •

mass versus radius inversion : Stable hydrostatic solutions with SM EoSs have
${M_{\star}}^{\prime}(R_{\star})<0$. This is not true for all stable IGY
stars: some hybrid stars with moderate coupling (fig. 7) and some hard stars
(all hard drops) with strong coupling (fig. 8) have
${M_{\star}}^{\prime}(R_{\star})>0$.

  * •

maximum mass / minimum radius dichotomy : For SM EoSs stable cold hydrostatic
configurations have a maximum mass, $M_{\star}^{\mathrm{mx}}$. Since
${M_{\star}}^{\prime}(R_{\star})<0$ this maximal star also has the minimum
radius. In contrast for IGY stars with couplings
$\mathfrak{g}_{\mathrm{n}}\gtrsim 5.4$ the minimum radius configuration is not
the maximal star (fig.s 7 & 8).

  * •

radius gap : For SM EoSs there is a minimum radius below which there are no
stable configurations and stable configurations extend from this minimal
radius up to an instability gap between neutron stars and white dwarfs. In
contrast for IGY stars there is a gap in radii in the neutron star range
($\sim 10$ km) for $\mathfrak{g}_{\mathrm{n}}\gtrsim 4$.

  * •

mass degeneracy : For SM EoSs there is single configuration with a given mass
and km scale radius. In contrast for IGY stars there are are two stable
solutions with different radii over some mass range if
$\mathfrak{g}_{\mathrm{n}}\gtrsim 4$.

  * •

small stars : SM EoSs have stable solutions with a minimum radius of several
km. In contrast hard stars with strong coupling (fig. 8) have stable solutions
(hard drops) extending down to “microscopic” radii:
$R_{\star}\sim{m_{\phi}}^{-1}$.

  * •

large stars : SM EoSs have stable solutions with a maximum mass of a few
$M_{\odot}$ with a radius of $\sim 10\,\mathrm{km}$. In contrast hard stars
with $\mathfrak{g}_{\mathrm{n}}>30$ have stable solutions with much larger
masses and these maximal stars have much larger radii
($M_{\star}\approx\frac{\mathfrak{g}_{\mathrm{n}}}{12}\,M_{\odot}$,
$R_{\star}\approx 0.4\,\mathfrak{g}_{\mathrm{n}}\,\mathrm{km}$, see §6.4 and
fig.s 8 & 10).

As more neutron star radius measurements become available each of these
anomalies will provide qualitative probe of nucleon Yukawa coupling. An
exception might be the small star / hard drops phenomena. Whether hard drops
would actually form in the course of normal stellar evolution is uncertain and
even if they did their small mass and radius would make them difficult to
detect. If no such anomalies are found only weak Yukawa coupling will be
allowed.

Other anomalies in addition to those related to mass and radius anomalies
include

  * •

thin crusts : Stable hydrostatic solutions with SM EoSs contain superfluid
neutrons with a crust of normal matter which is $\sim 1\,$km thick. While we
have not addressed the issue of superfluidity in IGY we note that some hybrid
stars (see §6.5 & fig. 11) and all hard stars (see §5.2.2 & fig. 5) are nearly
uniform up to an extremely rapidly falling density near the surface. Thus we
would expect any crust to be very thin in these stars. The outer region of a
hard star is “microscopic” with thickness $\sim{m_{\phi}}^{-1}$. Pulsar
glitches are believed to originate from rotational slippage between the
superfluid and normal matter in the crust. This crust must have sufficient
moment of inertia to produce observed glitches which is unlikely to be the
case in hard and some hybrid stars (see ref. [11] for this argument applied to
strange stars). Thus pulsar timing also places limits on scalar coupling to
nucleons.

  * •

large neutrino bursts : The energy radiated as a result of collapse to a SM
neutron star is limited by the gravitational binding energy of the neutron
star. In contrast an IGY neutron star has additional binding energy due to the
attractive $\phi$ mediated force between nucleons which will lead to much
larger radiated energy if $\mathfrak{g}_{\mathrm{n}}\gtrsim 10$:
$E_{\mathrm{rad}}^{\mathrm{mx}}\approx 8.5\times
10^{53}\,\mathrm{erg}\,\times\,\frac{\mathfrak{g}_{\mathrm{n}}}{\mathfrak{g}_{\mathrm{vac}}}\,\left(\sqrt{\frac{\mathfrak{g}_{\mathrm{n}}}{\mathfrak{g}_{\mathrm{vac}}}}-1\right)$
(see §6.6 & fig. 12c). This will lead to larger neutrinos bursts produced
during collapse and cooling which would be observed by neutrino telescopes.

##  7 Generality of IGY Phenomenology

While the IGY EoS is an idealization which does not provide a realistic model
of matter throughout the entirety of any known astrophysical object the
phenomenology it reveals can be considered generic in the sense that the
unusual phenomenology it exhibits in stellar structure is purely a function of
a simple feature in the EoS, i.e. in the $P(\mathcal{E})$ curve. This feature
is what one may call a dip in the EoS which softens with increasing density
(as SM EoSs do) but then stiffens again. In SM stars one finds such a dip
between the densities where electrons become relativistic and where neutron
pressure starts to dominate which leads to the radius gap between neutron
stars and white dwarfs. In IGY an effective Yukawa coupling,
$\mathfrak{g}_{\mathrm{n}}\gtrsim 4$ adds an additional dip of sufficient
magnitude to cause a 2nd radius gap on the $R_{\star}\sim 10\,$km scale which
widens with larger $\mathfrak{g}_{\mathrm{n}}$. For
$\mathfrak{g}_{\mathrm{n}}\gtrsim 5$ the dip is large enough so that
$\frac{\partial P}{\partial\mathcal{E}}<0$ in the dip which leads to a density
gap between soft and hard matter allowing for hybrid stars which can have much
smaller radii. For $\mathfrak{g}_{\mathrm{n}}\gtrsim 5.7$ the even larger dip
has $P<0$, producing two classes of stars: hard and soft. Hard stars have a
minimum density and can be extremely small and are quite different from SM
stars.

An EoS dip occurs naturally when SM particles are coupled to each other by a
light scalar (force carrier). The light mass of the scalar allows the
superposition of the exterior $\phi$-field of many SM particles, the analog of
the electric field of charged particles, to coherently add to a large $\phi$
polarization when the SM particle density is large. Stellar remnants with
their large particle density is a natural place for this polarization to
manifest. The scalar polarization, $\tilde{\phi}$, will increase with density
and is accompanied by increasing negative pressure, $-V(\tilde{\phi})$, and
softening of $P(\mathcal{E})$. As the SM particles become relativistic the
$\phi$ polarization begins to saturate and so does $-V(\tilde{\phi})$ so the
EoS stiffens at higher densities. This is a result of two effects: 1) the
Yukawa coupled “$\phi$ charge” per particle is smaller for relativistic
particles: $\tilde{n}=n\,\langle\gamma^{-1}\rangle$ from eq. 4.22; 2) the
Yukawa coupling decreases the fermion kinematic mass as the polarization
increases causing the particles to become relativistic at lower densities:
$\tilde{m}_{f}=m_{f}-g_{f}\,\tilde{\phi}$ from eq. 4.2. This softening and
then stiffening is the dip which leads to the phenomenology described in this
paper. It derives from the form of a Yukawa term in the Lagrangian. Adding
more realism to IGY should not change these basic features.

Since the baryonic chemical potential is
$\mu_{\mathrm{b}}=(\varepsilon+P)/n_{\mathrm{b}}$ (eq. 3.5) a dip in $P$ is
associated with a dip in $\mu_{\mathrm{b}}$. This decrease in
$\mu_{\mathrm{b}}$ is due to tighter binding of the baryons caused by baryonic
self-attraction mediated by $\phi$. More tightly bound matter causes stars
carrying the same baryon number to be smaller in size. A larger dip
corresponds to larger scalar binding energy, and if large enough one has stars
which are primarily bound by the scalar force rather than by gravity as in
hard stars.

###  7.1 White Dwarfs

A somewhat trivial example of the generality of the IGY model comes from white
dwarfs with a Yukawa coupling to the electrons (or more generally to leptons).
Since the pressure is primarily from degenerate free electrons and the density
is primarily from non-relativisitic nuclei one has

$\displaystyle\mathcal{E}\approx
n_{\mathrm{b}}\,m_{\mathrm{n}}+\mathcal{E}_{\mathrm{IGY}}\qquad P\approx
P_{\mathrm{IGY}}$ (7.1)

where one takes $f\rightarrow\mathrm{e}$ in the IGY EoS and the baryon number
density is $n_{\mathrm{b}}=\sum_{(Z,A)}A\,n_{(Z,A)}$. Here $Z$ is the atomic
number and $A$ the atomic mass number and $n_{(Z,A)}$ is the number density of
specific nuclei. Charge neutrality requires
$n_{\mathrm{e}}=\sum_{(Z,A)}Z\,n_{(Z,A)}$. The majority white dwarfs have low
mass and isotopically consist almost entirely of ${}^{12}_{\ 6}$C, ${}^{16}_{\
8}$O, ${}^{20}_{10}$Ne, ${}^{24}_{12}$Mg each with $\frac{Z}{A}=\frac{1}{2}$.
In this case $n_{\mathrm{b}}\approx 2\,n_{\mathrm{e}}$. More massive white
dwarfs with Fe cores have $\overline{Z/A}$ no more than 8% lower. Thus the IGY
model gives an accurate pressure for a white dwarf but the energy dependence
of $P(\mathcal{E})$ is rescaled relative to that in fig. 1; roughly shifted by
a factor of $\sim 4000$ when the electrons are non-relativistic. For fixed
$n_{\mathrm{e}}/n_{\mathrm{b}}$ one can show that the values of
$\mathfrak{g}_{\mathrm{crit}}$, $\mathfrak{g}_{0}$ and
$\mathfrak{g}_{\mathrm{vac}}$ are unchanged by this rescaling. So the dips in
the EoS which lead to anomalous phenomenology occur at the same
$\mathfrak{g}_{\mathrm{e}}$ as the $\mathfrak{g}_{\mathrm{n}}$ in the IGY
neutron star model. Thus one should expect very similar phenomenology for
white dwarfs as for IGY neutron stars but of course masses and radii will
differ numerically. This will be worked out in a future work and compared to
astronomical observations of white dwarfs.

##  8 Amelioration of IGY Deficiencies

IGY is not realistic for neutron stars in a number of ways:

  * •

IGY only models a single conserved fermion which we take to be the neutron for
neutron stars. This is only a rough approximation over a narrow range of
densities since protons will carry a significant fraction of the baryon number
especially at densities not much lower than nuclear densities.

  * •

The protons which are not present in IGY are accompanied by electrons and
muons which dominate the pressure at low densities. This leads to a more
complicated $P(\mathcal{E})$ at low densities and $M_{\star}(R_{\star})$ for
small masses than is reflected in the IGY model.

  * •

IGY treats the neutrons as free particles with no interactions whereas it is
known that nucleons experience significant nuclear interactions at densities
above the saturation density which is at or below the density scale of neutron
stars. Thus IGY is certainly inaccurate for neutron stars just as is the OV
model of neutron stars.

One can ameliorate some of these deficiencies with extensions of the IGY
model.

Much of the formalism for an extending IGY to ameliorate these deficiencies
have been given in this paper. §4.1-§4.4 include the energy density and
pressure of multiple free fermionic species, with or without Yukawa couplings;
as well as allowing for non-quadratic scalar potentials. These generalizations
were dispensed with in the sections that followed in order to focus on a
simple, specific model with a semi-analytic EoS which exhibits qualitatively
the EoS and stellar structure phenomenology of adding a Yukawa coupling.

Interactions among fermions were not considered here but are important for
neutron stars where nuclear interactions play a role. The form of
$\mathcal{E}_{f}$ and $P_{f}$ from IGY can be used only so long as the
fermions move freely and not bound in localized structures such as quarks in
nucleons or nucleons in a lattice but even so one must add the contribution of
the interaction energy to $\mathcal{E}$. The interaction contribution to $P$
is then given by the adiabaticity equation. There is currently a great deal of
uncertainty SM nuclear interactions and a model where light scalar is coupled
to quarks and/or gluon further complicates modeling since the nuclear force
carriers are themselves made of quarks and gluons. Constructing a more
realistic Yukawa coupled EoS for neutron stars will be left to future work.

##  9 Implications

A summary of the results of this paper has already been given in §1.1. Here we
try to draw more general conclusions about the prospect for Yukawa coupled
stars both for fitting data and excluding light scalars coupled to SM matter.
This paper gives the qualitative phenomenology of Yukawa coupled cold remnant
stars with only rough quantitative predictions. We have identified three mass
regimes for scalars heavy, light and ultralight which depends on the type of
star under consideration. This paper focuses on the intermediate light scalar
mass range which is valid throughout a neutron star if
$10^{-10}\,\mathrm{eV}\lesssim m_{\phi}\lesssim 10^{3}\,\mathrm{eV}$ (eq.
4.17) and throughout a white dwarf if $10^{-13}\,\mathrm{eV}\lesssim
m_{\phi}\lesssim 10^{3}\,\mathrm{eV}$ (eq. 4.18). In this mass range we find
that for a star whose pressure is dominated by a single fermion species $f$
(neutrons for neutron stars and electrons for white dwarfs) the equation of
state and thus the stellar structure depends on the effective coupling
$\mathfrak{g}_{f}=\frac{g_{f}\,m_{f}}{m_{\phi}}$. Thus this is only a one
parameter family of stars resulting from the two parameters: the Yukawa
coupling $g_{f}$ and the scalar mass $m_{\phi}$. In our simple model various
anomalies in the stellar structure (see §6.7 and preceding sections) occur for
$\mathfrak{g}_{f}\gtrsim 3$. The onset of anomalous stellar structure will
occur for somewhat different values of $\mathfrak{g}_{f}$ for EoS models with
greater realism but we expect correcting the inadequacies listed in §8 will
only change these values by $\mathcal{O}(1)$. If no anomalous behavior is
observed then one should expect very roughly that this limits
$\mathfrak{g}_{f}\lesssim\mathcal{O}(1)$ or $g_{f}\lesssim
m_{\phi}/m_{f}$.101010If $V(\phi)=\frac{1}{2}\,{m_{\phi}}^{2}\,\phi^{2}$ this
bound on $g_{\mathrm{n}}$ is much weaker than laboratory bounds for all values
of $m_{\phi}$ consistent with the LSA. However with an appropriate choice of
$V(\phi)$ one can make the effective scalar mass, $\tilde{m}_{\phi}$, at the
nuclear densities of neutron stars much smaller than that of terrestrial
laboratories; greatly weakening the laboratory bounds relative to the bound
from stellar structure.. To obtain much stronger bounds from stellar structure
would require precision measurements of stellar parameters.

For large $\mathfrak{g}_{f}$ ($\geq 5.8$ in our simple model) the attractive
binding force from the scalar field allows for a high density state of matter
(hard matter) which is more tightly bound than free fermions and all low
density matter is metastable to conversion to hard matter. There is no direct
path for conversion of low density matter so this instability would not
manifest unless low density matter comes into contact with pre-existing seeds
of hard matter. It is conceivable that the dark matter is made of drops of
hard matter, a possibility which will be explored in future work.

## Acknowledgement

We thank Gordan Baym, Reed Essick and Ivan Esteban for useful discussions as
well as Paddy Fox for comments. CG is grateful to the generous support from
Aspen Center for Physics (ACP) where part of the work was done. ACP is
supported by National Science Foundation grant PHY-1607611. AS and CG are
supported by Fermilab, operated by the Fermi Research Alliance, LLC under
Contract DE-AC02-07CH11359 with the U.S. Dept. of Energy.

##  Appendix A $\langle\bar{\psi}\psi\rangle$ in Terms of Density and Pressure

Here we re-express $\langle\bar{\psi}\psi\rangle$ in terms of the density and
pressure contributed by a free fermion species. We may ignore space-time
curvature in solving $\langle\bar{\psi}\psi\rangle$, when a negligible
fraction of particles have very low momentum: $p\,R\lesssim 1$ (momentum $p$,
curvature scale $R$). This is true for all the SM fermions in a neutron star.
For flat space-time quantum field theory the Dirac field operator, $\psi$, is

$\displaystyle\begin{split}\psi&=\int\frac{d^{3}{\bf
p}}{(2\pi)^{3}}\frac{1}{\sqrt{2E}}\sum_{s}\left(a_{\bf p}^{s}u^{s}({\bf
p})e^{-i{p\cdot x}}+{b_{\bf p}^{s}}^{\dagger}v^{s}({\bf p})e^{i{p\cdot
x}}\right)\\\ \bar{\psi}&=\int\frac{d^{3}{\bf
p}}{(2\pi)^{3}}\frac{1}{\sqrt{2E}}\sum_{s}\left({a_{\bf
p}^{s}}^{\dagger}\bar{u}^{s}({\bf p})e^{i{p\cdot x}}+{b_{\bf
p}^{s}}\bar{v}^{s}({\bf p})e^{-i{p\cdot x}}\right)\end{split}$ (A.1)

where $E=\sqrt{|{\bf p}|^{2}+m^{2}}$, $m$ is the particle mass, $a_{\bf
p}^{s}$ (${a_{\bf p}^{s}}^{\dagger}$) annihilates (creates) a particle with
momentum ${\bf p}$ and spin $s$, and $b_{\bf p}^{s}$ (${b_{\bf
p}^{s}}^{\dagger}$) creates (annihilates) an anti-particles with momentum
${\bf p}$ and spin $s$. The creation and annihilation operators satisfy the
anticommutation relations:

$\displaystyle\\{a^{s}_{\bf p},{a^{s^{\prime}}_{\bf
p^{\prime}}}^{\dagger}\\}=\\{b^{s}_{\bf p},{b^{s^{\prime}}_{\bf
p^{\prime}}}^{\dagger}\\}=(2\pi)^{3}\delta^{3}({\bf p}-{\bf
p^{\prime}})\,\delta^{s\,s^{\prime}}$ (A.2)

so

$\displaystyle\begin{split}\bar{\psi}\psi\supset&\sum_{ss^{\prime}}\int\frac{d^{3}{\bf
p}}{(2\pi)^{3}}\frac{d^{3}{\bf
p^{\prime}}}{(2\pi)^{3}}\frac{1}{2\sqrt{E\,E^{\prime}}}\times\\\
&\left({a_{\bf p}^{s}}^{\dagger}a_{\bf
p^{\prime}}^{s^{\prime}}\bar{u}^{s}({\bf p})u^{s^{\prime}}({\bf
p^{\prime}})e^{-i{\bf(p-p^{\prime})\cdot x}+i(E-E^{\prime})t}+{b_{\bf
p}^{s}}{b_{\bf p^{\prime}}^{s^{\prime}}}^{\dagger}\bar{v}^{s}({\bf
p})v^{s^{\prime}}({\bf p^{\prime}})e^{i{\bf(p-p^{\prime})\cdot
x}-i(E-E^{\prime})t}\right)\end{split}$ (A.3)

Apply $V^{-1}\int dV$ to $\bar{\psi}\psi$,

$\displaystyle\begin{split}-\frac{1}{V}\int
dV\,\bar{\psi}\psi\sim&-\frac{1}{V}\int\frac{d^{3}{\bf
p}}{(2\pi)^{3}}\frac{1}{2E_{\bf p}}\sum_{ss^{\prime}}\left({a_{\bf
p}^{s}}^{\dagger}a_{\bf p}^{s^{\prime}}\bar{u}^{s}({\bf p})u^{s^{\prime}}({\bf
p})-{b_{\bf p}^{s^{\prime}}}^{\dagger}{b_{\bf p}^{s}}\bar{v}^{s}({\bf
p})v^{s^{\prime}}({\bf p})\right)\\\ =&-\frac{1}{V}\int\frac{d^{3}{\bf
p}}{(2\pi)^{3}}\frac{2m}{2E_{\bf p}}\sum_{s}\left({a_{\bf
p}^{s}}^{\dagger}a_{\bf p}^{s^{\prime}}+{b_{\bf
p}^{s^{\prime}}}^{\dagger}{b_{\bf p}^{s}}\right)\\\
=&-\frac{1}{Vm}\int\frac{d^{3}{\bf p}}{(2\pi)^{3}}\left(E_{\bf p}-\frac{{\bf
p}^{2}}{E_{\bf p}}\right)\sum_{s}\left({a_{\bf p}^{s}}^{\dagger}a_{\bf
p}^{s^{\prime}}+{b_{\bf p}^{s^{\prime}}}^{\dagger}{b_{\bf
p}^{s}}\right)\end{split}$ (A.4)

where we have used $\bar{u}^{s}({\bf p})u^{s^{\prime}}({\bf
p})=-\bar{v}^{s}({\bf p})v^{s^{\prime}}({\bf p})=2m\delta^{ss^{\prime}}$ in
the second line. Therefore, $\langle\bar{\psi}{\psi}\rangle$ approximately
measures $(\mathcal{E}-3P)/m$, assuming an ideal Fermi gas. In the non-
relativistic limit, this is approximately equal to the number density of
fermion $n$.

##  Appendix B TOV Equations with Anisotropic Pressure

Consider a spherically symmetric metric,

$\displaystyle
ds^{2}=-e^{2\Phi}dt^{2}+\left(1-\frac{2Gm}{r}\right)^{-1}dr^{2}+r^{2}d\Omega^{2}$
(B.1)

where $\Phi=\Phi(r),m=m(r)$ under hydrostatic equilibrium. The Einstein tensor
is given by

$\displaystyle G_{\mu\nu}\equiv R_{\mu\nu}-\frac{1}{2}Rg_{\mu\nu}$ (B.2)

In the basis of $(\hat{t},\hat{r},\hat{\theta},\hat{\phi})$,

$\displaystyle\begin{split}{G^{t}}_{t}&=-\frac{2Gm^{\prime}}{r^{2}}\\\
{G^{r}}_{r}&=\frac{-2Gm+2r(r-2Gm)\Phi^{\prime}}{r^{3}}\\\
{G^{\theta}}_{\theta}&=\frac{(1+r\Phi^{\prime})\left(Gm-
Gm^{\prime}r+r(r-2Gm)\Phi^{\prime}\right)+r^{2}(r-2Gm)\Phi^{\prime\prime}}{r^{3}}\\\
{G^{\phi}}_{\phi}&={G^{\theta}}_{\theta}\end{split}$ (B.3)

The Einstein equations are

$\displaystyle{G^{\mu}}_{\nu}=8\pi G{T^{\mu}}_{\nu}$ (B.4)

An important property of the stress-energy tensor is that it satisfies the
continuity equation

$\displaystyle\nabla_{\mu}{T^{\mu}}_{\nu}=0$ (B.5)

We consider a system that is static and spherically symmetric, described by
the metric (B.1). Due to the symmetry, the $t$’s and $\phi$’s components of
the continuity equation are automatically satisfied. The $\theta$’s component
yields ${T^{\theta}}_{\theta}={T^{\phi}}_{\phi}.$ The $r$’s component yields

$\displaystyle
0=\partial_{r}{T^{r}}_{r}+\Phi^{\prime}\left({T^{r}}_{r}-{T^{t}}_{t}\right)+\frac{2}{r}\left({T^{r}}_{r}-{T^{\theta}}_{\theta}\right)$
(B.6)

If we identify
$T^{t}_{t}=-\mathcal{E},T^{r}_{r}=P_{\parallel},T^{\theta}_{\theta}=P_{\perp}$,
the continuity equation (B.6) becomes

$\displaystyle
0=\partial_{r}P_{\parallel}+\Phi^{\prime}\left(P_{\parallel}+\mathcal{E}\right)+\frac{2}{r}\left(P_{\parallel}-P_{\perp}\right)$
(B.7)

Now we can write down the TOV equations. Because of the continuity property of
the energy stress tensor, not all components of Einstein equations are
independent. Therefore, using (B.7), the Einstein equations (B.4) become

$\displaystyle\begin{split}\frac{dm}{dr}&=4\pi r^{2}\mathcal{E}\\\
\frac{dP_{\parallel}}{dr}&=-\frac{\mathcal{E}Gm}{r^{2}}\left(\frac{P_{\parallel}}{\mathcal{E}}+1\right)\left(\frac{4\pi
P_{\parallel}r^{3}}{m}+1\right)\left(1-\frac{2Gm}{r}\right)^{-1}-\frac{2}{r}\left(P_{\parallel}-P_{\perp}\right)\end{split}$
(B.8)

In the case of a perfect fluid, $P_{\parallel}=P_{\perp}$, it simplifies to

$\displaystyle\begin{split}\frac{dm}{dr}=&4\pi r^{2}\mathcal{E}\\\
\frac{dP}{dr}=&-\frac{\rho Gm}{r^{2}}\left(\frac{4\pi
Pr^{3}}{m}+1\right)\left(\frac{P}{\mathcal{E}}+1\right)\left(1-\frac{2Gm}{r}\right)^{-1}\end{split}$
(B.9)

##  Appendix C Ideal Degenerate Fermi Gas

To solve TOV, one needs to find the relation between $\mathcal{E}$ and
$P_{\parallel},P_{\perp}$, a.k.a. the equation of state (EoS). For an ideal
fermi gas in equilibrium, the distribution function is

$\displaystyle f(E)=\frac{1}{\textrm{exp}[(E-\mu)/kT]+1}$ (C.1)

For completely degenerate fermions ($T\to 0$), $\mu(T=0)$ is called the Fermi
energy $E_{F}$, and

$\displaystyle f(E)=\left\\{\begin{array}[]{cc}1&E\leq E_{F}\\\
0&E>E_{F}\end{array}\right.$ (C.4)

Define Fermi momentum $p_{F}$ as

$\displaystyle E_{F}\equiv\left(p_{F}^{2}+m^{2}\right)^{1/2}$ (C.5)

The number density of the degenerate ideal fermi gas is

$\displaystyle n=2\int_{0}^{p_{F}}4\pi p^{2}dp=\frac{1}{3\pi^{2}}m^{3}x^{3}$
(C.6)

where we use the “relativity parameter”, $x$, define as

$\displaystyle x\equiv p_{F}/m$ (C.7)

The energy density is

$\displaystyle\begin{split}\mathcal{E}=2\int_{0}^{p_{F}}\left(p^{2}+m^{2}\right)^{1/2}4\pi
p^{2}dp=m^{4}\xi(x)\\\
\xi(x)=\frac{1}{8\pi^{2}}\left\\{x(1+x^{2})^{1/2}(1+2x^{2})-\ln\left[x+(1+x^{2})^{1/2}\right]\right\\}\end{split}$
(C.8)

The pressure is

$\displaystyle\begin{split}P=\frac{1}{3}2\int_{0}^{p_{F}}\frac{p^{2}}{\left(p^{2}+m^{2}\right)^{1/2}}4\pi
p^{2}dp=m^{4}\varphi(x)\\\
\varphi(x)=\frac{1}{8\pi^{2}}\left\\{x(1+x^{2})^{1/2}(2x^{2}/3-1)+\ln\left[x+(1+x^{2})^{1/2}\right]\right\\}\end{split}$
(C.9)

Clearly, by scanning over x, one can numerically obtain an EoS
$\mathcal{E}=\mathcal{E}(P)$.

##  Appendix D Simpler Analytics for IGY

An algebraically and numerically simpler form of the IGY EoS of eq.s 4.5.1 is
given by

$\displaystyle\frac{\mathcal{E}_{\mathrm{IGY}}(x)}{m_{f}^{4}}$
$\displaystyle=\mathcal{A}(x)+\mathcal{B}(x)$
$\displaystyle\frac{n_{f}(x)}{m_{f}^{3}}$
$\displaystyle=r_{\mathrm{m}}(x)^{3}\,\frac{x^{3}}{3\,\pi^{2}}$
$\displaystyle\mathcal{A}(x)$
$\displaystyle=r_{\mathrm{m}}(x)^{4}\,\frac{x^{3}\sqrt{1+x^{2}}}{4\pi^{2}}$
$\displaystyle\frac{P_{\mathrm{IGY}}(x)}{m_{f}^{4}}$
$\displaystyle=\frac{\mathcal{A}(x)}{3}-\mathcal{B}(x)$
$\displaystyle\frac{\tilde{n}_{f}(x)}{m_{f}^{3}}$
$\displaystyle=\frac{1-r_{\mathrm{m}}(x)}{{\mathfrak{g}_{f}}^{2}}$
$\displaystyle\mathcal{B}(x)$