# Holographic cold dense matter constrained by neutron stars

Lin Zhang<EMAIL_ADDRESS>School of Nuclear Science and Technology,
University of Chinese Academy of Sciences,
Beijing, P.R.China 100049  Mei Huang<EMAIL_ADDRESS>School of Nuclear
Science and Technology, University of Chinese Academy of Sciences,
Beijing, P.R.China 100049

###### Abstract

The equation of state (EoS) for cold dense matter inside neutron stars is
investigated by using holographic QCD models in the framework of the Einstein-
Maxwell-dilaton (EMD) system and the improved Karch-Katz-Son-Stephanov (KKSS)
action for matter part. This method of describing holographic nuclear matter
in the EMD$+$KKSS framework is different from that by using the Dirac-Born-
Infeld (DBI) action and the Chern-Simons (CS) terms. Combining with the
Hebeler-Lattimer-Pethick-Schwenk (HLPS) intermediate equation of state (EoS),
the hybrid EoS inside the neutron stars is constructed. The obtained hybrid
EoS is located in the range that is defined by the low-density chiral
effective theory, the high-density perturbative QCD, and the polytropic
interpolations between them, and is constrained by the astrophysics
observations. The square of the sound velocity reaches a maximum value larger
than $0.8$ in the region of $2-5$ times the saturation baryon number density
and approaches the conformal limit at the high baryon density range. The mass-
radius relation and the tidal deformability of the neutron stars are in
agreement with astrophysical measurements. The possible maximum mass for the
neutron star is about $2.5M_{\odot}$ and the radius is about $12\mathrm{km}$
then. It is noticed that the holographic quark matter branch in the mass-
radius relation is always unstable and the holographic nuclear matter can
produce a stable branch. These results indicate that even in the core of the
NS, the matter is still in the confinement phase and the quark matter is not
favored.

##  I Introduction

Quantum chromodynamics (QCD) is believed to be the basic theory of the strong
interaction. However, the strong interaction is non-perturbative at low energy
scale. Thus many phenomena such as hadron and glueball spectra, chiral phase
transition and confinement-deconfinement phase transition, the equation of
state (EoS) of the strong interaction matter under extreme conditions (high
temperature, high baryon number density, rapid rotation, and strong magnetic
field) can’t be explained from the first principle until now. Therefore many
efforts have been paid on developing non-perturbative methods, such as lattice
QCD [1, 2, 3], chiral effective theory [4], Nambu–Jona-Lasinio (NJL) model [5,
6], QCD sum rule [7, 8, 9, 10, 11, 12, 13, 14], Dyson Schwinger equations
(DSE) and functional renormalization group (FRG) [15], and recently
holographic QCD method [16, 17, 18, 19, 20, 21, 22, 23, 24].

Exploring the equation of state (EoS) and phase structure of quantum
chromodynamics (QCD) at high density is boosted both by neutron star (NS) and
heavy ion collisions. Relativistic heavy ion collision (RHIC) can produce the
quark-gluon plasma (QGP) with high temperature and/or high baryon number
density, and strong magnetic field. Besides that, the QGP with rapid rotation
can also be produced in non-central heavy ion collisions. That is the most
extreme environment created in the laboratory. It has been theoretically
predicted that there exists a critical end point (CEP) for QCD at finite
baryon density. The search of the CEP has been one of the most important
targets for the beam energy scan program at RHIC as well as for the future
accelerator facilities at Facility for Antiproton and Ion Research (FAIR) in
Darmstadt and Nuclotron-based Ion Collider Facility (NICA) in Dubna.

In the nature, the early universe and the compact stars are also in extreme
conditions. All these provide good platforms to investigate the non-
perturbative features of strong interaction and the properties of the QCD
matter.

The neutron star (NS) is a kind of compact star, which is the remnant after a
massive super-giant star collapses. From August 1967 on, when the existence of
NSs was confirmed by the discovery of radio pulsars, more than 2700 radio
pulsars have been detected. It has been wondered for more than a half-century
what’s the internal structure of NSs. Whether quark matter exists inside the
core of NSs [25, 26]? Is there any exotic matter in NSs? Do strange stars
exist? We can’t answer these questions yet. There are just a few of the open
questions [27]. The strong interaction matter in the core of the NS is
extremely dense, the baryon number density $n_{B,\mathrm{NS}}$ of which can
reach as high as $5n_{S}-10n_{S}$, where $n_{S}$ is the nuclear saturation
density. Although the surface of the NS is hot and the temperature can reach
several hundred thousand Kelvin [28, 29], compared with the QCD energy scale
$\Lambda_{\mathrm{QCD}}$ this temperature is still very low and can be
neglected. Thus we can regard the matter in the NS as cold dense QCD matter.

On the experimental side, there are some observations of the NSs, including
the mass measurements through the general relativistic Shapiro delay [30] and
other methods [31, 32, 33, 34, 35, 36], the mass and radius measurements
through the detection of thermal emission from the surface of the stars [37,
38, 39, 40], and the radius, EoS, and tidal deformability measurements through
the gravitational wave signals of the binary systems [41, 42]. We will discuss
these observations in detail in Sec. II.

Recent developments in the field of observational astronomy, e.g. the
currently operating satellites, the Neutron star Interior Composition Explorer
(NICER), and the gravitational-wave laser interferometers (Advanced LIGO,
Virgo, and KAGRA) have provided tight constraints on neutron star masses and
radii, which in turn, put strict constraints on the equation of state of dense
matter.

The difficulties to predict the cold QCD matter from the theoretical side come
from the awkward energy scale of the QCD matter in NS $\Lambda_{\mathrm{NS}}$:
the baryon number density in the neutron star $n_{B,\mathrm{NS}}$ is high so
that $\Lambda_{\mathrm{NS}}$ exceeds the appropriate range where the
traditional nuclear physics theories validate; however it is still not high
enough to adapt the perturbative method.

Although it’s very hard to directly study the EoS of cold dense QCD matter
inside the NSs, we can still derive some useful information when we combine
the results of the traditional nuclear theories in the low energy scale and
the results of the perturbative QCD at the high energy scale with the
observation data. Under the constraints from the observation data, there are
two approaches to derive the allowed range of the EoS for the cold QCD matter
in the whole range of $n_{B}$: the “direct interpolation approach” [43, 44]
and the “causality and thermodynamic stability constraint approach” [45, 46].
In the former approach, combining the observation data, a generic family of
neutron-star-matter EoSs, which interpolate between theoretical results at low
and high baryon density, are obtained. In the latter one, utilizing the
thermodynamic stability and causality, it is found that $ab$-$initio$
calculations in QCD at high densities, as complementary to the astrophysical
observations, provide significant and nontrivial information about the EoS of
matter in the cores of the neutron star.

To directly calculate the NS interior composition, from the theoretical side,
normally Nambu-Jona-Lasinio (NJL) model [47, 48, 49, 50, 51, 52, 53] has been
used for quark matter. The holographic QCD model [54, 55, 56, 57, 58, 59, 60,
61, 62, 63, 64, 65, 66, 67, 68, 69] has been developed at finite baryon
density over recent years.

There are many efforts to describe the EoS of cold dense matter inside the
neutron star in the holographic frame. The main approach is to use the
background action in the holographic models to describe the holographic quark
matter. Then the phase equilibrium conditions can be used to connect the EoS
of the holographic quark matter and the EoS calculated in effective theories
adapted at low energy scale. This kind of hybrid EoSs can describe the NSs.
However, the problem is that the branch that corresponds to the holographic
quark matter in the mass-radius relation is unstable according to the Bardeen-
Thorne-Meltzer (BTM) stable criterion stated in Refs. [70, 71]. The BTM stable
criterion is stated in the following. As one moves along the $M$-$R$ curve,
where $M$ is the mass of the NS as one observes from outside and $R$ is its
radius, to make sure the value of the internal energy density at the center of
the star increases, all configurations remain stable until one reach the first
extreme value (maximum or minimum) of $M$. The fundamental mode of radial
oscillation becomes unstable at the first extreme value point. Between the
first extreme value point and the second one, all configurations contain one
unstable radial mode. At the second extreme value point, and also at each
succeeding extreme value point, the stability of one more mode of radial
oscillation changes. How the $M$-$R$ curve bends as it passes through the
extreme value point determines how the stability changes. One previously
stable radial mode becomes unstable at the critical point if the direction of
the bend is counter-clockwise (with increasing $\epsilon_{0}$). One previously
unstable mode becomes stable at the critical point if the direction of the
bend is clockwise.

The key point to produce the stable branch in the holographic method is to
consider the holographic nuclear matter. In Refs. [72, 60, 61, 62, 63, 64, 65,
66], the Dirac-Born-Infeld (DBI) action and the Chern-Simons (CS) terms has
been used to describe the baryons, where the solitons are regarded as baryons.
The authors in Refs. [60, 61, 62] used the V-QCD bottom-up model action [73,
74, 75] to describe the background. In Refs. [63, 64, 65], the authors chose
the Witten-Sakai-Sugimoto model [76, 20, 21] in top-down approach as the
background. In Ref. [66], the authors used the hard-wall AdS/QCD model action
to describe the background. They all obtained the stable branch that
corresponds to the holographic nuclear matter.

In this work, different from the topological description of the baryons, we
use the improved Karch-Katz-Son-Stephanov (KKSS) action [24, 77, 78, 79, 80]
to introduce the contribution of the hadrons and treat it as a probe. Then we
can calculate the EoS of the holographic nuclear matter. This article is
organized in the following. In Sec. II, we briefly summarize the observation
data of the mass and radius of the NSs. In Sec. III, we first introduce three
holographic models. Then, we solve the Einstein-Maxwell-dilaton (EMD) system
to derive the thermodynamics of the background and plot the $T-\mu$ phase
diagram for the confinement-deconfinement phase transition in the three
models. We also compare the deformed background metric in the three models
with that in the $AdS_{5}$ space-time. In Sec. IV we consider the
thermodynamic contribution made by the KKSS action and calculate the EoS of
the holographic nuclear matter. Combined with the Hebeler-Lattimer-Pethick-
Schwenk (HLPS) intermediate EoS, we then construct the hybrid EoS of the cold
strong interaction matter. In Sec. V, we use the hybrid EoSs to calculate some
properties of the neutron stars, such as mass-radius relation and tidal
deformability. We also compare our results with the observation data. Finally,
we give the conclusion and some discussions in Sec. VI.

##  II The measurement of the mass and radius of neutron stars

Here we give a summary of the astrophysical observations for NSs [81, 62].

PSR J1614-2230 is a recycled millisecond pulsar with a very low inferred
magnetic field, the magnetic induction intensity of which is $B=1.8\times
10^{8}\mathrm{G}$. It was discovered in a completed search of 56 unidentified
mid-latitude EGRET $\gamma$-ray sources for pulsations using the Parkes
telescope and the multibeam receiver in 2004 [82]. In 2010, the authors in
Ref. [31] presented radio timing observations of the binary millisecond pulsar
PSR J1614-2230. The detection of the general relativistic Shapiro delay [30]
implied the pulsar mass is $1.97\pm 0.04M_{\odot}$, which made it the heaviest
known NS at that time. In 2016, the authors in Ref. [32] analyzed 24 binary
radio pulsars in the North American Nanohertz Observatory for Gravitational
Waves (NANOGrav) nine-year data set and made 14 significant measurements of
the Shapiro delay. They presented the mass of PSR J1614-2230 is
$1.928_{-0.017}^{+0.017}M_{\odot}$ at $68.3\%$ credibility.

PSR J0348+0432 is a mildly partially recycled binary pulsar with
$P=39.1\mathrm{ms}$ and low magnetic field ($B_{\mathrm{surf}}=3.1\times
10^{9}\mathrm{G}$), in a $2.4\mathrm{hr}$ orbit with a low-mass companion,
which was detected in 2012 in a $350\mathrm{MHz}$ Drift-scan Survey using the
Robert C. Byrd Green Bank Telescope with the goal of finding new radio
pulsars, especially millisecond pulsars that can be timed to high precision
[83]. PSR J0348+0432 is only the second NS with a precisely determined mass
around $2M_{\odot}$, after PSR J1614-2230. In 2013, the authors in Ref. [33]
reported on radio-timing observations of the pulsar J0348+0432 and phase-
resolved optical spectroscopy of its white-dwarf companion and then derived a
NS mass in the range from $1.97M_{\odot}$ to $2.05M_{\odot}$ at $68.27\%$
confidence or $1.90M_{\odot}$ to $2.18M_{\odot}$ at $99.73\%$ confidence.

In 2014, 67 new radio pulsars including J0741+66
($P=0.00288570374(45)\mathrm{s}$) were discovered in the Green Bank Northern
Celestial Cap (GBNCC) survey, which surveys the entire sky visible to the
Robert C. Byrd Green Bank Radio Telescope (GBT) at 350 MHz for radio pulsars
and fast radio transients [84]. In 2018, the pulsar J0741+66 was renamed
J0740+6620 in Ref. [85]. In 2019, by combining data from the North American
Nanohertz Observatory for Gravitational Waves (NANOGrav) 12.5-yr data set with
orbital phase-specific observations using the Green Bank Telescope, the
authors in Ref. [34] measured the mass of the millisecond pulsars (MSP)
J0740+6620 to be $2.14_{-0.09}^{+0.10}M_{\odot}$ at $68.3\%$ credibility
interval or $2.14_{-0.18}^{+0.20}M_{\odot}$ at $95.4\%$ credibility interval.
In 2021, the authors in Ref. [35] used the data set that consisted of combined
pulse arrival-time measurements made with the $100\mathrm{m}$ Green Bank
Telescope and the Canadian Hydrogen Intensity Mapping Experiment telescope.
They presented the model-averaged result of the mass of PSR J0740+6620 is
$2.08_{-0.07}^{+0.07}M_{\odot}$ at $68.3\%$ credible interval, which is
reduced by about $2.8\%$ compared to the result given in Ref. [34]. This
model-averaged value remains the largest of all other precisely measured
pulsar masses determined with the Shapiro delay to date.

There is also a special class of MSPs called black widows (and redbacks), the
masses of which are higher. In black-widow systems, a millisecond pulsar is
accompanied by a low-mass, few $0.01M_{\odot}$ companion, which is bloated and
strongly irradiated by the pulsar, leading to outflows strong enough to
eclipse the pulsar signal for significant fractions of the orbit [36]. The
authors in Ref. [36] present evidence that the black-widow pulsar, PSR
B1957+20, has a high mass. Their best-fit pulsar mass is $2.40\pm
0.12M_{\odot}$. However, as it was stated in Ref. [81], there are many
difficulties to obtain accurate measurements from these ablated companions,
and the possibility of even more massive NSs than J0348+0432 is not yet as
robust as the results from radio timing.

These mass measurements provide a strict restriction on the EoS of the cold
QCD matter, which can’t be too “soft” and should lead to the maximum mass of
NSs is about $2$ times the solar mass. Otherwise, the high masses of the above
NSs cannot be explained. Many soft EoSs have already been ruled out through
these mass measurements.

The measurements of the radii of the NSs are more difficult than that of the
masses. In the past decade, several different techniques have been employed in
the NS radius measurements. Almost all of the methods that are currently used
base on the detection of thermal emission from the surface of the stars.
People either measure its apparent angular size or detect the effects of the
NS space-time on this emission to extract the information of the radius. These
methods can be divided into two kinds: spectroscopic measurements and timing
measurements. The ongoing NASA’s NS Interior Composition Explorer Mission
(NICER) [86] uses the latter method (pulse profile modeling) where a large
number of counts from a small number of rotation-powered pulsars which emit
thermally from their polar caps are collected. Combined with the refined
theoretical models that have been developed to analyze these waveforms, those
data will provide the measurement of radii and EoSs of the NSs. The results
have been published in Refs. [37, 38, 39, 40], in which they also use the data
from the XMM-Newton telescope.

In 2019, the authors in Ref. [37] derived that for the isolated
$205.53\mathrm{Hz}$ millisecond pulsar PSR J0030+0451, the estimated
equatorial circumferential radius is
$R_{\mathrm{e}}=13.02_{-1.06}^{+1.24}\mathrm{km}$ and the gravitational mass
is $M=1.44_{-0.14}^{+0.15}M_{\odot}$ at $68\%$ credible level. In the same
year, another group presented different values for PSR J0030+0451, inferred
mass $M$ and equatorial radius $R_{\mathrm{eq}}$ of which are
$1.34_{-0.16}^{+0.15}M_{\odot}$ and $12.71_{-1.19}^{+1.14}\mathrm{km}$,
respectively, where the credible interval bounds are approximately the $16\%$
and $84\%$ quantiles in marginal posterior mass, given relative to the median
[38].

In 2021, the authors in Ref. [39] found that the equatorial circumferential
radius of PSR J0740+6620 is $13.7_{-1.5}^{+2.6}\mathrm{km}$ at $68\%$
credibility based on fits of rotating hot spot patterns to NS Interior
Composition Explorer (NICER) and X-ray Multi-Mirror (XMM-Newton) X-ray
observations. On the same day, another group constrained the equatorial radius
and mass of PSR J0740+6620 to be $12.39_{-0.98}^{+1.30}\mathrm{km}$ and
$2.072_{-0.066}^{+0.067}M_{\odot}$, respectively, each reported as the
posterior credible interval bounded by the $16\%$ and $84\%$ quantiles [40].

There is another way to study the NSs, in which people observe the binary
merger events involving NSs. On August 17, 2017, the Advanced LIGO and
Advanced Virgo gravitational-wave detectors made their first observation of a
binary NS inspiral [87]. The signal, GW170817, was detected. Thereafter the
multi-messenger observations of this binary NS merger are finished by
telescopes and observatories ranging basically over the whole spectrum of
electromagnetic waves [88]. In 2018, under the hypothesis that both bodies in
the GW170817 event are NSs that are described by the same equation of state
and have spins within the range observed in Galactic binary NSs, the authors
in Refs. [41] employed two methods to give an analysis of the radii and EoS of
the NSs. From the LIGO and Virgo data alone and using equation-of-state-
insensitive relations between various macroscopic properties of the NSs, they
measured the two NS radii as $R_{1}=10.8_{-1.7}^{+2.0}\mathrm{km}$ for the
heavier one and $R_{2}=10.7_{-1.5}^{+2.1}\mathrm{km}$ for the lighter one at
the $90\%$ credible level. Using an efficient parametrization of the defining
function $P(\epsilon)$ of the equation of state itself and requiring that the
equation of state supports NSs with masses larger than $1.97M_{\odot}$, they
further derived $R_{1}=11.9_{-1.4}^{+1.4}\mathrm{km}$ for the heavier one and
$R_{2}=11.9_{-1.4}^{+1.4}\mathrm{km}$ for the lighter one at the $90\%$
credible level. They also obtained constraints on $P(\epsilon)$ at
supranuclear densities, with pressure at twice nuclear saturation density
measured at $3.5_{-1.7}^{+2.7}\times
10^{34}\,\mathrm{dyn}\,{\mathrm{cm}}^{-2}$ at the $90\%$ level. The tidal
deformability of a $1.4M_{\odot}$-NS could be estimated to be
$\Lambda_{1.4}=190_{-120}^{+390}$ at the $90\%$ level when a common but
unknown EoS is imposed. This measurement value of $\Lambda_{1.4}$ is very
interesting: contrary to the constraint from the mass measurement, it suggests
that “soft” EoSs such as APR4 are favored over “stiff” EoSs such as H4 or MS1.
So for an EoS that satisfies both the mass measurement and the tidal
deformability measurement, it should be stiff but not too stiff.

First identified in data from two detectors, LIGO Livingston and Virgo on 2019
August 14, GW190814 is the signal of a compact binary coalescence with the
most unequal mass ratio yet measured with gravitational waves:
$q=0.112_{-0.009}^{+0.008}$ [42]. The primary component of GW190814 is
conclusively a black-hole (BH) with mass $m_{1}=23.2_{-1.0}^{+1.1}M_{\odot}$.
Due to the lack of measurable tidal deformations and the absence of an
electromagnetic counterpart, the secondary component may be an NS or a BH.
However, comparing with the maximum NS mass predicted by studies of GW170817’s
remnant indicate it’s an NS by current knowledge of the NS equation of state,
and by electromagnetic observations of NSs in binary systems [42]. If it’s
true, this source with the mass of $2.59_{-0.09}^{+0.08}M_{\odot}$ is the
heaviest NS ever observed in a double compact-object system. Then the
corresponding constraints on the radius and tidal deformability of a canonical
$1.4M_{\odot}$-NS are $R_{1.4}=12.9_{-0.7}^{+0.8}\mathrm{km}$ and
$\Lambda_{1.4}=616_{-158}^{+273}$ [42].

Finally, we conclude all the data in this section in Tab. 1 and Tab. 2.

pulsar |  mass $M$ ($M_{\odot}$) |  radius $R$ ($\mathrm{km}$)
---|---|---
J1614-2230 |  $1.97\pm 0.04$ [30] | –
|  $1.928_{-0.017}^{+0.017}$ [32] |
J0348+0432 |  $1.97-2.05$ [33] | –
J0740+6620 |  $2.14_{-0.09}^{+0.10}$ [34] | –
|  $2.08_{-0.07}^{+0.07}$ [35] |
B1957+20 |  $2.40\pm 0.12$ [36] | –
J0030+0451 |  $1.44_{-0.14}^{+0.15}$ [37] |  $13.02_{-1.06}^{+1.24}$ [37]
|  $1.34_{-0.16}^{+0.15}$ [38] |  $12.71_{-1.19}^{+1.14}$ [38]
J0740+6620 |  |  $13.7_{-1.5}^{+2.6}$ [39]
|  $2.072_{-0.066}^{+0.067}$ [40] |  $12.39_{-0.98}^{+1.30}$ [40]
Table 1: The mass and radius for several discovered pulsars. gravitational-wave signal |  mass $M$ ($M_{\odot}$) |  radius $R$ ($\mathrm{km}$) |  tidal deformability $\Lambda_{1.4}$
---|---|---|---
GW170817 | – |  $R_{1}=10.8_{-1.7}^{+2.0}$, $R_{2}=10.7_{-1.5}^{+2.1}$ [41] |  $190_{-120}^{+390}$ [41]
|  |  $R_{1}=11.9_{-1.4}^{+1.4}$, $R_{2}=11.9_{-1.4}^{+1.4}$ [41] |
GW190814 |  $m_{1}=23.2_{-1.0}^{+1.1}$, $m_{2}=2.59_{-0.09}^{+0.08}$ [42] | – |  $616_{-158}^{+273}$ [42]
Table 2: Some discovery parameters and derived parameters for two binary
compact object systems from the gravitational-wave detection.

##  III The holographic framework and the background of the holographic models

Firstly we introduce the framework of the Einstein-Maxwell-dilaton (EMD)
system, which comes back to the gravity-dilaton coupling system at zero
chemical potential. As the background action, the EMD system is widely used in
the holographic QCD models. It can be used to describe the QCD confinement-
deconfinement phase diagram [89, 90], the effect of the magnetic field on the
QCD phase transition [91], the rotating effect on deconfinement phase
transition [92], the thermodynamic properties of QCD and the location of the
critical endpoint (CEP) [93], the critical behavior near the CEP [94, 95], and
the spectra glueballs and corresponding EoS [96].

Including hadrons excitation, the total action of the $5$-dimensional
holographic QCD model is

$\displaystyle S_{\text{total }}^{s}=S_{b}^{s}+S_{m}^{s},$ (1)

where $S_{b}^{s}$ is the action for the background in the string frame, and
$S_{m}^{s}$ is the action that describes hadrons in the string frame.

The EMD action in the string frame is denoted by $S_{b}^{s}$:

$\displaystyle S_{b}^{s}=\frac{1}{2\kappa_{5}^{2}}$ $\displaystyle\int$
$\displaystyle\mathrm{d}^{5}x\sqrt{-g^{s}}e^{-2\Phi}\Big{[}R^{s}+4{g^{s}}^{MN}\partial_{M}\Phi\partial_{N}\Phi$
$\displaystyle-V^{s}(\Phi)-\frac{h(\Phi)}{4}e^{\frac{4\Phi}{3}}{g^{s}}^{M\widetilde{M}}{g^{s}}^{N\widetilde{N}}F_{MN}F_{\widetilde{M}\widetilde{N}}\Big{]},$
(2)

where $s$ denotes the string frame, $\kappa_{5}^{2}=8\pi G_{5}$, $G_{5}$ is
the $5$-dimensional Newton constant. $g^{s}$ is the determinant of the metric
in the string frame: $g^{s}=\det\left(g_{MN}\right)$. The metric tensor in the
string frame can be extracted from

$\displaystyle
ds^{2}=\frac{L^{2}e^{2A_{s}(z)}}{z^{2}}\bigg{(}-f(z)dt^{2}+\frac{dz^{2}}{f(z)}+dy_{1}^{2}$
$\displaystyle\qquad\qquad\qquad\qquad+dy_{3}^{2}+dy_{3}^{2}\bigg{)},$ (3)

where $L$ is the curvature radius of the asymptotic $AdS_{5}$ space-time. For
simplicity and without loss of generality, we assume $L=1$ in the following
calculations. $R^{s}$ is the Ricci curvature scalar in the string frame. The
$5$-dimensional scalar field $\Phi(z)$ is the dilaton field that depends only
on the coordinate $z$. $F_{MN}$ is the field strength of the $U(1)$ gauge
field $A_{M}$:

$\displaystyle F_{MN}=\partial_{M}A_{N}-\partial_{N}A_{M}.$ (4)

The $5$-dimensional vector field $A_{M}$ is dual to the baryon number current.
The function $h(\Phi)$ describes the coupling strength of $A_{M}$ in the
theory, and the function $V^{s}(\Phi)$ represents the potential of the dilaton
field in the string frame. $h(\Phi)$ and $V^{s}(\Phi)$ are the functions that
depend only on the value of $\Phi$.

###  III.1 The Einstein-Maxwell-dilaton system

As discussed in Ref. [97], it is more convenient to choose the string frame
when solving the vacuum expectation value of the loop operator, and it is more
convenient to choose the Einstein frame to work out the gravity solution and
to study the equation of state. So we apply the following Weyl transformation
[98, 99] to Eq. (2):

$\displaystyle g_{MN}^{s}=\mathrm{e}^{\frac{4}{3}\Phi}g_{MN}^{E},$ (5)

where $g_{MN}^{E}$ is the metric tensor in the Einstein frame. The capital
letter ’E’ denotes the Einstein frame. Then, Eq. (2) becomes

$\displaystyle S^{E}=\frac{1}{2\kappa_{5}^{2}}$ $\displaystyle\int$
$\displaystyle\mathrm{d}^{5}x\sqrt{-g^{E}}\Big{[}R^{E}-\frac{4}{3}{g^{E}}^{MN}\partial_{M}\Phi\partial_{N}\Phi$
$\displaystyle-V^{E}(\Phi)-\frac{h(\Phi)}{4}{g^{E}}^{M\widetilde{M}}{g^{E}}^{N\widetilde{N}}F_{MN}F_{\widetilde{M}\widetilde{N}}\Big{]},$
(6)

where the function $V^{E}=\mathrm{e}^{\frac{4}{3}\Phi}V^{s}$.

The coefficient of the kinetic term of the dilaton field $\Phi$ is
$\frac{4}{3}$ in Eq. (6), which is not canonical. So we define a new dilaton
field $\phi$:

$\displaystyle\phi=\sqrt{\frac{8}{3}}\Phi.$ (7)

Now the action Eq. (6) becomes

$\displaystyle S^{E}=\int\mathrm{d}^{5}{\mathcal{L}}^{E},$ (8)

where ${\mathcal{L}}^{E}$ is the Lagrangian density in the Einstein frame:

$\displaystyle{\mathcal{L}}^{E}=\frac{1}{2\kappa_{5}^{2}}$
$\displaystyle\sqrt{-g^{E}}\Big{[}R^{E}-\frac{1}{2}{g^{E}}^{MN}\left(\partial_{M}\phi\right)\left(\partial_{N}\phi\right)$
$\displaystyle-
V_{\phi}(\phi)-\frac{h_{\phi}(\phi)}{4}{g^{E}}^{M\widetilde{M}}{g^{E}}^{N\widetilde{N}}F_{MN}F_{\widetilde{M}\widetilde{N}}\Big{]}.$
(9)

The function $V_{\phi}(\phi)=V^{E}(\Phi)$, and $h_{\phi}(\phi)=h(\Phi)$.
According to Eqs. (3), (5) and (7), we then derive the line element of the
space-time in the Einstein frame:

$\displaystyle
ds^{2}=\frac{L^{2}e^{2A_{E}(z)}}{z^{2}}\bigg{(}-f(z)dt^{2}+\frac{dz^{2}}{f(z)}+dy_{1}^{2}$
$\displaystyle\qquad\qquad\qquad\qquad+dy_{3}^{2}+dy_{3}^{2}\bigg{)},$ (10)

where

$\displaystyle A_{E}(z)=A_{s}(z)-\sqrt{\frac{1}{6}}\phi(z).$ (11)

Using the variation method, we can derive the Einstein field equations and the
equations of motion of $A_{M}$ and $\phi$ from Eq. (8)-(9):

$\displaystyle R_{MN}^{E}-\frac{1}{2}g_{MN}^{E}R^{E}-T_{MN}=0,$
$\displaystyle\nabla_{M}\left[h_{\phi}(\phi)F^{MN}\right]=0,$
$\displaystyle\partial_{M}\left[\sqrt{-g}\partial^{M}\phi\right]-\sqrt{-g}\left(\frac{\mathrm{d}V_{\phi}(\phi)}{\mathrm{d}\phi}+\frac{F^{2}}{4}\frac{\mathrm{d}h_{\phi}(\phi)}{\mathrm{d}\phi}\right)=0,$
(12)

where $T_{MN}$ is the energy-momentum tensor:

$\displaystyle T_{MN}$ $\displaystyle=$
$\displaystyle\frac{1}{2}\bigg{[}\left(\partial_{M}\phi\right)\left(\partial_{N}\phi\right)-\frac{1}{2}g_{MN}^{E}{g^{E}}^{P\widetilde{P}}\left(\partial_{P}\phi\right)\left(\partial_{\widetilde{P}}\phi\right)$
$\displaystyle-
g_{MN}^{E}V_{\phi}(\phi)\bigg{]}+\frac{h_{\phi}(\phi)}{2}\bigg{(}{g^{E}}^{P\widetilde{P}}F_{MP}F_{N\widetilde{P}}$
$\displaystyle-\frac{1}{4}g_{MN}^{E}{g^{E}}^{P\widetilde{P}}{g^{E}}^{Q\widetilde{Q}}F_{PQ}F_{\widetilde{P}\widetilde{Q}}\bigg{)}.$
(13)

It’s safe to suppose all the components of the vector field $A_{M}(z)$ are
zero except the $t$-component $A_{t}(z)$. After substituting the metric Eq.
(10) into the equations of motion (EoMs) Eq. (12), we then get the EoMs for
the components of the fields:

$\displaystyle
A_{t}^{\prime\prime}+A_{t}^{\prime}\left(-\frac{1}{z}+\frac{{h_{\phi}}^{\prime}}{h_{\phi}}+{A_{E}}^{\prime}\right)=0,$
(14) $\displaystyle
f^{\prime\prime}+f^{\prime}\left(-\frac{3}{z}+3{A_{E}}^{\prime}\right)-\frac{e^{-2{A_{E}}}A_{t}^{\prime
2}z^{2}h_{\phi}}{L^{2}}=0,$ (15) $\displaystyle
A_{E}^{\prime\prime}+\frac{f^{\prime\prime}}{6f}+A_{E}^{\prime}\left(-\frac{6}{z}+\frac{3f^{\prime}}{2f}\right)-\frac{1}{z}\left(-\frac{4}{z}+\frac{3f^{\prime}}{2f}\right)$
$\displaystyle+3{A_{E}}^{\prime
2}+\frac{L^{2}e^{2{A_{E}}}V_{\phi}}{3z^{2}f}=0,$ (16) $\displaystyle
A_{E}^{\prime\prime}-A_{E}^{\prime}\left(-\frac{2}{z}+A_{E}^{\prime}\right)+\frac{\phi^{\prime
2}}{6}=0,$ (17)
$\displaystyle\phi^{\prime\prime}+\phi^{\prime}\left(-\frac{3}{z}+\frac{f^{\prime}}{f}+3A_{E}^{\prime}\right)-\frac{L^{2}e^{2{A_{E}}}}{z^{2}f}\frac{\mathrm{d}V_{\phi}(\phi)}{\mathrm{d}\phi}$
$\displaystyle+\frac{z^{2}e^{-2{A_{E}}}A_{t}^{\prime
2}}{2L^{2}f}\frac{\mathrm{d}h_{\phi}(\phi)}{\mathrm{d}\phi}=0.$ (18)

Only 4 equations are independent in the above 5 equations. So we choose Eq.
(18) as a constraint and use it to check the solutions of the EoMs.

###  III.2 The matter part action in the holographic model

As mentioned earlier, $S_{m}^{s}$ in Eq. (1) is the action that describe
hadrons in the string frame. Here we use the improved KKSS action:

$\displaystyle
S_{m}^{s}=-\beta\int\mathrm{d}^{5}x\sqrt{-g^{s}}\mathrm{e}^{-\Phi}\operatorname{Tr}\Big{\\{}\left|D_{M}X\right|^{2}+V_{X}^{s}(X_{M})$
$\displaystyle\qquad\qquad\qquad\qquad\qquad\qquad\quad+\frac{1}{4g_{5}^{2}}\left(F_{L}^{2}+F_{R}^{2}\right)\Big{\\}}$
$\displaystyle\quad\qquad+S_{\text{baryons}}^{s},$ (19)

where

$\displaystyle D_{M}X=\nabla_{M}X-\mathrm{i}g_{c}A_{M}X$ (20)

is the covariant derivative of the $5$-dimensional scalar field $X$,

$\displaystyle
V_{X}^{s}(X_{M})={M_{X}}_{5}^{2}\left|X\right|^{2}+\lambda_{4}\left|X\right|^{4}$
(21)

is the potential of the $5$-dimensional scalar field $X$. The parameter
$\beta$ describe the coupling strength between the matter part and the
background in the action.
$D^{M}X=\partial^{M}X-\mathrm{i}A_{L}^{M}X+\mathrm{i}XA_{R}^{M}$,
$A_{L}^{M}=A_{L}^{a,M}t_{L}^{a}$, $A_{R}^{M}=A_{R}^{a,M}t_{R}^{a}$,
$t_{L}^{a}$ and $t_{R}^{a}$ are the generators of the $SU(2)_{L}$ and
$SU(2)_{R}$ group, respectively. $g_{5}^{2}=\frac{12\pi^{2}}{N_{c}}=4\pi^{2}$.
$S_{\text{baryons}}$ describes the baryons. The original KKSS action is
proposed in Ref. [24]. It is modified in Refs. [77, 78, 79, 80] to describe
the chiral phase transition and meson spectra.

According to the AdS/CFT dictionary [23], the bulk scalar field $X$ and the
chiral gauge fields $A_{L,R}^{M}$ are dual to the relevant QCD operators at
the ultraviolet (UV) boundary $z=0$. The bulk scalar field $X$ can be
decomposed as

$\displaystyle
X=\left(\frac{\chi(z)}{2}+S(x,z)\right)\mathrm{e}^{2\mathrm{i}\pi(x,z)},$ (22)

where $\pi(x,z)=\pi^{a}(x,z)t^{a}$ is the pseudoscalar meson field and
$S(x,z)$ is the scalar meson field. The field $\chi(z)$ is related to the
vacuum expectation value (VEV):

$\displaystyle\langle X\rangle=\frac{\chi}{2}I_{2},$ (23)

where $I_{2}$ is the $2\times 2$ identity matrix. The $5$-dimensional bulk
gauge field $A_{L,R}^{M}$ can be recombined into the vector field $V^{M}$ and
the axial-vector field $A^{M}$:

$\displaystyle V^{M}=\frac{1}{2}(A_{L}^{M}+A_{R}^{M}),$ $\displaystyle
A^{M}=\frac{1}{2}(A_{L}^{M}-A_{R}^{M}).$ (24)

The field-strength tensor for the vector field and the axial-vector field are

$\displaystyle F_{V}^{MN}=\frac{1}{2}\left(F_{L}^{MN}+F_{R}^{MN}\right)$
$\displaystyle\qquad=\partial^{M}V^{N}-\partial^{N}V^{M}-\mathrm{i}\left[V^{M},V^{N}\right]-\mathrm{i}\left[A^{M},A^{N}\right],$
(25) $\displaystyle F_{A}^{MN}=\frac{1}{2}\left(F_{L}^{MN}-F_{R}^{MN}\right)$
$\displaystyle\qquad=\partial^{M}A^{N}-\partial^{N}A^{M}-\mathrm{i}\left[V^{M},A^{N}\right]-\mathrm{i}\left[A^{M},V^{N}\right].$
(26)

According to mass-dimension relationship $M^{2}=(\Delta-p)(\Delta+p-4)$ and
$\Delta=3$, $p=0$, the $5$-dimensional mass square of the bulk field $X$ is
$M_{X}^{2}=-3$.

The action $S_{m}^{s}$ that describes the matter part can be decomposed as

$\displaystyle
S_{m}^{s}=S_{\chi}^{s}+S_{\text{mesons}}^{s}+S_{\text{baryons}}^{s},$ (27)

where

$\displaystyle
S_{\chi}=-\beta\int\mathrm{d}^{5}x\sqrt{-g^{s}}\mathrm{e}^{-\Phi}\bigg{\\{}\frac{1}{2}\left|\partial_{M}\chi-\mathrm{i}g_{c}A_{M}\chi\right|^{2}$
$\displaystyle\qquad\qquad\qquad\qquad\qquad\qquad+V_{\chi}^{s}(\chi)\bigg{\\}}$
(28)

is the action for the VEV $\chi$. The $5$-dimensional fields that dual to the
mesons and baryons tower are perturbations compared with the field $\chi$. The
contribution made by $S_{\text{mesons}}^{s}$ and $S_{\text{baryons}}^{s}$ to
the thermodynamics is believed to be small compared with that from
$S_{\chi}^{s}$. Thus we will neglect $S_{\text{mesons}}^{s}$ and
$S_{\text{baryons}}^{s}$ in the calculation in this article.

The EoM of the field $\chi(z)$ can be obtained from the action Eq. (19):

$\displaystyle\chi^{\prime\prime}(z)+\left(3{A_{E}}^{\prime}(z)+\frac{\sqrt{6}}{4}\phi^{\prime}(z)-\frac{3}{z}+\frac{f^{\prime}(z)}{f(z)}\right)\chi^{\prime}(z)$
$\displaystyle-\frac{1}{f(z)}\frac{1}{z^{2}}\mathrm{e}^{2A_{E}(z)+\frac{\sqrt{6}}{3}\phi(z)}$
$\displaystyle\times\bigg{[}M_{X}^{2}\chi(z)-g_{c}^{2}z^{2}\mathrm{e}^{-2A_{E}(z)-\frac{\sqrt{6}}{3}\phi(z)}\frac{1}{f(z)}{A_{t}(z)}^{2}\chi(z)$
$\displaystyle\qquad+\frac{1}{2}\lambda_{4}{\chi(z)}^{3}\bigg{]}=0.$ (29)

Here we use the Eq. (7) and Eq. (11).

The asymptotic expansion of the field $\chi(z)$ in the UV region is

$\displaystyle\chi(z\rightarrow 0)=m_{q}\zeta
z+\frac{\sigma}{\zeta}z^{3}+\cdots,$ (30)

where $m_{q}=3.5\mathrm{MeV}$ [100] is the current quark mass, $\sigma$ is the
chiral condensate, $\zeta=\frac{\sqrt{N_{c}}}{2\pi}=\frac{\sqrt{3}}{2\pi}$ is
the normalization constant.

###  III.3 Three different holographic models

In this work, we take three different holographic models.

####  III.3.1 Model I

Model I is a variant of Gubser model [101]. The potential for the dilation
field $\phi(z)$ and the function $h_{\phi}(\phi)$ are taken from Ref. [93]:

$\displaystyle
V_{\phi}(\phi)=-12\cosh(c_{1}\phi)+(6c_{1}^{2}-\frac{3}{2}){\phi}^{2}+c_{2}{\phi}^{6},$
(31) $\displaystyle
h_{\phi}(\phi)=\frac{1}{1+c_{3}}\text{sech}(c_{4}{\phi}^{3})+\frac{c_{3}}{1+c_{3}}{\mathrm{e}}^{-c_{5}\phi},$
(32)

where $c_{1}=0.71$, $c_{2}=0.0037$, $c_{3}=1.935$, $c_{4}=0.085$, and
$c_{5}=30$. The $5$-dimensional Newtonian constant is $G_{5}=0.42$. The
characteristic energy scale [102] of the EMD system is
$\Lambda=1.085\mathrm{GeV}$. The Eqs. (31)-(32) and the values of the above
parameters are taken from Ref. [93], which can give the thermodynamic
properties in good agreement with the lattice simulations for the system of
$2+1+1$ dynamical quarks with physical masses [103].

Besides those parameters, $\beta=0.834954$, $g_{c}=0.3$, $\lambda_{4}=16$. The
critical value of baryon chemical potential is $\mu_{P0}=977\mathrm{MeV}$. The
meaning of this parameter will be interpreted in Subsec. IV.1.

####  III.3.2 Model II

In model II, the potential for the dilation field $\phi(z)$ is

$\displaystyle
V_{\phi}(\phi)=-12-\left(\frac{3}{2}-\frac{2}{3}c_{1}\right){\phi}^{2}$
$\displaystyle\qquad\qquad-4c_{1}(1+{\phi}^{2})^{\frac{c_{2}}{2}}\left[\sinh(\frac{\phi}{\sqrt{6}})\right]^{2},$
(33)

and the function $h_{\phi}(\phi)$ is

$\displaystyle h_{\phi}(\phi)=\text{sech}(c_{3}\phi),$ (34)

where $c_{1}=1$, $c_{2}=\frac{1}{2}$, $c_{3}=0.7$. The function
$V_{\phi}(\phi)$ are taken from Ref. [104]. The $5$-dimensional Newtonian
constant is $G_{5}=1$. The characteristic energy scale [102] of the EMD system
is $\Lambda=1\mathrm{GeV}$.

Besides those parameters, $\beta=1$, $g_{c}=0$, $\lambda_{4}=16$. The critical
value of baryon chemical potential is $\mu_{P0}=1085\mathrm{MeV}$.

####  III.3.3 Model III

In model III, the dilation field $\phi(z)$ is

$\displaystyle\phi(z)=c_{1}z^{2},$ (35)

and the function $h_{\phi}$ is

$\displaystyle h_{\phi}(z)=\mathrm{e}^{c_{2}z^{2}-A_{E}(z)},$ (36)

where $c_{1}=1.536{\mathrm{GeV}}^{2}$ and $c_{2}=1.16{\mathrm{GeV}}^{2}$. The
$5$-dimensional Newtonian constant is $G_{5}=17.930$. The value of the
parameter $c_{1}$ is taken from Ref. [96], which can produce good results for
the glueballs/oddballs spectra compared with those of the lattice QCD
simulation, as “Model III, IV(1)” in Table III in Ref. [96]. The form of the
function $h_{\phi}(z)$ is inspired by Refs. [105, 106].

Besides those parameters, $\beta=0.7808$, $g_{c}=0$, $\lambda_{4}=16$. The
critical value of baryon chemical potential is $\mu_{P0}=977\mathrm{MeV}$.

###  III.4  $T-\mu$ phase diagram for the confinement-deconfinement phase
transition

We can get the black-hole solution from the EoMs Eq. (14)-(18). The
temperature of the black-hole is

$\displaystyle T=\left|\frac{f^{\prime}(z_{h})}{4\pi}\right|,$ (37)

where $z_{h}$ is the location of the horizon for the black-hole. The entropy
density of the black-hole is

$\displaystyle
s_{b}=\frac{\mathrm{e}^{3A_{E}(z_{h})}}{\frac{\kappa_{5}^{2}}{2\pi}z_{h}^{3}}.$
(38)

According to the AdS/CFT dictionary, the temperature and the entropy density
of the QCD matter are $T$ and $s_{b}$, respectively.

According to the AdS/CFT dictionary, the chemical potential is

$\displaystyle\mu=A_{t}(z=0),$ (39)

and the baryon number density is

$\displaystyle n_{b}=\left|\lim_{z\rightarrow
0}{\frac{\partial{{\mathcal{L}}^{E}}}{\partial\left(\partial_{z}A_{t}\right)}}\right|$
$\displaystyle\quad=-\frac{1}{2{\kappa_{5}}^{2}}\lim_{z\rightarrow
0}\left[\frac{{\mathrm{e}}^{A_{E}(z)}}{z}h_{\phi}(\phi)\frac{\mathrm{d}}{\mathrm{d}z}A_{t}(z)\right],$
(40)

where ${\mathcal{L}}^{E}$ is the Lagrangian density in the Einstein frame in
Eq. (9). From the Euler-Lagrange equation of the field $A_{t}(z)$:

$\displaystyle\partial_{M}{\frac{\partial{{\mathcal{L}}^{E}}}{\partial\left(\partial_{M}A_{t}\right)}}-\frac{\partial{{\mathcal{L}}^{E}}}{\partial
A_{t}}=0,$ (41)

we then derive that

$\displaystyle\partial_{z}\left[\frac{{\mathrm{e}}^{A_{E}(z)}}{z}h_{\phi}(\phi)\frac{\mathrm{d}}{\mathrm{d}z}A_{t}(z)\right]=0.$
(42)

Eq. (42) is actually the EoM of $A_{t}$ Eq. (14). Thus, we get the conserved
Gauss charge associated with the field $A_{t}$:

$\displaystyle
Q_{G}=\frac{{\mathrm{e}}^{A_{E}(z)}}{z}h_{\phi}(\phi)\frac{\mathrm{d}}{\mathrm{d}z}A_{t}(z).$
(43)

Then we derive

$\displaystyle n_{b}=-\frac{1}{2{\kappa_{5}}^{2}}Q_{G}.$ (44)

Sometimes the following expression of $n_{b}$ which is related to the
quantities at the black-hole horizon is more convenient in the calculation:

$\displaystyle
n_{b}=-\frac{1}{2\kappa_{5}^{2}}\frac{\mathrm{e}^{A_{E}(z_{h})}}{z_{h}}h_{\phi}(\phi=\phi(z_{h})){A_{t}}^{\prime}(z_{h}).$
(45)

The internal energy density and the free energy density are

$\displaystyle\epsilon_{b}=Ts_{b}-P_{b}+\mu n_{b},$
$\displaystyle\mathcal{F}_{b}=-P_{b}=\epsilon_{b}-Ts_{b}-\mu n_{b}.$ (46)

The differential relation of the above equations are

$\displaystyle\mathrm{d}\epsilon_{b}=T\mathrm{d}s_{b}+\mu\mathrm{d}n_{b},$
$\displaystyle\mathrm{d}\mathcal{F}_{b}=-\mathrm{d}P_{b}=-s_{b}\mathrm{d}T-n_{b}\mathrm{d}\mu.$
(47)

From Eq. (47), we can calculate the pressure $P$ and the internal energy
density $\epsilon$. The trace anomaly is

$\displaystyle I_{b}(T,\mu)=\epsilon_{b}(T,\mu)-3P_{b}(T,\mu)$
$\displaystyle\qquad\qquad\quad=T{s_{b}}(T,\mu)+\mu{n_{b}}(T,\mu)-4P_{b}(T,\mu).$
(48)

The square of the speed of sound is

$\displaystyle{{c_{s}}_{b}}^{2}=\frac{\mathrm{d}P_{b}}{\mathrm{d}{\epsilon_{b}}}.$
(49)

Here all the thermodynamic quantities are labeled with “b”, which denotes
“background”. The corresponding phase transition is the confinement-
deconfinement phase transition.

Figure 1: The $T-\mu$ phase diagram of the confinement-deconfinement phase
transition in models I, II, and III. The $T_{c}$ is the critical temperature
of the confinement-deconfinement phase transition, $\mu$ is the baryon
chemical potential. The blue full line represents that the confinement-
deconfinement phase transition is first order transition, and the tan dashed
line represents that it’s the cross-over behavior. In model I and II, the
critical temperature of the cross-over behavior is defined as the temperature
where the value of the speed of sound is minimal. In model III, the critical
temperature of the cross-over behavior is defined as the temperature where
$\lvert\frac{\partial}{\partial T}{z_{h}}(T,\,\mu)\rvert$ is maximal. The red
point is the critical end point (CEP), the position of which is
${(T,\,\mu)}_{\mathrm{CEP}}$. In model I,
${(T,\,\mu)}_{\mathrm{CEP}}=(105.032\mathrm{MeV},\,555.52\mathrm{MeV})$. In
model II,
${(T,\,\mu)}_{\mathrm{CEP}}=(101.478\mathrm{MeV},\,474.5\mathrm{MeV})$. In
model III,
${(T,\,\mu)}_{\mathrm{CEP}}=(283.359\mathrm{MeV},\,525\mathrm{MeV})$.

The $T-\mu$ phase diagrams for the confinement-deconfinement phase transition
in model I, II, and III are shown in in Fig. 1. $T_{c}$ is the critical
temperature of the confinement-deconfinement phase transition, $\mu$ is the
baryon chemical potential. The blue full line represents that the confinement-
deconfinement phase transition is first order transition, and the tan dashed
line represents that it’s the cross-over behavior.

It is noticed that the confinement-deconfinement phase transition is the
cross-over behavior at low baryon chemical potential range and a first order
phase transition at high baryon chemical potential range in model I. The CEP
is located at
${(T,\,\mu)}_{\mathrm{CEP}}=(105.032\mathrm{MeV},\,555.52\mathrm{MeV})$, which
is close to the FRG result in Ref. [107] and that of the rPNJL model [108].
The confinement-deconfinement phase transition is a first order phase
transition at low baryon chemical potential range and the cross-over behavior
at high baryon chemical potential range in model II and model III, which
contrasts with that in model I. The CEP is located at
${(T,\,\mu)}_{\mathrm{CEP}}=(101.478\mathrm{MeV},\,474.5\mathrm{MeV})$ in
model II and
${(T,\,\mu)}_{\mathrm{CEP}}=(283.359\mathrm{MeV},\,525\mathrm{MeV})$ in model
III.

When $\mu=0$, the critical temperature is $142.227\mathrm{MeV}$,
$146.803\mathrm{MeV}$, and $343.466\mathrm{MeV}$ in model I, II and III,
respectively. When $T=0$, the critical baryon chemical potential is about
$1.22\mathrm{GeV}$, $0.81\mathrm{GeV}$, and $2.765\mathrm{GeV}$ in model I, II
and III, respectively. In Subsec. IV.1 we will see that these critical
chemical potentials will be changed when we consider the contribution made by
the matter part action.

###  III.5 Metric deformation in the three models

To understand the system described by model I, model II, and model III, before
we calculate the EoS of the quark matter, we compare the backgrounds in the
three models with the ${\mathrm{AdS}}_{5}$ space-time.

Figure 2: The function $\frac{{\mathrm{e}}^{A_{s}(z)}}{z}$ in the deformed
metric in three models at $\mu=1.2\mathrm{GeV}$, $\mu=1.4\mathrm{GeV}$,
$\mu=1.6\mathrm{GeV}$, and $\mu=1.8\mathrm{GeV}$ at $T=0$. The steel blue
line, the goldenrod line, and the red line represent the function
$\frac{{\mathrm{e}}^{A_{s}(z)}}{z}$ in model I, model II, and model III,
respectively. The black line is that in the ${\mathrm{AdS}}_{5}$ space-time,
which is $\frac{1}{z}$.

In Fig. 2, we plot the function $\frac{{\mathrm{e}}^{A_{s}(z)}}{z}$ in the
deformed metric in three models at four different values of baryon chemical
potential: $\mu=1.2\mathrm{GeV}$, $\mu=1.4\mathrm{GeV}$,
$\mu=1.6\mathrm{GeV}$, and $\mu=1.8\mathrm{GeV}$ at $T=0$. In the $AdS_{5}$
space-time, $\frac{{\mathrm{e}}^{A_{s}(z)}}{z}=\frac{1}{z}$, which is shown by
the black solid line. So the deviation of the lines in Fig. 2 from the black
lines shows the deformation of the metric in the three holographic models.

We can see from Fig. 2 that the function $\frac{{\mathrm{e}}^{A_{s}(z)}}{z}$
in the three models are all located above the $AdS_{5}$ case. As mentioned in
Subsec. III.3, under the parameter values in model I, the EMD system can
derive the thermodynamics that is consistent with the lattice QCD result of
the system of $2+1+1$ dynamical quarks with physical masses [93, 103]. Thus,
we think the background of model I describes the $2+1+1$ flavor system. The
parameters in model III can produce good results of the glueballs/oddballs
spectra compared with those of the lattice QCD simulation [96]. Thus, the
background of model III is considered to describe the pure gluon system. Fig.
2 shows that the pure gluon system is most deformed in the infrared (IR)
range.

In Subsec. V.1 we will see that the astrophysical observations of neutron star
prefer model I and model III than model II. So we surmise that the background,
the function $\frac{{\mathrm{e}}^{A_{s}(z)}}{z}$ of which is located between
the steel blue line and the red line in Fig. 2, can produce neutron star
properties that are consistent with the astrophysical observations.

##  IV Equation of state for the cold QCD matter

###  IV.1 Nuclear matter

In this subsection, we introduce how to describe nuclear matter in the
holographic model. The total pressure consists of the contribution made by the
background part and that made by the matter part:

$\displaystyle P_{\mathrm{tot}}=P_{b}+P_{m}.$ (50)

As discussed in Subsec. III.4, the background part of the pressure $P_{b}$ can
be calculated from the confined solutions of the EoMs Eq. (14)-(18). But the
matter part of the pressure $P_{\text{m}}$ is obtained through calculating the
on-shell Euclidean action ${S_{m}}_{\text{Euclidean}}^{\text{on-shell}}$:

$\displaystyle P_{m}=-{\mathcal{F}}_{m}=-{S_{m}}_{\text{Euclidean}}^{\text{on-
shell}}.$ (51)

When we calculate the on-shell Euclidean action for the matter part, we
neglect the contribution from the $S_{\text{mesons}}$ and
$S_{\text{baryons}}$:

$\displaystyle{S_{m}}_{\text{Euclidean}}^{\text{on-
shell}}\simeq{S_{\chi}}_{\text{Euclidean}}^{\text{on-shell}}.$ (52)

As mentioned in Subsec. III.2, the pressure density from $5$-dimensional
fields that dual to the mesons and baryons tower are perturbations compared
with that from the field $\chi$. So we think the contribution to the Euclidean
on-shell action from $S_{\text{mesons}}$ and $S_{\text{baryons}}$ are
negligible.

Now we use model I as the example to illustrate how to describe nuclear matter
and then to calculate the EoS of nuclear matter. From Fig. 1, we can see that
there is a intersection point between the $T_{c}-\mu$ curve and the horizontal
axis, the position of which is denoted as $(T=0,\,\mu_{c}(T=0)\simeq
1.22\mathrm{GeV})$. When $T=0$ and $\mu<\mu_{c}(T=0)$, the background is the
thermal-gas solutions, the line element for which in Einstein frame is

$\displaystyle
ds^{2}=\frac{L^{2}e^{2A_{E}(z)}}{z^{2}}\left(-dt^{2}+dz^{2}+dy_{1}^{2}+dy_{3}^{2}+dy_{3}^{2}\right).$
(53)

In this case, the values of $s_{b}$, $n_{b}$, $P_{b}$, $\epsilon_{b}$ are all
$0$. Thus, from Eqs. (50)-(52) the total pressure is

$\displaystyle P_{\mathrm{tot}}\simeq-{S_{\chi}}_{\text{Euclidean}}^{\text{on-
shell}}.$ (54)

When $T=0$ and $\mu>\mu_{c}(T=0)$, although the thermal-gas solutions also
exist, the black-hole solutions are more stable, and the corresponding line
element in the Einstein frame is Eq. (10).

In principle, the calculation should be implemented in the case of $T=0$.
However, it’s very difficult to calculate the solutions in strict zero
temperature. Thus, we calculate the EoS of cold strong interaction matter in a
very low temperature, such as $T=0.1085\mathrm{MeV}$ for quark matter. When we
calculate the total pressure $P_{\mathrm{tot}}$, we choose a critical value of
baryon chemical potential, $\mu_{P0}=977\mathrm{MeV}$, so that
$P_{\mathrm{tot}}(\mu=\mu_{P0})=0$. This value is slightly different from the
$939\mathrm{MeV}$ in Sec. “Baryons from a homogeneous bulk gauge field” in
Ref. [60].

Figure 3: The pressure $P$ of the cold QCD matter as a function of the baryon
chemical potential $\mu$. The olive green line, the silver lake blue line, and
the red line represent the HLPS intermediate EoS [109], the holographic result
which contains contribution made by the background part and that made by the
matter part, and the holographic result which contains only the contribution
made by the background part, respectively. The full lines are the
thermodynamic stable parts, and the dashed lines are the unstable parts.
Additionally, we also plot the dotted purple line which represents the
contribution made by the matter part.

We then plot the function $P_{\mathrm{tot}}$ as the silver lake blue line in
Fig. 3. There are also other three different colors in Fig. 3: the olive green
line represents the HLPS intermediate case [109], the red line, which is
calculated in Subsec. III.4, represents the holographic result which contains
only the contribution made by the background part, and the dotted purple line
represents the contribution made by the matter part. There is an intersection
between the olive green line and the silver lake blue line, where the baryon
chemical potential is $1.084\mathrm{GeV}$. The intersection represents a
second-order phase transition, where $P$ and $\frac{\partial
P}{\partial_{\mu}}\rvert_{T=0}$ is continuous, but
$\frac{\partial^{2}P}{\partial_{\mu}^{2}}\rvert_{T=0}$ are different at both
side of this critical values of $\mu$.

Similarly, the intersection between the silver lake blue line and the red
line, where $\mu=1.723\mathrm{GeV}$, represents a first order transition. The
value of $\frac{\partial P}{\partial_{\mu}}\rvert_{T=0}$ are different at both
side of this critical values of $\mu$. We regard this phase transition as the
confinement-deconfinement phase transition in QCD. It is noticed that after
considering the contribution made by the matter part action, the critical
chemical potential $\mu_{c}=1.723\mathrm{GeV}$ is different from that in
Subsec. III.4. Thus, the silver lake blue line is the EoS for the holographic
nuclear matter and the red line is the EoS for the holographic quark matter.

We note that the silver lake blue line is full line only when
$1.084\mathrm{GeV}<\mu<1.723\mathrm{GeV}$. It means the free energy, in this
case, is minimal between the three cases (the olive green line, the silver
lake blue line, and the red line) and it is thermodynamic stable when
$1.084\mathrm{GeV}<\mu<1.723\mathrm{GeV}$. When $\mu<1.084\mathrm{GeV}$, the
nuclear matter with the HLPS intermediate EoS is stable. When
$\mu>1.723\mathrm{GeV}$, the holographic quark matter is stable.

After we get the pressure of the nuclear matter $P_{NM}=P_{\mathrm{tot}}$, we
can calculate the entropy and the baryon number density of the nuclear matter
using the differential relation:

$\displaystyle\mathrm{d}P_{NM}(T,\mu)=s_{NM}\mathrm{d}T+n_{NM}\mathrm{d}\mu.$
(55)

Then we can derive all other thermodynamic quantities of the nuclear matter.
That’s why all the thermodynamic quantities here are labeled with “NM”, which
denotes “nuclear matter”.

###  IV.2 Hybrid equation of state for neutron stars

We show the final results of the hybrid EoS for cold QCD matter from model I,
II and III in Fig. 4, and the corresponding square of the speed of sound in
Fig. 5 and Fig. 6.

Figure 4:  The hybrid EoS for the cold QCD matter in the pressure-internal
energy density plane calculated in model I, model II, and model III. The green
band, which is taken from Ref. [43], is the “allowed” range of the EOS for
cold QCD matter. It is defined by the low-density chiral effective theory, the
high-density perturbative QCD, and the polytropic interpolations between them,
and is constrained by the astrophysics observations. In model I, there are
three parts in the hybrid EoS. The olive green, the silver lake blue, and the
red segment represent the HLPS intermediate EoS [109], holographic nuclear
matter EoS, and holographic quark matter EoS, respectively. In model II, there
are two parts in the hybrid EoS. The olive green and the red segment represent
the HLPS intermediate EoS [109] and holographic quark matter EoS,
respectively. In model III, there are two parts in the hybrid EoS. The olive
green and the silver lake blue segment represent the HLPS intermediate EoS
[109], and holographic nuclear matter EoS, respectively.

Figure 5:  The square of the speed of sound $c_{s}^{2}$ as a function of the
baryon chemical potential $\mu$ for the cold QCD matter calculated in model I,
model II, and model III. In model I, there are three parts in the panel. The
olive green, the silver lake blue, and the red segment correspond to the HLPS
intermediate EoS [109], holographic nuclear matter EoS, and holographic quark
matter EoS, respectively. In model II, there are two parts in the panel. The
olive green and the red segment correspond to the HLPS intermediate EoS [109],
and holographic quark matter EoS, respectively. In model III, there are two
parts in the panel. The olive green and the silver lake blue segment
correspond to the HLPS intermediate EoS [109], and holographic nuclear matter
EoS, respectively.

Figure 6:  The square of the speed of sound $c_{s}^{2}$ as a function of the
scaled baryon number density $\frac{n_{B}}{n_{s}}$ for the cold QCD matter
calculated in model I, model II, and model III. In model I, there are three
parts in the panel. The olive green, the silver lake blue, and the red segment
correspond to the HLPS intermediate EoS [109], holographic nuclear matter EoS,
and holographic quark matter EoS, respectively. In model II, there are two
parts in the panel. The olive green and the red segment correspond to the HLPS
intermediate EoS [109], and holographic quark matter EoS, respectively. In
model III, there are two parts in the panel. The olive green and the silver
lake blue segment correspond to the HLPS intermediate EoS [109], and
holographic nuclear matter EoS, respectively.

The hybrid EoS of model I includes three parts, i.e., the HLPS intermediate
EoS [109], the holographic nuclear matter EoS, and the holographic quark
matter EoS. The hybrid EoS of model II includes the HLPS intermediate EoS and
the holographic quark matter EoS. The hybrid EoS of model III includes the
HLPS intermediate EoS and the holographic nuclear matter EoS.

It is observed that the EoSs in the three models all fall in the green band,
which is taken from Ref. [43], is the “allowed” range of the EOS for cold QCD
matter. It is defined by the low-density chiral effective theory, the high-
density perturbative QCD, and the polytropic interpolations between them, and
is constrained by the astrophysics observations.

In Fig. 5, the corresponding square of the speed of sound is shown as a
function of the chemical potential. We can observe that the square of the
speed of sound increases from a small value at zero chemical potential to a
relatively large value more than $0.8$, then decreases to a smaller value at
the high baryon chemical potential range. The violation of sound velocity
bound in neutron stars has been discussed in much literature, e.g. Refs. [110,
111]. The sound velocity in model I and model III reaches maximum in
holographic nuclear matter. However, the maximum sound velocity in model II
lies in the HLPS intermediate part.

In Fig. 6, we show the square of the speed of sound as a function of the
scaled baryon number density $\frac{n_{B}}{n_{s}}$, where $n_{s}$ is the
nuclear saturation baryon number density.

From Figs. 4-6, we clearly observe that in model I, there are two phase
transition. The first one is a second order transition. At this phase
transition, the baryon number density $n_{B}=2.221n_{s}$, the baryon chemical
potential $\mu=1.084\mathrm{GeV}$, the internal energy density
$\epsilon=349\mathrm{MeV}/{{\mathrm{fm}}^{3}}$, and the pressure
$P=35.962\mathrm{MeV}/{{\mathrm{fm}}^{3}}$. This phase transition connects the
HLPS intermediate nuclear matter with the holographic nuclear matter.

The HLPS intermediate nuclear matter and the holographic nuclear matter may be
in the same QCD deconfinement phase. So they should be joined together
smoothly. However, in practice calculation, the two EoSs come from different
models and we have tried our best to use a second order phase transition to
join them together. The square of the speed of sound, $c_{s}^{2}$, jumps from
$0.373$ to $0.854$ through this second order transition. Another possible
interpretation for the holographic nuclear matter is that it’s in the
quarkyonic phase [112, 90] which is in confinement with chiral symmetry
restoration. We don’t investigate the quarkyonic phase further here and leave
it in future work.

The second one is a first order transition, which is interpreted as the
confinement-deconfinement phase transition in Subsec. IV.1 and is located at
$\mu=1.723\mathrm{GeV}$ and $P=349.276\mathrm{MeV}/{{\mathrm{fm}}^{3}}$. The
value of $n_{B}$ jumps from $3.999n_{s}$ to $7.328n_{s}$, the value of
$c_{s}^{2}$ jumps from $0.710$ to $0.288$, and the value of $\epsilon$ jumps
from $753.260\mathrm{MeV}/{{\mathrm{fm}}^{3}}$ to
$1671.826\mathrm{MeV}/{{\mathrm{fm}}^{3}}$ through this first order
transition. At $\mu=1.890\mathrm{GeV}$, the baryon number density $n_{B}$
reaches $10n_{s}$, which is a very high value, and $c_{s}^{2}=0.304$. We will
see in Subsec. V.1 that this first order transition leads to the branch that
corresponds to the holographic quark matter is unstable in the mass-radius
figure of the NS.

In model II, there is a first order phase transition, which is interpreted as
the confinement-deconfinement phase transition and located at
$\mu=1.781\mathrm{GeV}$ and $P=427.250\mathrm{MeV}/{{\mathrm{fm}}^{3}}$. The
value of $n_{B}$ jumps from $4.602n_{s}$ to $8.163n_{s}$, the value of
$c_{s}^{2}$ jumps from $0.815$ to $0.283$, and the value of $\epsilon$ jumps
from $884.477\mathrm{MeV}/{{\mathrm{fm}}^{3}}$ to
$1899.401\mathrm{MeV}/{{\mathrm{fm}}^{3}}$ through this first order
transition. At $\mu=1.888\mathrm{GeV}$, the baryon number density $n_{B}$
reaches $10n_{s}$ and $c_{s}^{2}=0.290$. Being similar to the case in model I,
this first order phase transition here also leads to the unstable branch that
corresponds to the holographic quark matter.

In model III, there is a second order phase transition. The parameter values
of this phase transition are $n_{B}=2.223n_{s}$, $\mu=1.085\mathrm{GeV}$,
$\epsilon=349.7\mathrm{MeV}/{{\mathrm{fm}}^{3}}$, and
$P=36.083\mathrm{MeV}/{{\mathrm{fm}}^{3}}$. The square of the speed of sound,
$c_{s}^{2}$, jumps from $0.374$ to $0.842$ through this second order
transition. At $\mu=2.645\mathrm{GeV}$, the baryon number density $n_{B}$
reaches $10n_{s}$ and $c_{s}^{2}=0.399$.

##  V The neutron star

###  V.1 Mass-radius relation

The Tolman-Oppenheimer-Volkoff (TOV) equation [113, 114, 115] describes the
structure of a spherically symmetric and isotropic object in static
gravitational equilibrium:

$\displaystyle\frac{\mathrm{d}}{\mathrm{d}r}P(r)+\left[P(r)+\epsilon(r)\right]\frac{G\left[m(r)+4\pi
r^{3}P(r)\right]}{r^{2}\left[1-2\frac{Gm(r)}{r}\right]}=0,$ (56)

where $r$ is the radial coordinate, $P(r)$ is the pressure, and $\epsilon(r)$
is the internal energy density of the material at the radius $r$. $G$ is the
Newtonian constant. $m(r)$ is the total mass inside the sphere with radius $r$
measured by the gravitational field felt by a distant observer. It satisfies

$\displaystyle\frac{\mathrm{d}}{\mathrm{d}r}m(r)-4\pi r^{2}\epsilon(r)=0,$
(57) $\displaystyle m(r=0)=0.$ (58)

The metric inside the object can be obtained from the element:

$\displaystyle{\mathrm{d}s_{\mathrm{in}}}^{2}=\mathrm{e}^{\nu(r)}{\mathrm{d}t}^{2}-\frac{{\mathrm{d}r}^{2}}{1-2\frac{Gm(r)}{r}}$
$\displaystyle\qquad\quad-r^{2}\left[{\mathrm{d}\theta_{1}}^{2}+\left(\sin{\theta_{1}}\right)^{2}{\mathrm{d}\theta_{2}}^{2}\right],$
(59)

where $\nu(r)$ obeys the following equation:

$\displaystyle\frac{\mathrm{d}}{\mathrm{d}r}\nu(r)=-\left(\frac{2}{P(r)+\epsilon(r)}\right)\frac{\mathrm{d}}{\mathrm{d}r}P(r).$
(60)

If the sphere is in the vacuum, the metric outside it is the Schwarzschild
metric:

$\displaystyle{\mathrm{d}s_{\mathrm{out}}}^{2}=\left[1-2\frac{GM}{r}\right]{\mathrm{d}t}^{2}-\frac{{\mathrm{d}r}^{2}}{1-2\frac{GM}{r}}$
$\displaystyle\qquad\quad-r^{2}\left[{\mathrm{d}\theta_{1}}^{2}+\left(\sin{\theta_{1}}\right)^{2}{\mathrm{d}\theta_{2}}^{2}\right],$
(61)

where $M=m(r=r_{B})$ is the total mass of the object measured by the
gravitational field felt by a distant observer, $r_{B}$ is the radius
coordinate of the boundary of the object.

To make sure the metric is continuous at the boundary $r=r_{B}$, we impose the
boundary condition

$\displaystyle\mathrm{e}^{\nu(r=r_{B})}=1-2\frac{GM}{r}.$ (62)

Two other boundary conditions are the zero pressure condition at the boundary
of the object

$\displaystyle P(r=r_{B})=0,$ (63)

and the given internal energy density at $r=0$

$\displaystyle\epsilon(r=0)=\epsilon_{0}.$ (64)

The NS consists of cold QCD matter. We only consider the spherical symmetric
and isotropic neutron star in the vacuum without rotation. The EoS of the cold
QCD matter can be written in the form

$\displaystyle\mathcal{C}\left(P(r),\epsilon(r)\right)=0.$ (65)

Figure 7:  The mass-radius relation of the neutron star calculated in model I,
model II, and model III. In model I, there are three parts in the M-R curve.
The olive green, the silver lake blue, and the red segment correspond to the
HLPS intermediate EoS [109], holographic nuclear matter, and holographic quark
matter, respectively. In model II, there are two parts in the M-R curve. The
olive green and the red segment correspond to the HLPS intermediate EoS [109],
and holographic quark matter, respectively. In model III, there are two parts
in the M-R curve. The olive green and the silver lake blue segment correspond
to the HLPS intermediate EoS [109], and holographic nuclear matter,
respectively.

From Fig 7, we observe that in model I, the radius for the $1.4M_{\odot}$-NS
is $12.320\mathrm{km}$. The max possible radius for the NS is
$12.833\mathrm{km}$, which correspond a $2.32M_{\odot}$-NS. The max possible
mass for the NS is $2.630M_{\odot}$, which correspond a
$12.559\mathrm{km}$-NS.

In model II, the radius for the $1.4M_{\odot}$-NS is $12.287\mathrm{km}$. The
max possible radius for the NS is $12.475\mathrm{km}$, which correspond a
$1.98M_{\odot}$-NS. The max possible mass for the NS is $2.466M_{\odot}$,
which correspond a $11.902\mathrm{km}$-NS.

In model III, the radius for the $1.4M_{\odot}$-NS is $12.319\mathrm{km}$. The
max possible radius for the NS is $12.795\mathrm{km}$, which correspond a
$2.27M_{\odot}$-NS. The max possible mass for the NS is $2.598M_{\odot}$,
which correspond a $12.24\mathrm{km}$-NS.

In conclusion, the max possible mass for the neutron star is about
$2.5M_{\odot}$ and the radius is about $12\mathrm{km}$ then. It is noticed
that the branch that corresponds to the holographic quark matter is always
unstable, and the holographic nuclear matter can produce the stable branch.
These results indicate that even in the core of the NS, the matter is still in
the confinement phase and the quark matter is not favored. This is consistent
with the conclusion in Ref. [116], in which the authors combined the input
from the gauge-gravity duality and the input from various $ab$-$initio$
methods for the nuclear matter at low density and analyzed families of hybrid
equations of state of cold QCD matter. Their prediction is that all neutron
stars are fully hadronic without quark matter cores.

###  V.2 Tidal deformability

At last, we investigate the tidal deformability of the NSs. We consider a
static, spherically symmetric object in a static external quadrupolar tidal
field ${\mathcal{E}}_{ij}$ [117, 118, 119]. Then the tidal deformability
$\lambda$ is defined to linear order as

$\displaystyle Q_{ij}=-\lambda{\mathcal{E}}_{ij},$ (66)

where $Q_{ij}$ is the quadrupole moment. The dimensional tidal deformability
is also related to the $l=2$ dimensionless tidal Love number $k_{2}$:

$\displaystyle k_{2}=\frac{3}{2}\lambda GR^{-5}.$ (67)

The asymptotic expansion of the $tt$ component of the full metric outside the
object at a large distance $r$ from the object is

$\displaystyle{g_{\mathrm{out}}}_{tt}=1-2\frac{GM}{r}-\frac{3GQ_{ij}}{r^{3}}\frac{x^{i}}{r}\frac{x^{j}}{r}+\cdots$
$\displaystyle\qquad\quad+{\mathcal{E}}_{ij}r^{2}\frac{x^{i}}{r}\frac{x^{j}}{r}+\cdots.$
(68)

The right hand side (r.h.s.) of Eq. (68) includes two leading-order terms of
the perturbative expansions. The second one describes an external tidal field
increasing with $r^{2}$ and the first one describes the resulting tidal
distortion decreasing with $r^{-3}$.

Introducing a linear $l=2$ perturbation on the object leads to a static, even-
parity perturbation on the metric [118, 119] in the Regge-Wheeler gauge [120]:

$\displaystyle{\mathrm{d}s_{\mathrm{perturb}}}^{2}={\mathrm{e}}^{\nu(r)}\left[1+H(r)Y_{20}(\theta_{1},\theta_{2})\right]{\mathrm{d}t}^{2}$
$\displaystyle\qquad\quad-\frac{\left[1-H(r)Y_{20}(\theta_{1},\theta_{2})\right]}{1-2\frac{Gm(r)}{r}}{\mathrm{d}r}^{2}$
$\displaystyle\qquad\quad-r^{2}\left[1-K(r)Y_{20}(\theta_{1},\theta_{2})\right]\left[{\mathrm{d}\theta_{1}}^{2}+\left(\sin{\theta_{1}}\right)^{2}{\mathrm{d}\theta_{2}}^{2}\right],$
(69)

where

$\displaystyle\frac{\mathrm{d}}{\mathrm{d}r}K(r)=\frac{\mathrm{d}}{\mathrm{d}r}H(r)+H(r)\frac{\mathrm{d}}{\mathrm{d}r}\nu(r).$
(70)

The corresponding perturbations of the stress-energy tensor of the perfect
fluid are

$\displaystyle\delta
T_{0}^{0}=-\delta\epsilon(r)Y_{20}(\theta_{1},\theta_{2})$ (71)
$\displaystyle\delta T_{i}^{i}=\delta P(r)Y_{20}(\theta_{1},\theta_{2}).$ (72)

The function $H(r)$ obeys the following differential equation:

$\displaystyle\Bigg{\\{}-\frac{6}{r^{2}\left(1-2\frac{Gm(r)}{r}\right)}-\frac{1}{2}\left[\frac{\mathrm{d}}{\mathrm{d}r}\nu(r)\right]^{2}+\frac{{\mathrm{d}}^{2}}{\mathrm{d}r^{2}}\nu(r)$
$\displaystyle+\frac{3}{r}\frac{G\left[r\frac{\mathrm{d}}{\mathrm{d}r}m(r)-m(r)\right]}{r^{2}\left[1-2\frac{Gm(r)}{r}\right]}+\frac{7}{2r}\frac{\mathrm{d}}{\mathrm{d}r}\nu(r)$
$\displaystyle-\frac{\mathrm{d}}{\mathrm{d}r}\nu(r)\frac{G\left[r\frac{\mathrm{d}}{\mathrm{d}r}m(r)-m(r)\right]}{r^{2}\left[1-2\frac{Gm(r)}{r}\right]}$
$\displaystyle+\frac{1}{2r}\frac{\frac{\mathrm{d}}{\mathrm{d}r}\epsilon(r)}{\frac{\mathrm{d}}{\mathrm{d}r}P(r)}\bigg{\\{}\frac{\mathrm{d}}{\mathrm{d}r}\nu(r)+\frac{2G\left[r\frac{\mathrm{d}}{\mathrm{d}r}m(r)-m(r)\right]}{r^{2}\left[1-2\frac{Gm(r)}{r}\right]}\bigg{\\}}\Bigg{\\}}H(r)$
$\displaystyle+\frac{1}{2}\left\\{\frac{4}{r}+\frac{\mathrm{d}}{\mathrm{d}r}\nu(r)-\frac{2G\left[r\frac{\mathrm{d}}{\mathrm{d}r}m(r)-m(r)\right]}{r^{2}\left[1-2\frac{Gm(r)}{r}\right]}\right\\}\frac{\mathrm{d}}{\mathrm{d}r}H(r)$
$\displaystyle+\frac{{\mathrm{d}}^{2}}{\mathrm{d}r^{2}}H(r)=0.$ (73)

Being similar to the method of solving the moment of inertia in the slow
rotation approximation [121], we can calculate the tidal perturbation for a
given EoS $\mathcal{C}\left(P(r),\epsilon(r)\right)=0$ [119].

$\displaystyle\beta(r)-\frac{\mathrm{d}}{\mathrm{d}r}H(r)=0,$ (74)
$\displaystyle\frac{\mathrm{d}}{\mathrm{d}r}\beta(r)-\frac{2H(r)}{\left[1-2\frac{Gm(r)}{r}\right]}\Bigg{\\{}-2\pi
G\bigg{[}5\epsilon(r)+9P(r)$
$\displaystyle\qquad+\frac{\frac{\mathrm{d}}{\mathrm{d}r}\epsilon(r)}{\frac{\mathrm{d}}{\mathrm{d}r}P(r)}\left[P(r)+\epsilon(r)\right]\bigg{]}+\frac{3}{r^{2}}$
$\displaystyle\qquad+\frac{2}{1-2\frac{Gm(r)}{r}}G^{2}\left[\frac{M}{r^{2}}+4\pi
rP(r)\right]^{2}\Bigg{\\}}$
$\displaystyle\quad-\frac{2\beta(r)}{r\left[1-2\frac{Gm(r)}{r}\right]}\Bigg{\\{}-1+\frac{Gm(r)}{r}$
$\displaystyle\qquad+2\pi Gr^{2}\left[\epsilon(r)-P(r)\right]\Bigg{\\}}=0.$
(75)

The boundary conditions are

$\displaystyle H(r\rightarrow 0)\rightarrow{a_{0}}r^{2},$ (76)
$\displaystyle\beta(r\rightarrow 0)\rightarrow 2a_{0}r.$ (77)

The constant $a_{0}$ determines the extent of deformation of the object and
will be canceled in the calculation of the Love number. Eq. (74)-(77) can be
solved together with the TOV system Eq. (56), Eq. (57), Eq. (60), Eq. (65),
Eq. (58), Eq. (62), Eq. (63), and Eq. (64).

Figure 8: The tidal deformability as a function of the mass ratio between the
neutron star and the sun calculated in model I, model II, and model III. The
black dot with error bar is the tidal deformability of a $1.4M_{\odot}$-NS,
estimated in Ref. [41]: $\Lambda_{1.4}=190_{-120}^{+390}$. In model I, our
calculation result for the tidal deformability of a $1.4M_{\odot}$-NS is
$\Lambda_{1.4}=484$. In model II, $\Lambda_{1.4}=472$. In model III,
$\Lambda_{1.4}=484$. Our results of the tidal deformability in all three
models are consistent with this experiment data.

The $l=2$ Love number $k_{2}$ is [119]

$\displaystyle
k_{2}=\frac{8C^{5}}{5}\left(1-2C\right)^{2}\left[2+2C(y-1)-y\right]$
$\displaystyle\quad\times\bigg{\\{}2C\left[6-3y+3C(5y-8)\right]$
$\displaystyle\qquad+4C^{3}\left[13-11y+C\left(3y-2\right)+2C^{2}\left(1+y\right)\right]$
$\displaystyle\qquad+3\left(1-2C\right)^{2}\left[2-y+2C\left(y-1\right)\right]\ln{\left(1-2C\right)}\bigg{\\}}^{-1},$
(78)

where

$\displaystyle y=\frac{r_{B}\beta(r_{B})}{H(r_{B})},$ (79)

and

$\displaystyle C=\frac{GM}{r_{B}}$ (80)

is the compactness of the object. From Eq. (67), we then define a
dimensionless tidal deformability $\Lambda$ [64]:

$\displaystyle\Lambda=\frac{2}{3}\frac{k_{2}}{C^{5}}.$ (81)

The tidal deformability as a function of the mass ratio between the neutron
star and the sun calculated in model I, model II, and model III are shown in
Fig. 8. The black dot with error bar is the tidal deformability of a
$1.4M_{\odot}$-NS, estimated in Ref. [41]: $\Lambda_{1.4}=190_{-120}^{+390}$.
In model I, our calculation result for the tidal deformability of a
$1.4M_{\odot}$-NS is $\Lambda_{1.4}=484$. In model II, $\Lambda_{1.4}=472$. In
model III, $\Lambda_{1.4}=484$. Our results of the tidal deformability in all
three models are consistent with this experiment data.

In Ref. [42], the authors gave that the radius of a canonical $1.4M_{\odot}$
NS is $R_{1.4}=12.9_{-0.7}^{+0.8}\mathrm{km}$ and the tidal deformability is
$\Lambda_{1.4}=616_{-158}^{+273}$ by assuming that the secondary of GW190814
proves to be a NS. Because the secondary mass of GW190814 is
$m_{2}=2.59_{-0.09}^{+0.08}M_{\odot}$, that assumption requires that the max
possible mass of the NS should be no less than $m_{2}$. Model I and model III
satisfy this requirement. Actually, we find the results of $R_{1.4}$ and
$\Lambda_{1.4}$ given by the three models are all consistent with those from
Ref. [42].

Refs. [122, 123, 124] also discussed the sound velocity and tidal
deformability. For example, the authors in Ref. [122] used the pseudo-
conformal model [125, 126, 127] and presented that $\Lambda_{1.4}=785$ for
$n_{1/2}=2.0$, and $\Lambda_{1.4}=652$ for $n_{1/2}=3.0$ and $n_{1/2}=4.0$,
where $n_{1/2}$ denoted the density at which the topology change takes place.

##  VI Conclusion and discussion

We construct the hybrid equations of state (EoSs) of the cold strong
interaction matter using three holographic models in this work. In the low
baryon number density $n_{B}$ range, the HLPS intermediate EoS is used to
describe the nuclear matter. At the intermediate $n_{B}$ range, we use the
holographic QCD to calculate the EoS. The EoS from the holographic QCD is
connected with the HLPS EoS through the phase equilibrium condition. The key
point is to describe the nuclear matter in the holographic frame, the EoS of
which can be connected with the HLPS EoS through a second order phase
transition. It is found that the branch that corresponds to the holographic
nuclear matter is stable in the mass-radius relation of the NS.

The holographic action consists of two parts: the Einstein-Maxwell-dilaton
(EMD) system for the background and the KKSS action for the matter part. The
EMD system is solved to get the background metric which is used to calculate
the EoS. For holographic nuclear matter, the thermodynamic contribution made
by the KKSS action is also included to calculate the EoS. From the previous
study in Ref. [90], most probably the holographic nuclear matter is in the
quarkyonic phase [112] in which the chiral symmetry is restored but the quarks
and gluons are still confined. We will study the quarkyonic phase and its
relationship with the NS in our holographic models in the future. There are
already some articles [72, 60, 61, 62, 63, 64, 65, 66] to describe the baryons
by introducing the Dirac-Born-Infeld (DBI) action and the Chern-Simons (CS)
terms, where the solitons are regarded as baryons. However, the method to
calculate the EoS of the holographic nuclear matter here is new and hadn’t
been used within the holographic QCD yet.

The hybrid EoSs we constructed are located in the “allowed” range constrained
by interpolating the results of nuclear physics and the results of
perturbative QCD with the constraint of the astrophysical observations. Then
we calculate the properties of the neutron star, including the mass-radius
relation and the tidal deformability.

We use three different holographic models here. In model I, the radius for the
$1.4M_{\odot}$-NS is $12.320\mathrm{km}$. The max possible radius for the NS
is $12.833\mathrm{km}$, which correspond a $2.32M_{\odot}$-NS. The max
possible mass for the NS is $2.630M_{\odot}$, which correspond a
$12.559\mathrm{km}$-NS. In model II, the radius for the $1.4M_{\odot}$-NS is
$12.287\mathrm{km}$. The max possible radius for the NS is
$12.475\mathrm{km}$, which correspond a $1.98M_{\odot}$-NS. The max possible
mass for the NS is $2.466M_{\odot}$, which correspond a
$11.902\mathrm{km}$-NS. In model III, the radius for the $1.4M_{\odot}$-NS is
$12.319\mathrm{km}$. The max possible radius for the NS is
$12.795\mathrm{km}$, which correspond a $2.27M_{\odot}$-NS. The max possible
mass for the NS is $2.598M_{\odot}$, which correspond a $12.24\mathrm{km}$-NS.

The branch that corresponds to the holographic quark matter is always
unstable, and the holographic nuclear matter can produce the stable branch.
These results indicate that even in the core of the NS, the matter is still in
the confinement phase and the quark matter is not favored. This is consistent
with the conclusion in Ref. [116], in which the authors combined the input
from the gauge-gravity duality and the input from various $ab$-$initio$
methods for the nuclear matter at low density and analyzed families of hybrid
equations of state of cold QCD matter. Their prediction is that all neutron
stars are fully hadronic without quark matter cores.

In model I, our calculation result for the tidal deformability of a
$1.4M_{\odot}$-NS is $\Lambda_{1.4}=484$. In model II, $\Lambda_{1.4}=472$. In
model III, $\Lambda_{1.4}=484$. Our results of the tidal deformability in all
three models are consistent with that estimated in Ref. [41]. In Ref. [42],
the authors gave that the radius of a canonical $1.4M_{\odot}$ NS is
$R_{1.4}=12.9_{-0.7}^{+0.8}\mathrm{km}$ and the tidal deformability is
$\Lambda_{1.4}=616_{-158}^{+273}$ by assuming that the secondary of GW190814
proves to be a NS. Because the secondary mass of GW190814 is
$m_{2}=2.59_{-0.09}^{+0.08}M_{\odot}$, that assumption requires that the max
possible mass of the NS should be no less than $m_{2}$. Model I and model III
satisfy this requirement. Actually, we find the results of $R_{1.4}$ and
$\Lambda_{1.4}$ given by the three models are all consistent with those from
Ref. [42].

For future studies, we will investigate the relation between the chiral phase
transition and confinement-deconfinement phase transition, and figure out
whether the stable holographic nuclear matter is in the quarkyonic phase.
Currently, the EMD system and the KKSS action are solved separately, which
means we use the “probe approximation” here. The KKSS action is regarded as a
probe and its back-reaction to the background EMD system is neglected. We will
investigate this holographic model beyond the probe approximation and solve it
consistently in the future. Also, as we stated in the introduction, the
“allowed” range of the EoS we used here is derived by the “direct
interpolation approach” [43, 44]. However the “causality and thermodynamic
stability constraint approach” proposed in Refs. [45, 46] will give more
strict constraints on the EoS of cold QCD matter. We will consider these
constraints in future works.

###### Acknowledgements.

We thank Minghua Wei and Yidian Chen for helpful discussions. This work is
supported in part by the National Natural Science Foundation of China (NSFC)
Grant Nos. 11725523 and 11735007, the Strategic Priority Research Program of
Chinese Academy of Sciences under Grant Nos XDB34030000 and XDPB15, the start-
up funding from University of Chinese Academy of Sciences(UCAS), and the
Fundamental Research Funds for the Central Universities.

## References

  * [1] John B. Kogut.  An Introduction to Lattice Gauge Theory and Spin Systems.  Rev. Mod. Phys., 51:659, 1979.
  * [2] John B. Kogut.  A Review of the Lattice Gauge Theory Approach to Quantum Chromodynamics.  Rev. Mod. Phys., 55:775, 1983.
  * [3] Heng-Tong Ding, Frithjof Karsch, and Swagato Mukherjee.  Thermodynamics of strong-interaction matter from Lattice QCD.  Int. J. Mod. Phys. E, 24(10):1530007, 2015.
  * [4] Steven Weinberg.  Phenomenological Lagrangians.  Physica A, 96(1-2):327–340, 1979.
  * [5] Yoichiro Nambu and G. Jona-Lasinio.  Dynamical Model of Elementary Particles Based on an Analogy with Superconductivity. 1.  Phys. Rev., 122:345–358, 1961.
  * [6] Yoichiro Nambu and G. Jona-Lasinio.  DYNAMICAL MODEL OF ELEMENTARY PARTICLES BASED ON AN ANALOGY WITH SUPERCONDUCTIVITY. II.  Phys. Rev., 124:246–254, 1961.
  * [7] C. A. Dominguez and N. Paver.  LOCAL DUALITY CONSTRAINTS ON SCALAR GLUONIUM.  Z. Phys. C, 31:591, 1986.
  * [8] C. A. Dominguez and N. Paver.  Tensor Gluonium Spectrum in QCD.  Z. Phys. C, 32:391, 1986.
  * [9] J. I. Latorre, Stephan Narison, and S. Paban.  0++ TRIGLUONIUM SUM RULES.  Phys. Lett. B, 191:437–441, 1987.
  * [10] Stephan Narison and G. Veneziano.  QCD Tests of $G$ (1.6) = Glueball.  Int. J. Mod. Phys. A, 4:2751, 1989.
  * [11] Stephan Narison.  Masses, decays and mixings of gluonia in QCD.  Nucl. Phys. B, 509:312–356, 1998.
  * [12] Stephan Narison.  Masses, decays and mixings of gluonia in QCD.  Nucl. Phys. B Proc. Suppl., 64:210–219, 1998.
  * [13] Tao Huang, Hong-Ying Jin, and Ai-Lin Zhang.  Determination of the scalar glueball mass in QCD sum rules.  Phys. Rev. D, 59:034026, 1999.
  * [14] Stephan Narison.  Light scalar mesons in QCD.  Nucl. Phys. B Proc. Suppl., 186:306–311, 2009.
  * [15] Wei-jie Fu.  QCD at finite temperature and density within the fRG approach: An overview.  5 2022.
  * [16] Juan Martin Maldacena.  The Large N limit of superconformal field theories and supergravity.  Adv. Theor. Math. Phys., 2:231–252, 1998.
  * [17] S. S. Gubser, Igor R. Klebanov, and Alexander M. Polyakov.  Gauge theory correlators from noncritical string theory.  Phys. Lett. B, 428:105–114, 1998.
  * [18] Edward Witten.  Anti-de Sitter space and holography.  Adv. Theor. Math. Phys., 2:253–291, 1998.
  * [19] Martin Kruczenski, David Mateos, Robert C. Myers, and David J. Winters.  Towards a holographic dual of large N(c) QCD.  JHEP, 05:041, 2004.
  * [20] Tadakatsu Sakai and Shigeki Sugimoto.  Low energy hadron physics in holographic QCD.  Prog. Theor. Phys., 113:843–882, 2005.
  * [21] Tadakatsu Sakai and Shigeki Sugimoto.  More on a holographic dual of QCD.  Prog. Theor. Phys., 114:1083–1118, 2005.
  * [22] Leandro Da Rold and Alex Pomarol.  Chiral symmetry breaking from five dimensional spaces.  Nucl. Phys. B, 721:79–97, 2005.
  * [23] Joshua Erlich, Emanuel Katz, Dam T. Son, and Mikhail A. Stephanov.  QCD and a holographic model of hadrons.  Phys. Rev. Lett., 95:261602, 2005.
  * [24] Andreas Karch, Emanuel Katz, Dam T. Son, and Mikhail A. Stephanov.  Linear confinement and AdS/QCD.  Phys. Rev. D, 74:015005, 2006.
  * [25] J. M. Lattimer and M. Prakash.  The physics of neutron stars.  Science, 304:536–542, 2004.
  * [26] Eemeli Annala, Tyler Gorda, Aleksi Kurkela, Joonas Nättilä, and Aleksi Vuorinen.  Evidence for quark-matter cores in massive neutron stars.  Nature Phys., 16(9):907–910, 2020.
  * [27] Werner Becker et al.  Neutron stars and pulsars, volume 357.  Springer, 2009.
  * [28] Pawel Haensel, Aleksander Yu Potekhin, and Dmitry G Yakovlev.  Neutron Stars 1.  Springer, 2007.
  * [29] Bulent Kiziltan.  Reassessing the fundamentals: On the evolution, ages and masses of neutron stars.  University of California, Santa Cruz, 2010.
  * [30] Irwin I. Shapiro.  Fourth Test of General Relativity.  Phys. Rev. Lett., 13:789–791, 1964.
  * [31] Paul Demorest, Tim Pennucci, Scott Ransom, Mallory Roberts, and Jason Hessels.  Shapiro Delay Measurement of A Two Solar Mass Neutron Star.  Nature, 467:1081–1083, 2010.
  * [32] Emmanuel Fonseca et al.  The NANOGrav Nine-year Data Set: Mass and Geometric Measurements of Binary Millisecond Pulsars.  Astrophys. J., 832(2):167, 2016.
  * [33] John Antoniadis et al.  A Massive Pulsar in a Compact Relativistic Binary.  Science, 340:6131, 2013.
  * [34] H. T. Cromartie et al.  Relativistic Shapiro delay measurements of an extremely massive millisecond pulsar.  Nature Astron., 4(1):72–76, 2019.
  * [35] E. Fonseca et al.  Refined Mass and Geometric Measurements of the High-mass PSR J0740+6620.  Astrophys. J. Lett., 915(1):L12, 2021.
  * [36] M. H. van Kerkwijk, R. Breton, and S. R. Kulkarni.  Evidence for a Massive Neutron Star from a Radial-Velocity Study of the Companion to the Black Widow Pulsar PSR B1957+20.  Astrophys. J., 728:95, 2011.
  * [37] M. C. Miller et al.  PSR J0030+0451 Mass and Radius from $NICER$ Data and Implications for the Properties of Neutron Star Matter.  Astrophys. J. Lett., 887(1):L24, 2019.
  * [38] Thomas E. Riley et al.  A $NICER$ View of PSR J0030+0451: Millisecond Pulsar Parameter Estimation.  Astrophys. J. Lett., 887(1):L21, 2019.
  * [39] M. C. Miller et al.  The Radius of PSR J0740+6620 from NICER and XMM-Newton Data.  Astrophys. J. Lett., 918(2):L28, 2021.
  * [40] Thomas E. Riley et al.  A NICER View of the Massive Pulsar PSR J0740+6620 Informed by Radio Timing and XMM-Newton Spectroscopy.  Astrophys. J. Lett., 918(2):L27, 2021.
  * [41] B. P. Abbott et al.  GW170817: Measurements of neutron star radii and equation of state.  Phys. Rev. Lett., 121(16):161101, 2018.
  * [42] R. Abbott et al.  GW190814: Gravitational Waves from the Coalescence of a 23 Solar Mass Black Hole with a 2.6 Solar Mass Compact Object.  Astrophys. J. Lett., 896(2):L44, 2020.
  * [43] Eemeli Annala, Tyler Gorda, Aleksi Kurkela, and Aleksi Vuorinen.  Gravitational-wave constraints on the neutron-star-matter Equation of State.  Phys. Rev. Lett., 120(17):172703, 2018.
  * [44] Eemeli Annala, Tyler Gorda, Evangelia Katerini, Aleksi Kurkela, Joonas Nättilä, Vasileios Paschalidis, and Aleksi Vuorinen.  Multimessenger Constraints for Ultradense Matter.  Phys. Rev. X, 12(1):011058, 2022.
  * [45] Oleg Komoltsev and Aleksi Kurkela.  How Perturbative QCD Constrains the Equation of State at Neutron-Star Densities.  Phys. Rev. Lett., 128(20):202701, 2022.
  * [46] Tyler Gorda, Oleg Komoltsev, and Aleksi Kurkela.  Ab-initio QCD calculations impact the inference of the neutron-star-matter equation of state.  4 2022.
  * [47] Klaus Schertler, Stefan Leupold, and Jurgen Schaffner-Bielich.  Neutron stars and quark phases in the NJL model.  Phys. Rev. C, 60:025801, 1999.
  * [48] S. Lawley, Wolfgang Bentz, and Anthony William Thomas.  Neutron star properties from an NJL model modified to simulate confinement.  Nucl. Phys. B Proc. Suppl., 141:29–33, 2005.
  * [49] Takahide Kambe, Tetsuya Katayama, and Koichi Saito.  Equation of State for Neutron Star Matter with NJL Model and Dirac–Brueckner–Hartree–Fock Approximation.  JPS Conf. Proc., 14:020807, 2017.
  * [50] German Malfatti, Milva Orsaria, Gustavo A. Contrera, and Fridolin Weber.  Quark-hadron Phase Transition in Proto-Neutron Stars Cores based on a Non-local NJL Model.  Int. J. Mod. Phys. Conf. Ser., 45:1760039, 2017.
  * [51] T. Yazdizadeh and G. H. Bordbar.  The structure of cold neutron star with a quark core within the MIT and NJL models.  Iran. J. Sci. Technol. A, 43(5):2691–2698, 2019.
  * [52] Takuya Minamikawa, Toru Kojo, and Masayasu Harada.  Chiral condensates for neutron stars in hadron-quark crossover: From a parity doublet nucleon model to a Nambu–Jona-Lasinio quark model.  Phys. Rev. C, 104(6):065201, 2021.
  * [53] Takuya Minamikawa, Toru Kojo, and Masayasu Harada.  Chiral condensates for neutron stars in hadron-quark crossover; from a parity doublet nucleon model to an NJL quark model.  2 2022.
  * [54] Jan de Boer, Kyriakos Papadodimas, and Erik Verlinde.  Holographic Neutron Stars.  JHEP, 10:020, 2010.
  * [55] Kazuo Ghoroku, Kouki Kubo, Motoi Tachibana, Tomoki Taminato, and Fumihiko Toyoda.  Holographic cold nuclear matter as dilute instanton gas.  Phys. Rev. D, 87(6):066006, 2013.
  * [56] Kazuo Ghoroku, Kouki Kubo, Motoi Tachibana, and Fumihiko Toyoda.  Holographic cold nuclear matter and neutron star.  Int. J. Mod. Phys. A, 29:1450060, 2014.
  * [57] Carlos Hoyos, David Rodríguez Fernández, Niko Jokela, and Aleksi Vuorinen.  Holographic quark matter and neutron stars.  Phys. Rev. Lett., 117(3):032501, 2016.
  * [58] Niko Jokela, Matti Järvinen, and Jere Remes.  Holographic QCD in the Veneziano limit and neutron stars.  JHEP, 03:041, 2019.
  * [59] Luis A. H. Mamani, Cesar V. Flores, and Vilson T. Zanchin.  Phase diagram and compact stars in a holographic QCD model.  Phys. Rev. D, 102(6):066006, 2020.
  * [60] Takaaki Ishii, Matti Järvinen, and Govert Nijs.  Cool baryon and quark matter in holographic QCD.  JHEP, 07:003, 2019.
  * [61] Niko Jokela, Matti Järvinen, Govert Nijs, and Jere Remes.  Unified weak and strong coupling framework for nuclear matter and neutron stars.  Phys. Rev. D, 103(8):086004, 2021.
  * [62] Matti Järvinen.  Holographic modeling of nuclear matter and neutron stars.  Eur. Phys. J. C, 82(4):282, 2022.
  * [63] Nicolas Kovensky and Andreas Schmitt.  Isospin asymmetry in holographic baryonic matter.  SciPost Phys., 11(2):029, 2021.
  * [64] Nicolas Kovensky, Aaron Poole, and Andreas Schmitt.  Building a realistic neutron star from holography.  Phys. Rev. D, 105(3):034022, 2022.
  * [65] Nicolas Kovensky, Aaron Poole, and Andreas Schmitt.  Predictions for neutron stars from holographic nuclear matter.  SciPost Phys. Proc., 6:019, 2022.
  * [66] Lorenzo Bartolini, Sven Bjarke Gudnason, Josef Leutgeb, and Anton Rebhan.  Neutron stars and phase diagram in a hard-wall AdS/QCD model.  Phys. Rev. D, 105(12):126014, 2022.
  * [67] Sitthichai Pinkanjanarod and Piyabut Burikham.  Massive neutron stars with holographic multiquark cores.  Eur. Phys. J. C, 81(8):705, 2021.
  * [68] Sitthichai Pinkanjanarod, Piyabut Burikham, and Supakchai Ponglertsakul.  Tidal deformation and radial pulsations of neutron star with holographic multiquark core.  Eur. Phys. J. C, 82(2):141, 2022.
  * [69] Piyabut Burikham, Sitthichai Pinkanjanarod, and Supakchai Ponglertsakul.  Slowly rotating neutron star with holographic multiquark core: I-Love-Q relations.  Phys. Rev. D, 105(10):104018, 2022.
  * [70] David W Meltzer and Kip S Thorne.  Normal modes of radial pulsation of stars at the end point of thermonuclear evolution.  The Astrophysical Journal, 145:514, 1966.
  * [71] James M Bardeen, Kip S Thorne, and David W Meltzer.  A catalogue of methods for studying the normal modes of radial pulsation of general-relativistic stellar models.  The Astrophysical Journal, 145:505, 1966.
  * [72] Hiroyuki Hata, Tadakatsu Sakai, Shigeki Sugimoto, and Shinichiro Yamato.  Baryons from instantons in holographic QCD.  Prog. Theor. Phys., 117:1157, 2007.
  * [73] Matti Jarvinen and Elias Kiritsis.  Holographic Models for QCD in the Veneziano Limit.  JHEP, 03:002, 2012.
  * [74] U. Gursoy and E. Kiritsis.  Exploring improved holographic theories for QCD: Part I.  JHEP, 02:032, 2008.
  * [75] U. Gursoy, E. Kiritsis, and F. Nitti.  Exploring improved holographic theories for QCD: Part II.  JHEP, 02:019, 2008.
  * [76] Edward Witten.  Anti-de Sitter space, thermal phase transition, and confinement in gauge theories.  Adv. Theor. Math. Phys., 2:505–532, 1998.
  * [77] Zhen Fang, Yue-Liang Wu, and Lin Zhang.  Chiral phase transition and meson spectrum in improved soft-wall AdS/QCD.  Phys. Lett. B, 762:86–95, 2016.
  * [78] Zhen Fang, Yue-Liang Wu, and Lin Zhang.  Chiral Phase Transition with 2+1 quark flavors in an improved soft-wall AdS/QCD Model.  Phys. Rev. D, 98(11):114003, 2018.
  * [79] Zhen Fang, Yue-Liang Wu, and Lin Zhang.  Chiral phase transition and QCD phase diagram from AdS/QCD.  Phys. Rev. D, 99(3):034028, 2019.
  * [80] Zhen Fang, Yue-Liang Wu, and Lin Zhang.  Octet meson spectra and chiral phase diagram in the improved soft-wall AdS/QCD model.  Phys. Rev. D, 100(5):054008, 2019.
  * [81] Feryal Özel and Paulo Freire.  Masses, Radii, and the Equation of State of Neutron Stars.  Ann. Rev. Astron. Astrophys., 54:401–440, 2016.
  * [82] Jason Hessels, Scott Ransom, Mallory Roberts, Victoria Kaspi, Margaret Livingstone, Cindy Tam, and Fronefield Crawford.  Three new binary pulsars discovered with Parkes.  ASP Conf. Ser., 328:395, 2005.
  * [83] Ryan S. Lynch et al.  The Green Bank Telescope 350 MHz Drift-scan Survey II: Data Analysis and the Timing of 10 New Pulsars, Including a Relativistic Binary.  Astrophys. J., 763:81, 2013.
  * [84] K. Stovall et al.  The Green Bank Northern Celestial Cap Pulsar Survey - I: Survey Description, Data Analysis, and Initial Results.  Astrophys. J., 791(1):67, 2014.
  * [85] Ryan S. Lynch et al.  The Green Bank North Celestial Cap Pulsar Survey III: 45 New Pulsar Timing Solutions.  Astrophys. J., 859(2):93, 2018.
  * [86] Space Telescopes and Instrumentation 2012: Ultraviolet to Gamma Ray, volume 8443 of Society of Photo-Optical Instrumentation Engineers (SPIE) Conference Series, September 2012.
  * [87] B. P. Abbott et al.  GW170817: Observation of Gravitational Waves from a Binary Neutron Star Inspiral.  Phys. Rev. Lett., 119(16):161101, 2017.
  * [88] B. P. Abbott et al.  Multi-messenger Observations of a Binary Neutron Star Merger.  Astrophys. J. Lett., 848(2):L12, 2017.
  * [89] Rong-Gen Cai, Song He, and Danning Li.  A hQCD model and its phase diagram in Einstein-Maxwell-Dilaton system.  JHEP, 03:033, 2012.
  * [90] Xun Chen, Danning Li, Defu Hou, and Mei Huang.  Quarkyonic phase from quenched dynamical holographic QCD model.  JHEP, 03:073, 2020.
  * [91] Song He, Yi Yang, and Pei-Hung Yuan.  Analytic Study of Magnetic Catalysis in Holographic QCD.  4 2020.
  * [92] Xun Chen, Lin Zhang, Danning Li, Defu Hou, and Mei Huang.  Gluodynamics and deconfinement phase transition under rotation from holography.  JHEP, 07:132, 2021.
  * [93] Rong-Gen Cai, Song He, Li Li, and Yuan-Xu Wang.  Probing QCD critical point and induced gravitational wave by black hole physics.  1 2022.
  * [94] Zhibin Li, Yidian Chen, Danning Li, and Mei Huang.  Locating the QCD critical end point through the peaked baryon number susceptibilities along the freeze-out line.  Chin. Phys. C, 42(1):013103, 2018.
  * [95] Xun Chen, Danning Li, and Mei Huang.  Criticality of QCD in a holographic QCD model with critical end point.  Chin. Phys. C, 43(2):023105, 2019.
  * [96] Lin Zhang, Chutian Chen, Yidian Chen, and Mei Huang.  Spectra of glueballs and oddballs and the equation of state from holographic QCD.  Phys. Rev. D, 105(2):026020, 2022.
  * [97] Danning Li, Song He, Mei Huang, and Qi-Shu Yan.  Thermodynamics of deformed AdS5 model with a positive/negative quadratic correction in graviton-dilaton system.  JHEP, 09:041, 2011.
  * [98] Hermann Weyl and Henry L Brose.  Raum–zeit–materie [space–time–matter].  Lectures on General Relativity (in German)(Springer, Berlin, 1921), 1921.
  * [99] H Weyl and J Ehlers.  Space, time, matter: Lectures on general relativity.  Berlin, Germany: Springer, 1993.
  * [100] M. Tanabashi, K. Hagiwara, K. Hikasa, et al.  Review of particle physics.  Phys. Rev. D, 98:030001, Aug 2018.
  * [101] Oliver DeWolfe, Steven S. Gubser, and Christopher Rosen.  A holographic critical point.  Phys. Rev. D, 83:086005, 2011.
  * [102] Renato Critelli, Jorge Noronha, Jacquelyn Noronha-Hostler, Israel Portillo, Claudia Ratti, and Romulo Rougemont.  Critical point in the phase diagram of primordial quark-gluon matter from black hole physics.  Phys. Rev. D, 96(9):096026, 2017.
  * [103] S. Borsányi, Z. Fodor, J. N. Guenther, R. Kara, S. D. Katz, P. Parotto, A. Pásztor, C. Ratti, and K. K. Szabó.  Lattice QCD equation of state at finite chemical potential from an alternative expansion scheme.  Phys. Rev. Lett., 126(23):232001, 2021.
  * [104] Yuta Hamada, Elias Kiritsis, and Francesco Nitti.  Holographic Theories at Finite $\theta$-Angle, CP-Violation, Glueball Spectra and Strong-Coupling Instabilities.  Fortsch. Phys., 69(2):2000111, 2021.
  * [105] Song He, Shang-Yu Wu, Yi Yang, and Pei-Hung Yuan.  Phase Structure in a Dynamical Soft-Wall Holographic QCD Model.  JHEP, 04:093, 2013.
  * [106] David Dudal and Subhash Mahapatra.  Thermal entropy of a quark-antiquark pair above and below deconfinement from a dynamical holographic QCD model.  Phys. Rev. D, 96(12):126010, 2017.
  * [107] Wei-jie Fu, Jan M. Pawlowski, and Fabian Rennecke.  QCD phase structure at finite temperature and density.  Phys. Rev. D, 101(5):054032, 2020.
  * [108] Zhibin Li, Kun Xu, Xinyang Wang, and Mei Huang.  The kurtosis of net baryon number fluctuations from a realistic Polyakov–Nambu–Jona-Lasinio model along the experimental freeze-out line.  Eur. Phys. J. C, 79(3):245, 2019.
  * [109] K. Hebeler, J. M. Lattimer, C. J. Pethick, and A. Schwenk.  Equation of state and neutron star properties constrained by nuclear physics and observation.  Astrophys. J., 773:11, 2013.
  * [110] Paulo Bedaque and Andrew W. Steiner.  Sound velocity bound and neutron stars.  Phys. Rev. Lett., 114(3):031103, 2015.
  * [111] Len Brandes, Wolfram Weise, and Norbert Kaiser.  Inference of the sound speed and related properties of neutron stars.  8 2022.
  * [112] Larry McLerran and Sanjay Reddy.  Quarkyonic Matter and Neutron Stars.  Phys. Rev. Lett., 122(12):122701, 2019.
  * [113] Richard C. Tolman.  Effect of imhomogeneity on cosmological models.  Proc. Nat. Acad. Sci., 20:169–176, 1934.
  * [114] Richard C. Tolman.  Static solutions of Einstein’s field equations for spheres of fluid.  Phys. Rev., 55:364–373, 1939.
  * [115] J. R. Oppenheimer and G. M. Volkoff.  On massive neutron cores.  Phys. Rev., 55:374–381, 1939.
  * [116] Niko Jokela, Matti Järvinen, and Jere Remes.  Holographic QCD in the NICER era.  Phys. Rev. D, 105(8):086005, 2022.
  * [117] Eanna E. Flanagan and Tanja Hinderer.  Constraining neutron star tidal Love numbers with gravitational wave detectors.  Phys. Rev. D, 77:021502, 2008.
  * [118] Tanja Hinderer.  Tidal Love numbers of neutron stars.  Astrophys. J., 677:1216–1220, 2008.
  * [119] Tanja Hinderer, Benjamin D. Lackey, Ryan N. Lang, and Jocelyn S. Read.  Tidal deformability of neutron stars with realistic equations of state and their gravitational wave signatures in binary inspiral.  Phys. Rev. D, 81:123016, 2010.
  * [120] Tullio Regge and John A. Wheeler.  Stability of a Schwarzschild singularity.  Phys. Rev., 108:1063–1069, 1957.
  * [121] James B. Hartle and Kip S. Thorne.  Slowly Rotating Relativistic Stars. III. Static Criterion for Stability.  Astrophys. J. , 158:719, November 1969.
  * [122] Yong-Liang Ma and Mannque Rho.  Sound velocity and tidal deformability in compact stars.  Phys. Rev. D, 100(11):114003, 2019.
  * [123] A. Kanakis-Pegios, P. S. Koliogiannis, and Ch. C. Moustakidis.  Speed of sound constraints from tidal deformability of neutron stars.  Phys. Rev. C, 102(5):055801, 2020.
  * [124] Alkiviadis Kanakis-Pegios, Polychronis Koliogiannis, and Charalampos Moustakidis.  Constraints on the speed of sound of dense nuclear matter through the tidal deformability of neutron stars.  EPJ Web Conf., 252:05005, 2021.
  * [125] Won-Gi Paeng, Thomas T. S. Kuo, Hyun Kyu Lee, Yong-Liang Ma, and Mannque Rho.  Scale-invariant hidden local symmetry, topology change, and dense baryonic matter. II.  Phys. Rev. D, 96(1):014031, 2017.
  * [126] Yong-Liang Ma, Hyun Kyu Lee, Won-Gi Paeng, and Mannque Rho.  Pseudoconformal equation of state in compact-star matter from topology change and hidden symmetries of QCD.  Sci. China Phys. Mech. Astron., 62(11):112011, 2019.
  * [127] Yong-Liang Ma and Mannque Rho.  Pseudoconformal structure in dense baryonic matter.  Phys. Rev. D, 99(1):014034, 2019.