# A Short Proof of a Conjecture of Erdős
Proved by Moreira, Richter and Robertson

Bernard Host

###### Abstract

We give a short proof of a sumset conjecture of Erdős, recently proved by
Moreira, Richter and Robertson: every subset of the integers of positive
density contains the sum of two infinite sets. The proof is written in the
framework of classical ergodic theory.

title = A Short Proof of a Conjecture of Erdős Proved by Moreira, Richter and
Robertson author = Bernard Host, plaintextauthor = Bernard Host, runningtitle
= A Short Proof of a Conjecture of Erdős, keywords = Sumset, Erdős Conjecture,
Positive density, Følner sequence, year=2019, number=19, received=3 May 2019,
published=3 December 2019, doi=10.19086/da.11129,

[classification=text]

##  1 Introduction

In this note we give a short proof of a result of Moreira, Richter and
Robertson, conjectured by Erdős:

######  Theorem 1 ([MRR]).

If $A\subset{\mathbb{N}}$ has positive density then there exist infinite
subsets $B$ and $C$ of ${\mathbb{N}}$ with $A\supset B+C$.

Here and below, by the _density_ of a subset $A$ of ${\mathbb{N}}$ we mean its
upper uniform density ${\mathrm{d}}^{*}(A)$ defined by

${\mathrm{d}}^{*}(A):=\lim_{N\to\infty}\sup_{M\in{\mathbb{N}}}\frac{|A\cap[M,M+N)|}{N}.$

The proof extends to the case of countable abelian groups up to notational
changes. We believe that a generalization of the method to countable amenable
groups is possible, by using the maximal isometric factor in place of the
Kronecker factor.

The combinatorial part of our prof is closely related to the corresponding
part of [MRR]: Theorem 3 is the dynamical counterpart of [MRR, Theorem 2.2].
But the whole proof of Theorem 1 presented here is very different of [MRR], as
it is written inside the framework of classical ergodic theory.

In [MRR] the authors ask whether any set of positive density contains the sum
of three infinite sets. A possible strategy would be to show that it contains
the sum of a set of positive density and an infinite set, and then to use
Theorem 1. However, this method does not work:

######  Proposition 2.

There exists a set of positive density not containing any sum of a set of
positive density and an infinite set.

Section 2 contains some preliminaries. Theorem 1 is proved in Section 3 and
Proposition 2 in Section 4.

##  2 Preliminaries

In this paper we use the convention that ${\mathbb{N}}=\\{0,1,2,\dots\\}$.

###  2.1 Følner sequences

We make a constant use of Følner sequences but, since we are dealing only with
subsets of the integers, we never need general Følner sequences and we can
assume that a Følner sequence is a sequence
${\boldsymbol{\Phi}}=(\Phi_{N})_{N\geq 1}$ of intervals of ${\mathbb{N}}$
whose lengths tend to $\infty$. If $A\subset{\mathbb{N}}$ we write

$\displaystyle\overline{{\mathrm{d}}}_{\boldsymbol{\Phi}}(A)=\limsup_{N\to\infty}\frac{|A\cap\Phi_{N}]}{|\Phi_{N}|};$
$\displaystyle
d_{\boldsymbol{\Phi}}(A)=\lim_{N\to\infty}\frac{|A\cap\Phi_{N}]}{|\Phi_{N}|}\text{
if this limit exists.}$

Therefore,

${\mathrm{d}}^{*}(A)=\sup_{\boldsymbol{\Phi}}{\mathrm{d}}_{\boldsymbol{\Phi}}(A)$

where the supremum is taken for all Følner sequences ${\boldsymbol{\Phi}}$
such that ${\mathrm{d}}_{\boldsymbol{\Phi}}(A)$ exists, and the supremum is
attained.

###  2.2 Topological dynamical systems

A topological dynamical system $(X,T)$ is a compact metric space endowed with
a map $T\colon X\to X$, continuous and onto. We write $d_{X}$ for the distance
on $X$.

Let $x_{0}\in X$, ${\boldsymbol{\Phi}}$ a Følner sequence and $\mu$ a
probability measure on $X$. We say that $x_{0}$ is _generic for $\mu$ along
${\boldsymbol{\Phi}}$_ if

$\frac{1}{|\Phi_{N}|}\sum_{n\in\Phi_{N}}\delta_{T^{n}x_{0}}\to\mu\text{
weakly* as }N\to\infty$

where $\delta_{x}$ is the Dirac mass at $x$. This is equivalent to

$\text{for every
}f\in{\mathcal{C}}(X),\quad\frac{1}{|\Phi_{N}|}\sum_{n\in\Phi_{N}}f(T^{n}x_{0})\to\int
f\,d\mu\text{ as }N\to\infty.$

In this case, $\mu$ is supported on the closed orbit of $x_{0}$ under $T$ and
is invariant under $T$. For $U\subset X$ we have

$\displaystyle\overline{{\mathrm{d}}}_{\boldsymbol{\Phi}}\bigl{(}\\{(n\colon
T^{n}x_{0}\in U\\}\bigr{)}$ $\displaystyle\geq\mu(U)\text{ if $U$ is open;}$
$\displaystyle{\mathrm{d}}_{\boldsymbol{\Phi}}\bigl{(}\\{(n\colon
T^{n}x_{0}\in U\\}\bigr{)}$ $\displaystyle=\mu(U)\text{ if $U$ is clopen (that
is, open and closed).}$

We remark that for every $x_{0}\in X$, every Følner sequence admits a
subsequence along which $x_{0}$ is generic for some measure $\mu$.

We say that $x_{0}$ is generic for $\mu$ without mention of a Følner sequence
if it is generic for the usual Følner sequence $\bigl{(}[0,N)\bigr{)}_{N\geq
1}$.

###  2.3 The key combinatorial result

The next theorem is the dynamical counterpart of [MRR, Theorem 2.2].

######  Theorem 3.

Let $(X,T)$ be a topological dynamical system, $x_{0}\in X$, $E$ a clopen
subset of $X$ and

$A=\\{n\in{\mathbb{N}}\colon T^{n}x_{0}\in E\\}.$

Let $X\times X$ be endowed with the transformation $T\times T$. Let $x_{1}\in
X$ and let $\nu$ be a measure on $X\times X$, for which the point
$(x_{0},x_{1})$ is generic along some Følner sequence
${\boldsymbol{\Phi}}=(\Phi_{N})_{N\in{\mathbb{N}}}$. Assume that there exists
$\varepsilon>0$ and a sequence $(m_{i})$ of integers tending to $\infty$ with

$T^{m_{i}}x_{0}\to x_{1}\ \text{ and }\ \nu(T^{-m_{i}}E\times
E)\geq\varepsilon\text{ for every }i.$ (1)

Then there exist infinite subsets $B$, $C$ of ${\mathbb{N}}$ with $A\supset
B+C$.

This theorem can be deduced from [MRR, Proposition 2.5] but we give a proof
for completeness. The construction gives that the set $B$ ban be taken as
included in any given infinite subset of $\\{m_{i}\colon i\in{\mathbb{N}}\\}$;
the set $C$ can be taken as included in any infinite given subset of
$\\{\ell\in{\mathbb{N}}\colon T^{\ell}x_{1}\in E\\}$.

###### Proof.

We use the following result of Bergelson:

######  Lemma 4 ([1]).

Let $(Y,{\mathcal{Y}},\lambda)$ be a probability space, $\varepsilon>0$, and
$(B_{n})_{n\geq 1}$ a sequence in ${\mathcal{Y}}$ with
$\lambda(B_{n})\geq\varepsilon$ for every $n$. Then this sequence admits a
subsequence $(B_{n_{j}})_{j\geq 1}$ such that

$\text{ for every }k\geq
1,\quad\lambda\Bigl{(}\bigcap_{j=1}^{k}B_{n_{j}}\Bigr{)}>0.$

Applying this lemma to the probability $\nu$ on $X\times X$ and to the sets
$T^{-m_{i}}E\times E$, we obtain that the sequence $(m_{i})_{i\geq 1}$ admits
a subsequence, still written $(m_{i})$, such that

$\text{for every }r\geq 1,\quad\nu\Bigl{(}\bigcap_{i=1}^{r}(T^{-m_{i}}E\times
E)\Bigr{)}>0.$ (2)

We write

$L=\\{\ell\in{\mathbb{N}}\colon T^{\ell}x_{1}\in E\\}.$

For $n\in{\mathbb{N}}$, we have $(T\times T)^{n}(x_{0},x_{1})\in
T^{-m_{i}}E\times E$ if and only if $n\in(A-m_{i})\cap L$. By definition of
$\nu$ and (2),

$\text{for every }r\geq
1,\quad{\mathrm{d}}_{{\boldsymbol{\Phi}}}\Bigl{(}\bigcap_{i=1}^{r}(A-m_{i})\cap
L\Bigr{)}>0.$ (3)

and in particular this set is infinite.

We build by induction strictly increasing sequences $(b_{n})_{n\geq 1}$ and
$(c_{n})_{n\geq 1}$ with

$b_{n}\in\\{m_{i}\colon i\geq 1\\}\text{ and }c_{n}\in L\text{ for every }n;\
b_{i}+c_{j}\in A\text{ for all }i,j\geq 1.$

Let $n\geq 0$. If $n=0$ we chose an arbitrary $c_{1}\in L$. Assume that $n\geq
1$ and that $b_{1},\dots,b_{n}$ and $c_{1},\dots,c_{n}$ are built with
$b_{i}+c_{j}\in A$ for all $i,j\leq n$. By (3) applied with $r=b_{n}$ we can
chose $c_{n+1}>c_{n}$ with

$c_{n+1}\in\bigcap_{i=1}^{n}(A-b_{i})\cap L$

and we have $b_{i}+c_{n+1}\in A$ for $i\leq n$.

In both cases, for $j\leq n+1$ we have $c_{j}\in L$ that is $T^{c_{j}}x_{1}\in
E$. Since $E$ is open and since $T^{m_{i}}x_{0}\to x_{1}$, for every
sufficiently large $i$ we have $T^{m_{i}+c_{j}}x_{0}\in E$ for every $j\leq
n+1$, that is $m_{i}+c_{j}\in A$. We chose $i$ with this property and moreover
$m_{i}>b_{n}$ if $n\geq 1$, and we set $b_{n+1}=m_{i}$. We have
$b_{n+1}+c_{j}\in A$ for $1\leq j\leq n+1$ and we are done. ∎

##  3 Proof of Theorem 1

In this Section, $A$ is a subset of ${\mathbb{N}}$ with
$\delta={\mathrm{d}}^{*}(A)>0$, and we build the different objects appearing
in the statement of Theorem 3 and check that they satisfy the required
conditions. Theorem 1 then follows from Theorem 3.

###  3.1 Construction of a first system

Let $\\{0,1\\}^{\mathbb{N}}$ be endowed with the product topology and with the
shift $T$. Elements of $\\{0,1\\}^{\mathbb{N}}$ are written as
$x=(x(n))_{n\in{\mathbb{N}}}$. We consider ${\boldsymbol{1}}_{A}$ as an
element of $\\{0,1\\}^{\mathbb{Z}}$ that we write $x_{0}$, and define

  * (i)

$X$ is the closed orbit of $x_{0}$ under $T$.

Let $E$ be the cylinder set $E=\\{x\in X\colon x(0)=1\\}$. We have:

  * (ii)

$E$ is a clopen subset of $X$ and $A=\\{n\in{\mathbb{N}}\colon T^{n}x_{0}\in
E\\}$.

Let ${\boldsymbol{\Phi}}$ be a Følner sequence with
${\mathrm{d}}_{{\boldsymbol{\Phi}}}(A)=\delta$. Replacing
${\boldsymbol{\Phi}}$ by a subsequence we can assume that

  * (iii)

$x_{0}$ is generic along ${\boldsymbol{\Phi}}$ for some measure $\mu$ on $X$

and since $E$ is clopen we have

  * (iv)

$\mu(E)\geq\delta$.

###  3.2 Construction of an ergodic measure

Let $\mu=\int_{\Omega}\mu_{\omega}\,dP(\omega)$ be the ergodic decomposition
of $\mu$ under $T$. Since $\int_{\Omega}\mu_{\omega}(E)\,dP(\omega)=\delta$,
we have $P(\\{\omega\in\Omega\colon\mu^{\prime}_{\omega}(E)\geq\delta\\})>0$.
Therefore there exists an ergodic measure $\mu^{\prime}$ with
$\mu^{\prime}(E)\geq\delta$.

By the pointwise ergodic theorem, the measure $\mu^{\prime}$ admits a generic
point $x_{1}\in X$ in the usual sense, meaning that for every
$f\in{\mathcal{C}}(X)$ we have

$\frac{1}{N}\sum_{n=0}^{N-1}f(T^{n}x_{1})\to\int f\,d\mu\text{ as
}N\to\infty.$

Since the orbit of $x_{0}$ is dense and $T$ is continuous there exists a
sequence $(m_{N})_{N\geq 1}$ of integers such that

$\text{for every }N\text{ and }0\leq n<N,\
d_{X}(T^{m_{N}+n}x_{0},T^{n}x_{1})<1/N.$

For every $f\in{\mathcal{C}}(X)$, by uniform continuity of $f$ we have

$\frac{1}{N}\sum_{n=0}^{N-1}f(T^{m_{N}+n}x_{0})\to\int f\,d\mu\text{ as
}N\to\infty.$

This means that $x_{0}$ is generic of $\mu^{\prime}$ along the Følner sequence
${\boldsymbol{\Phi}}^{\prime}=(\Phi^{\prime}_{N})_{N\geq 1}$ where
$\Phi^{\prime}_{N}=[m_{N},m_{N}+N)$.

Substituting ${\boldsymbol{\Phi}}^{\prime}$ for ${\boldsymbol{\Phi}}$ and
$\mu^{\prime}$ for $\mu$, the properties (i)-(iv) remain valid and moreover

  * (v)

$\mu$ is ergodic under $T$.

###  3.3 Using the Kronecker factor

We recall that the Kronecker factor $(Z,m_{Z},S)$ of an ergodic system
$(X,\mu,T)$ is a compact abelian group $Z$ endowed with its Haar measure
$m_{Z}$, with a translation $S$ and with a factor map $\pi\colon X\to Z$,
characterized by

$\text{for all }\phi,\psi\in L^{2}(\mu)\text{ with
}{\mathbb{E}}_{\mu}(\psi\mid Z)=0\\\ \text{for every
}\varepsilon>0,\quad{\mathrm{d}}^{*}\Bigl{(}\Bigl{\\{}n\in{\mathbb{N}}\colon\Bigl{|}\int
T^{n}\psi\cdot\phi\,d\mu\Bigr{|}\geq\varepsilon\Bigr{\\}}\Bigr{)}=0.$ (4)

The transformation $S$ has the form

$\text{for every }z\in Z,\quad Sz=z+\alpha$ (5)

where $\alpha$ is a fixed element of $Z$.

It should be noticed that, even when $\mu$ is an invariant ergodic measure on
a topological dynamical system $(X,T)$, the factor map $\pi\colon X\to Z$ is
not continuous in general. The next Proposition explains how to modify the
given system and obtain continuity. It is stated in [2, 3] for a distal
system, and thus applies in particular to the isometric system $(Z,S)$.

######  Proposition 5 ([2, Proposition 6.1], [3, Proposition 24.3]).

Let $(X,T)$ be a topological dynamical system, $x_{0}\in X$, and $\mu$ an
ergodic invariant probability measure supported on the closed orbit of $x_{0}$
under $T$. Let $(Z,m_{Z},S)$ be the Kronecker factor of $(X,\mu,T)$, with
factor map $\pi\colon X\to Z$. Let $X\times Z$ be endowed with the
transformation $T\times S$. Let $\widetilde{\mu}$ be the measure on $X\times
Z$, image of $\mu$ under the map $x\mapsto(x,\pi(x))\colon X\to X\times Z$.
Then there exists a Følner sequence $\widetilde{\boldsymbol{\Phi}}$ and a
point $z_{0}$ of $Z$ such that $(x_{0},z_{0})$ is generic for
$\widetilde{\mu}$ along $\widetilde{\boldsymbol{\Phi}}$.

The conclusion means that, for every $f\in{\mathcal{C}}(X)$ and every
$h\in{\mathcal{C}}(Z)$ we have

$\int_{X}f\cdot
h\circ\pi\,d\mu=\lim_{N\to\infty}\frac{1}{|\widetilde{\Phi}_{N}|}\sum_{n\in\widetilde{\Phi}_{N}}f(T^{n}x_{0})\,h(S^{n}z_{0}).$

The first projection $p\colon(X\times Z,T\times T)\to(X,T)$ is a factor map in
the topological sense. Moreover, $p\colon(X\times Z,\widetilde{\mu},T\times
T)\to(X,\mu,T)$ is an isomorphism in the ergodic sense. In particular,
$(X\times Z,\widetilde{\mu},T\times T)$ is ergodic, its Kronecker factor is
$(Z,m_{Z},T)$, the factor map is the second projection and thus is continuous.
The key point in the Proposition is the genericity of $(x_{0},z_{0})$.

We apply this Proposition to the system $(X,\mu,T)$ and the point $x_{0}$
introduced above and let $(Z,m_{Z},S)$, $\pi$, $\widetilde{\mu}$ and $z_{0}$
and $\widetilde{\boldsymbol{\Phi}}$ be given by this Proposition. Let
$\widetilde{x_{0}}=(x_{0},z_{0})$ and let $\widetilde{X}$ be the closed orbit
of $\widetilde{x_{0}}$ in $X\times Z$. Since the point $\widetilde{x_{0}}$ is
generic for $\widetilde{\mu}$ along $\widetilde{\boldsymbol{\Phi}}$,
$\widetilde{\mu}$ is supported on $\widetilde{X}$. Let
$\widetilde{p}\colon\widetilde{X}\to X$ and
$\widetilde{\pi}\colon\widetilde{X}\to Z$ be the restriction to
$\widetilde{X}$ of the two natural projections $X\times Z\to X$ and $X\times
Z\to Z$ respectively, and let $\widetilde{E}=\widetilde{p}^{-1}(E)$. Then
$\widetilde{E}$ is a clopen subset of $\widetilde{X}$,
$A=\\{n\colon\widetilde{T}^{n}\widetilde{x_{0}}\in\widetilde{E}\\}$ and
${\mathrm{d}}_{\widetilde{\boldsymbol{\Phi}}}(A)=\widetilde{\mu}(\widetilde{E})=\mu(E)\geq\delta$.

Substituting $\widetilde{x_{0}}$ for $x_{0}$, $\widetilde{X}$ for $X$,
$\widetilde{\boldsymbol{\Phi}}$ for ${\boldsymbol{\Phi}}$, $\widetilde{\mu}$
for $\mu$ and $\widetilde{\pi}$ for $\pi$, the properties (i)-(v) above remain
valid and moreover,

  * (vi)

The factor map $\pi\colon X\to Z$ from $X$ to the Kronecker factor
$(Z,m_{Z},S)$ of $(X,\mu,T)$ is continuous.

We define

$z_{0}=\pi(x_{0}).$

###  3.4 Choosing a point $x_{1}$

We copy the next definition and results from [4].

######  Definition ([4]).

A Følner sequence ${\boldsymbol{\Phi}}=(\Phi_{N})$ is _tempered_ if there
exists a constant $C>0$ such that

$\Bigl{|}\bigcup_{k<N}(\Phi_{N}-\Phi_{k})\Bigr{|}\leq C|\Phi_{N}|$

for every $N$. (Here $\Phi_{N}-\Phi_{k}=\\{a-b\colon a\in\Phi_{N},\
b\in\Phi_{k}\\}$.)

######  Lemma 6 ([4]).

Every Følner sequence admits a tempered subsequence.

######  Theorem 7 ([4]).

If $(X,\mu,T)$ is a system and ${\boldsymbol{\Phi}}$ is a tempered Følner
sequence in ${\mathbb{N}}$ then for every $f\in L^{1}(\mu)$ the limit

$\lim_{N\to\infty}\frac{1}{|\Phi_{N}|}\sum_{n\in\Phi_{N}}f(T^{n}x)$

exists of $\mu$-almost every $x$. If $\mu$ is ergodic, then the limit is the
constant $\int f\,d\mu$.

######  Corollary 8.

Let $(X,T)$ be a topological dynamical system, $\mu$ an ergodic measure on $X$
and ${\boldsymbol{\Phi}}$ a tempered Følner sequence. Then $\mu$-almost every
$x\in X$ is generic for $\mu$ along ${\boldsymbol{\Phi}}$.

###### Proof.

It suffices to apply Theorem 7 to all functions in a countable subset of
${\mathcal{C}}(X)$, dense in ${\mathcal{C}}(X)$ under the uniform norm. ∎

Returning to our construction, we use Lemma 6 and replace the Følner sequence
with a tempered subsequence, still written ${\boldsymbol{\Phi}}$. We write

$X_{1}=\\{x\in X\colon x\text{ is generic for $\mu$ along
}{\boldsymbol{\Phi}}\\}.$

By Corollary 8 we have $\mu(X_{1})=1$. We define also

$X_{2}\text{ is the topological support of }\mu.$

We have $\mu(X_{2})=1$. Since the image of $\mu$ under $\pi$ is equal to
$m_{Z}$ we have $\pi(X_{2})=Z$.

We write

$\phi={\mathbb{E}}_{\mu}({\boldsymbol{1}}_{E}\mid Z)\text{ and
}\psi={\boldsymbol{1}}_{E}-\phi\circ\pi.$

We chose an open neighborhood $V$ of $0$ in $Z$ such that

$\int_{Z}|\phi(z)-\phi(z+\xi)|\,dm_{Z}(z)<\delta^{2}/4\text{ for every }\xi\in
V$ (6)

and define

$X_{3}=\pi^{-1}(V+z_{0}).$

Note that for every $n$ we have $T^{n}x_{0}\in X_{3}$ if and only if
$n\alpha\in V$.

$X_{3}$ is an open neighborhood of $x_{0}$ and, since $m_{Z}(V+z_{0})>0$, we
have $\mu(X_{3})>0$ and we can chose

$x_{1}\in X_{1}\cap X_{2}\cap X_{3}.$

We write

$z_{1}=\pi(x_{1})\text{ and }\beta=z_{1}-z_{0}.$

Since $x_{1}\in X_{3}$ we have $z_{1}\in V+z_{0}$ and

$\beta\in V.$

###  3.5 The joining $\nu$

Replacing the Følner sequence ${\boldsymbol{\Phi}}$ by a subsequence we can
assume that the point $(x_{0},x_{1})$ is generic along ${\boldsymbol{\Phi}}$
for a measure $\nu$ on $X\times X$. This measure is invariant under $T\times
T$; since $x_{0}$ is generic for $\mu$, the image of $\nu$ under the first
projection is equal to $\mu$; since $x_{1}\in X_{1}$, $x_{1}$ is generic for
$\mu$ too, and the second projection of $\nu$ is equal to $\mu$. We have that
$(X\times X,\nu,T\times T)$ is a joining of $(X,\mu,T)$ with itself.

Let $\lambda$ be the image of $\nu$ under the projection $\pi\times\pi\colon
X\times X\to Z\times Z$. Then $(Z\times Z,\lambda,S\times S)$ is a joining of
$(Z,m_{Z},S)$ with itself. Since $(x_{0},x_{1})$ is generic for $\nu$ along
${\boldsymbol{\Phi}}$, $(z_{0},z_{1})$ is generic for $\lambda$ along the same
Følner sequence. Therefore, $\lambda$ is supported on a single closed orbit
under $S\times S$ and, since $S\times S$ is a translation on the compact
abelian group $Z\times Z$, the measure $\lambda$ is ergodic under $S\times S$,
and for $\lambda$-almost every $(z,z^{\prime})$ we have
$z^{\prime}-z=z_{1}-z_{0}=\beta$. It follows that for all
$h,h^{\prime}\in{\mathcal{C}}(Z)$ we have

$\int_{Z\times
Z^{\prime}}h(z)h^{\prime}(z^{\prime})\,d\lambda(z,z^{\prime})=\int_{Z}h(z)h^{\prime}(z+\beta)\,dm_{Z}(z).$
(7)

By density the same relation holds for $h,h^{\prime}\in L^{2}(m_{Z})$.

Let $\theta\in L^{2}(\mu)$ be the image of ${\boldsymbol{1}}_{E}$ under the
Markov operator associated to the joining $\nu$. This function is
characterized by

$\text{for every }h\in L^{2}(\mu),\ \int_{X\times
X}h\otimes{\boldsymbol{1}}_{E}\,d\nu=\int_{X}h\cdot\theta\,d\mu.$

For every $n$ we have

$\nu(T^{-n}E\times E)-\int_{X\times
X}(S^{n}\phi\circ\pi)\otimes{\boldsymbol{1}}_{E}\,d\nu=\int_{X\times
X}T^{n}\psi\otimes{\boldsymbol{1}}_{E}\,d\nu\\\
=\int_{X}T^{n}\psi\cdot\,\theta\,d\mu.$

Since ${\mathbb{E}}_{\mu}(\psi\mid Z)=0$, by the characterization (4) of the
Kronecker factor we have

${\mathrm{d}}^{*}\Bigl{(}\Bigl{\\{}n\in{\mathbb{N}}\colon\Bigl{|}\nu(T^{-n}E\times
E)-\int_{X\times
X}(S^{n}\phi\circ\pi)\otimes{\boldsymbol{1}}_{E}\,d\nu\Bigr{|}\geq\delta^{2}/8\Bigr{\\}}\Bigr{)}=0.$
(8)

For every $n$, by definition of the measure $\lambda$ we have

$\int_{X\times
X}(S^{n}\phi\circ\pi)\otimes{\boldsymbol{1}}_{E}\,d\nu-\int_{X\times
X}(S^{n}\phi)\otimes\phi\,d\lambda\\\ =\int_{X\times
X}(S^{n}\phi\circ\pi)\otimes\psi\,d\nu=\int_{X\times
X}(\phi\circ\pi)\otimes(T^{-n}\psi)\,d\nu.$

Exchanging the role of the coordinates and proceeding as above we obtain

${\mathrm{d}}^{*}\Bigl{(}\Bigl{\\{}n\in{\mathbb{N}}\colon\Bigl{|}\int_{X\times
X}(S^{n}\phi\circ\pi)\otimes{\boldsymbol{1}}_{E}\,d\nu-\int_{Z\times
Z}(S^{n}\phi)\otimes\phi\,d\lambda\Bigr{|}\geq\delta^{2}/8\Bigr{\\}}\Bigr{)}=0.$
(9)

Combining (8) and (9) and using Formula (7) for $\lambda$ and formula (5) for
$S$ we obtain that

${\mathrm{d}}^{*}\Bigl{(}\Bigl{\\{}n\in{\mathbb{N}}\colon\Bigl{|}\nu(T^{-n}E\times
E)-\int_{Z}\phi(z+n\alpha)\phi(z+\beta)\,dm_{Z}(z)\Bigr{|}\geq\delta^{2}/4\Bigr{\\}}\Bigr{)}=0.$

We use now the fact that $\beta\in V$, the definition of this set and the
bound $\int\phi^{2}\,dm_{Z}\geq\delta^{2}$ and we obtain that the set

$\Lambda=\big{\\{}n\in{\mathbb{N}}\colon n\alpha\in V\text{ and
}\nu(T^{-n}E\times E)\leq\delta/4\bigr{\\}}$

satisfies

${\mathrm{d}}^{*}(\Lambda)=0.$

###  3.6 Conclusion

We check now that the point $x_{1}$ and the measure $\nu$ satisfy all the
hypothesis needed to apply Theorem 3. It remains only to build a good sequence
$(m_{i})$.

For every integer $i>0$ the set

$X_{3}\cap\\{x\in X\colon d_{X}(x,x_{1})<1/i\\}$

is an open neighborhood of $x_{1}$, and since $x_{1}$ belongs to the
topological support $X_{2}$ of $\mu$, the measure of this set is $>0$ and thus
the set

$F_{i}=\\{n\in{\mathbb{N}}\colon T^{n}x_{0}\in X_{3}\text{ and
}d_{X}(T^{n}x_{0},x_{1})<1/i\\}$

satisfies $\overline{{\mathrm{d}}}_{\boldsymbol{\Phi}}(F_{i})>0$. Since
$\overline{{\mathrm{d}}}_{\boldsymbol{\Phi}}(\Lambda)=0$,
$\overline{{\mathrm{d}}}_{\boldsymbol{\Phi}}(F_{i}\setminus\Lambda)>0$ and in
particular $F_{i}\setminus\Lambda$ is infinite. We thus can chose a sequence
of integers $(m_{i})$, tending to $\infty$, with $m_{i}\in
F_{i}\setminus\Lambda$ for every $i$. We have that

$T^{m_{i}}x_{0}\to x_{1}\text{ as }i\to\infty$

and, for every $i$,

$\displaystyle T^{m_{i}}x_{0}\in X_{3}\text{ and thus }m_{i}\alpha\in V;$
$\displaystyle m_{i}\notin\Lambda\text{ and thus }\nu(T^{-m_{i}}E\times
E)>\delta^{2}/4.$

All the hypothesis of Theorem 3 are satisfied with $\varepsilon=\delta^{2}/4$
and we are done. ∎

##  4 Proof of Proposition 2

Let $(X,T)$ be an uniquely ergodic topological dynamical system with invariant
measure $\mu$ and assume that the measure preserving system $(X,\mu,T)$ is
mixing. Let $U\subset X$ be an open subset with
$0<\mu(U)\leq\mu(\overline{U})<1$, $x_{0}\in X$ and

$A=\\{n\in{\mathbb{N}}\colon T^{n}x_{0}\in U\\}.$

Since $\mu(U)>0$ and by unique ergodicity, ${\mathrm{d}}^{*}(A)>0$. We show
that it does not contain a sum $B+C$ with ${\mathrm{d}}^{*}(B)>0$ and $C$
infinite.

Assume by contradiction that $B$ and $C$ exist withe these properties and let

$K=\overline{\\{T^{n}x_{0}\colon n\in B\\}}.$

By unique ergodicity,

$\mu(K)\geq{\mathrm{d}}^{*}(\\{n\in{\mathbb{N}}\colon T^{n}x_{0}\in
K\\})\geq{\mathrm{d}}^{*}(B)>0.$

On the other hand, for every $c\in C$ we have $T^{c}\\{T^{b}x_{0}\colon b\in
B\\}\subset U$ and thus

$\text{for every }c\in C,\quad T^{c}K\subset\overline{U}.$

Since $(X,\mu,T)$ is mixing and since $K$ and $X\setminus\overline{U})$ have
positive measure, for every sufficiently large $n$ we have
$\mu(T^{n}K\cap(X\setminus\overline{U}))>0$ and thus
$T^{n}K\not\subset\overline{U}$. We have a contradiction. ∎

## References

  * [1] Vitaly Bergelson.  Sets of recurrence of ${\mathbb{Z}}^{m}$-actions and properties of sets of differences in ${\mathbb{Z}}^{m}$.  J. London Math. Soc. (2), 31(2):295–-304, 1985.
  * [2] Bernard Host and Bryna Kra.  Uniformity seminorms on $\ell^{\infty}$ and applications.  J. Anal. Math., 108:219–276, 2009.
  * [3] Bernard Host and Bryna Kra.  Nilpotent Structures in Ergodic Theory. Mathematical Surveys and Monographs 236, American Mathematical Society, Providence, RI, 2018.
  * [4] Elon Lindenstrauss.  Pointwise convergence for amenable groups.  Invent. Math., 146(2):259–295, 2001.
  * [MRR] Joel Moreira, Florian Karl Richter and Donald Robertson  A proof of a sumset conjecture of Erdős.  Ann. of Math. (2), 89(2):605-–652, 2019.

[bh] Bernard Host
Professor
LAMA, Université Paris-Est Marne-la-Vallée, France.
bernard.hostu-pemfr


# Resolving phase transition properties of dense matter through tidal-excited
g-mode from inspiralling neutron stars

Zhiqiang Miao Department of Astronomy, Xiamen University, Xiamen, Fujian
361005, China  Enping Zhou Department of Astronomy, School of Physics,
Huazhong University of Science and Technology, Wuhan, 430074, China  Ang Li
Department of Astronomy, Xiamen University, Xiamen, Fujian 361005, China Ang
Li<EMAIL_ADDRESS>
###### Abstract

The investigation of the phase state of dense matter is hindered by
complications of first-principle nonperturbative quantum chromodynamics. By
performing the first consistent general-relativistic calculations of tidal-
excited g-mode of neutron stars with a first-order strong interaction phase
transition in the high-density core, we demonstrate that gravitational wave
signal during binary neutron star inspiral probes their innermost hadron-quark
transition and provides potent constraints from present and future
gravitational-wave detectors.

Gravitational waves (678); High energy astrophysics (739) Neutron star cores
(1107); Relativistic stars(1392)

##  1 Introduction

Despite decades of exploration, the phase state of cold, dense matter is still
a fundamental problem that remains unresolved. This problem pertains to low-
energy nonperturbative quantum chromodynamics (QCD) (Kurkela et al., 2010) and
is one of the key scientific goals in current and anticipated relativistic
heavy-ion collision experiments (Huth et al., 2022). It is widely believed
that this type of matter exists in the core of neutron stars (NSs) (Weber,
1999; Heiselberg & Pandharipande, 2000; Lattimer & Prakash, 2001). The
equation of state (EOS) of dense stellar matter has drawn much attention in
the era of gravitational wave (GW) astronomy. This is because the observations
of GW170817 led to the first measurement of NS tidal deformability (Abbott et
al., 2017a), which offered valuable insight into the phase state of dense
matter (Abbott et al., 2018).

The tidal deformability accounts for the quadrupole deformations of
inspiralling NSs, which are induced by the tidal force from each other, when
the quasi-equilibrium approximation can be utilized (i.e., the tidal field
changes slowly). While in the late stage of the inspiral phase, once the
binary’s orbital frequency approaches the frequency of a particular normal
mode of NSs, a resonance can occur. This resonance will transfer the orbital
energy to the stellar oscillation via a mechanism referred to as the dynamical
tide. This event will accelerate the orbital decay and result in a phase shift
in the GW waveform, which will provide insights into the NS inner structure
with the help of GW detection from inspiralling NSs.

Previous studies have focused on the interfacial (i) modes, which are excited
at the interface between the solid crust and fluid core (Tsang et al., 2012;
Pan et al., 2020; Zhu et al., 2023). These modes typically have frequencies of
$\sim 100\;\text{Hz}$ and contribute $\sim\mathcal{O}(0.1)$ phase shifts
during resonance, making them promising candidates for detection using
Advanced LIGO (aLIGO) (LIGO Scientific Collaboration et al., 2015) and 3rd
generation GW detectors such as Einstein Telescope (ET) (Maggiore et al.,
2020) and Cosmic Explorer (CE) (Reitze et al., 2019). Attempts have also been
made to infer NS EOS through tide-induced f-modes and g-modes, which arise
from the composition gradients within the NS’s interior. Studies (see e.g.,
Kokkotas & Schmidt, 1999; Ho & Lai, 1999; Reisenegger & Goldreich, 1994; Lai,
1994; Yu & Weinberg, 2017a; Xu & Lai, 2017) have shown that f-modes are
challenging to excite due to their high resonant frequencies, whereas g-modes
demonstrate negligible phase shifts because their couplings to tidal fields
are meager.

In this work, we investigate a unique type of g-modes, known as discontinuity
g-modes, which are expected to occur in NSs that possibly have a strong phase
transition in their high-density inner cores (Finn, 1987)111We study in the
context of slow phase transition, namely the timescale of hadron-quark
conversion is much longer than the typical oscillation timescale. An
alternative fast conversion would result in a g-mode with zero frequency
(Tonetto & Lugones, 2020) and the present analysis will not apply.. Our study
showcases how discontinuity g-modes’ dynamical tide can be leveraged to probe
the dense matter’s phase transition properties using the data obtained from
GW170817, along with possible observations from the next-generation GW
detectors.

##  2 Formalism and methods

###  2.1 Relativistic g-modes and tidal coupling

The existence of a quark core inside NSs (also called hybrid stars) has been
discussed with Newtonian formulations of perturbations and tidal couplings,
focusing on the i-mode at the interface of a crystalline quark core and a
fluid hadronic envelope (Lau & Yagi, 2021). With the solid quark core
condition relaxed, we here employ a full general relativistic (GR) formalism
as follows: 1) Constructing the background spherical hybrid stars in
hydrostatic equilibrium using Einstein’s field equations, i.e., solving the
Tolman-Oppenheimer-Volkoff (TOV) equations, and 2) Computing their
discontinuity g-modes by solving the linear pulsation equations in general
relativity (Lindblom & Detweiler, 1983; Detweiler & Lindblom, 1985). Our study
focuses solely on the $l=2$ g-modes, as they have the strongest coupling to
the tidal gravitational field.

Our approach prioritizes maintaining consistency in how we treat the
oscillation mode and the tidal overlap integral $Q$ within the framework of GR
(see Appendix .2 for a derivation of this form):

$Q=\frac{1}{MR^{2}}\int\sqrt{-g}e^{-2\Phi}d^{3}x(p+\varepsilon)\xi^{*\mu}\nabla_{\mu}(r^{2}Y_{2m})\
.$ (1)

Here, $g$ represents the determinant of the metric of the background spacetime
and $e^{2\Phi}=g_{tt}$ is the $(t,t)$ component of the metric. $Y_{2m}$ is the
$l=2$ spherical harmonic.
$\xi^{\mu}=(0,re^{-\lambda}WY_{2m},V\partial_{\theta}Y_{2m},V\sin^{-2}\theta\partial_{\phi}Y_{2m})e^{i\omega
t}$ is the Lagrangian displacement of the eigenmode with the normalization
$\langle\xi^{\mu}|\xi^{\mu}\rangle=MR^{2}$, where $e^{2\lambda}=g_{rr}$ is the
$(r,r)$ component of the metric and
$\langle\xi^{\mu}|\psi^{\mu}\rangle=\int\sqrt{-g}e^{-2\Phi}d^{3}x(p+\varepsilon)\xi^{*\mu}\psi_{\mu}$
is the relativistic inner product. $f$ is the frequency of the mode and
$\omega=2\pi f$ is the angular frequency. The present self-consistent approach
effectively mitigates the issue of “f-mode contamination” (Reisenegger, 1994)
in the calculation of $Q$ encountered in the hybrid (GR background + Newtonian
perturbation) approach commonly used previously. We provide the full
description of this relativistic formula in Appendix .2, where we’ve also
validated that our relativistic $Q$ satisfies both orthogonality and the sum
rule.

We determine the mass-radius ($M$-$R$) relationship of NSs through the
solution of the TOV equations, using the EOS - i.e., pressure ($p$) as a
function of energy density ($\varepsilon$) - as our basic input. To describe
the fiducial soft and stiff hadronic matter (HM), we employ the QMF (Zhu et
al., 2018) and NL3$\omega\rho$ (Horowitz & Piekarewicz, 2001) models, ensuring
their consistency with results from laboratory nuclear experiments. For
modeling the sharp first-order hadron-quark phase transition and representing
quark matter within the high-density cores of NSs, we utilize the well-
established generic constant-speed-of-sound (CSS) parameterization (Alford et
al., 2013). Previous studies have thoroughly discussed the weak density
dependence of the speed of sound in quark matter (see e.g., Alford et al.,
2013, 2015). This approach enables us to establish the EOS from the start of
the phase transition up to the maximum central density of the star by using
three phase-transition parameters: the transition pressure $p_{\rm t}$ (or the
transition density $\varepsilon_{\rm t}\equiv\varepsilon_{\rm HM}(p_{\rm
t})$), the transition strength $\Delta\varepsilon$ and the speed of sound in
quark matter $c_{s}$. The transition density and transition strength affect
the tidal overlap integral and its typical value remains at
$\sim\mathcal{O}(0.01)$.

It is worth noting that since the discontinuity g-modes are a result of the
discontinuity in the density profile of the star, we expect the frequency of
these modes to depend on both the transition strength and the transition
density. Based on our calculations, we observe that the g-mode frequency
typically falls within the range of $100\;\text{Hz}$ to $1500\;\text{Hz}$,
which is strongly dependent on the transition strength. The transition
strength varies between $\sim 5\,{\rm MeV/fm^{3}}$ and $\sim 300\,{\rm
MeV/fm^{3}}$. Our study also indicates that the speed of sound ($c_{s}$) has a
subdominant effect on mode frequency. Specifically, the difference in mode
frequency between the $c_{s}=1/\sqrt{3}$ case (corresponding to the conformal
limit) and the $c_{s}=1$ case (corresponding to the casual limit) is found to
be less than $\sim 10\%$. Thus, for the purpose of our analysis, we solely
focus on extremely stiff quark matter with a speed of sound $c_{s}=1$.

Figure 1: Total phase shift $|\delta\bar{\phi}|$ and resonant frequency
$\bar{f}$ for a $1.4\,{M}_{\odot}-1.4\,{M}_{\odot}$ system (left panel) and a
$1.8\,{M}_{\odot}-1.8\,{M}_{\odot}$ system (right panel), depending on the
transition densities (represented by $p_{\rm t}/\varepsilon_{\rm t}$). $f_{\rm
ISCO}\simeq 4400\;\text{Hz}/(M+M^{\prime})$ denotes the GW frequency
corresponding to the innermost stable circular orbit, at which we assume the
signal of inspiral phase shut off (Cutler & Flanagan, 1994). Also plotted are
the detectability thresholds of aLIGO, ET and CE, which are calculated from
their designed sensitivities by assuming the system is located at
$D_{L}=100\,{\rm Mpc}$.  Figure 2: Detectable regions of phase transition
parameters for a $1.4\,{M}_{\odot}-1.4\,{M}_{\odot}$ system (orange) and a
$1.8\,{M}_{\odot}-1.8\,{M}_{\odot}$ system (blue) on the 3rd generation
detectors. The red region depicts the parameter space ruled out by the mode
resonance analysis of GW170817 at 95% confidence level, and the grey regions
represent the portions eliminated by the maximum mass of observed pulsars. The
results of QMF+CSS and NL3$\omega\rho$+CSS are represented by dotted and
dashed lines, respectively. The critical line corresponds to
$\Delta\varepsilon/\varepsilon_{\rm t}=1/2+3p_{\rm t}/2\varepsilon_{\rm t}$,
dividing the parameter space into two parts. Below the critical line, there is
a hybrid star branch connected to the hadronic branch, while above the
critical line, there is no connected hybrid star branch (Alford et al., 2013).

###  2.2 Waveform correction

Throughout the inspiral process, the binary’s individual components will be
subject to external tidal forces from their respective companions. Upon the
tidal driving frequency nearing the g-mode frequency of the NSs, it becomes
possible to trigger internal stellar oscillations through resonance,
ultimately causing a phase shift in the GW waveform. With the mode frequency
and mode eigenfunction obtained in our GR framework, we then are able to
estimate the phase shift of the waveform as (Lai, 1994):

$\delta\phi\simeq-0.12\left(\frac{Q}{0.01}\right)^{2}\left(\frac{600\;\text{Hz}}{f}\right)^{2}\left(\frac{2q}{1+q}\right)M_{1.4}^{-4}R_{12}^{2}\
,$ (2)

where $q=M^{\prime}/M$ with $M^{\prime}$ representing the companion star mass.
$M_{1.4}=M/1.4\,{M}_{\odot}$ and $R_{12}=R/12\hbox{$\,{\rm km}$}$. The total
phase correction for a binary NS waveform $h(f)=\mathcal{A}e^{i\Psi(f)}$,
which takes into account the mode resonances of both components, is expressed
as (Yu & Weinberg, 2017b; Pan et al., 2020):

$\begin{split}\Delta\Psi(f)&=-\sum_{i}\delta\phi_{i}\left(1-\frac{f}{f_{i}}\right)\Theta(f-f_{i})\\\
&\approx-\delta\bar{\phi}\left(1-\frac{f}{\bar{f}}\right)\Theta(f-\bar{f})\
,\end{split}$ (3)

where $f_{i}$ is the g-mode frequency of the $i$-th star and
$\delta\bar{\phi}_{i}$ represents the corresponding phase shift due to
resonance. Here, $\Theta(\cdot)$ is the Heaviside step function. The second
line of Eq. (3) introduces the total phase shift
$\delta\bar{\phi}=\sum_{i}\delta\phi_{i}$ and the resonant frequency
$\bar{f}=\delta\bar{\phi}/\sum_{i}(\delta\phi_{i}/f_{i})$ to reduce the number
of parameters (Pan et al., 2020). For simplicity, $|\delta\bar{\phi}|$ will be
used as the waveform parameter from here on.

##  3 Results

###  3.1 Detectability

Fig. 1 displays the total phase shift $|\delta\bar{\phi}|$ against the
resonant frequency $\bar{f}$ for two equal-mass binary systems, along with the
detection thresholds of aLIGO and 3rd generation detectors. These thresholds
are computed using a Fisher matrix method (Cutler & Flanagan, 1994) and are
based on the designed sensitivities (Abbott et al., 2017b) assuming that the
binary systems are located at $D_{L}=100\,{\rm Mpc}$ (details in Appendix .5).
Generally, $|\delta\bar{\phi}|$ decreases with the increase of $\bar{f}$ since
$\delta\phi\propto f^{-2}$ [cf. Eq. (2)] (Lai, 1994). For a significant number
of hybrid star models, the total phase shifts are greater than the
detectability threshold of 3rd generation detectors. It is worth noting that
there is a threshold frequency $f_{\rm th}$, which is independent of the
hadronic EOS. Binaries with $\bar{f}\lesssim f_{\rm th}$ may be detectable by
3rd generation detectors, while those with $\bar{f}\gtrsim f_{\rm th}$ may not
be detectable. For a $1.4\,{M}_{\odot}-1.4\,{M}_{\odot}$ system, $f_{\rm
th}\simeq 1200\;\text{Hz}$, while for a $1.8\,{M}_{\odot}-1.8\,{M}_{\odot}$
system, $f_{\rm th}\simeq 1050\;\text{Hz}$. These threshold frequencies lead
to constraints on the phase transition parameters, which are illustrated in
the shaded regions (orange for $1.4\,{M}_{\odot}-1.4\,{M}_{\odot}$ system and
blue for $1.8\,{M}_{\odot}-1.8\,{M}_{\odot}$ system) in Fig. 2. It is evident
that the detection of a signal of mode resonance provides a new opportunity to
probe the properties of the hadron-quark phase transition, which can serve as
a complementary constraint to other observations, such as the $\sim
2\,{M}_{\odot}$ maximum mass constraint from the most massive pulsars (grey
regions in Fig. 2).

Figure 3: Constraints on phase transition densities (represented by $p_{\rm
t}/\varepsilon_{\rm t}$) from GW170817 data for QMF+CSS (upper panel) and
NL3$\omega\rho$+CSS (lower panel). Colored bands depict calculations for a
GW170817-like system with chirp mass $\mathcal{M}=1.186\,{M}_{\odot}$ and mass
ratio $q=0.7$–$1$. The posterior probability density for $\Delta|\bar{\phi}|$
and $\bar{f}$ inferred from the GW170817 data is shown in blue heatmap, with
the 95% and 68% boundaries being overlaid by dashed lines.

###  3.2 GW170817 constraints

We proceed by making use of the available GW170817 data and conduct a search
for mode resonance signals in full range of plausible resonant frequencies
$[50,1600]\;\text{Hz}$ (the upper limit corresponds to the shut-off frequency
of the waveform, $4400\;\text{Hz}/(M+M^{\prime})\approx 1600\;\text{Hz}$, for
the GW170817 system)222In principle, the g-mode frequency can be made
arbitrarily small if $\Delta\varepsilon$ is sufficiently tiny. However, given
the substantial noise levels at low frequencies ($\lesssim 10\;\text{Hz}$), we
have opted for a lower limit of $50\;\text{Hz}$ for the search. This choice
corresponds to a transition strength of $\Delta\varepsilon\simeq 1\,{\rm
MeV}$., with the application of waveform correction as described in Eq. (3).
This is similar to a previous search conducted in Pan et al. (2020) for
signals in the frequency range of $[30,300]\;\text{Hz}$. We employ a Bayesian
parameter estimation approach, following the data analysis framework presented
in Abbott et al. (2019). Our analysis utilizes PyMultiNest (Buchner et al.,
2014) implemented in Bilby (Ashton et al., 2019), enabling us to
simultaneously obtain the posterior distribution of parameters and the model
evidence $\mathcal{Z}$. By denoting $\mathcal{H}_{1}$ as the hypothesis with
mode resonance and $\mathcal{H}_{0}$ as the one without, we are able to
compute the Bayes factor from the evidence,
$\mathcal{B}^{1}_{0}=\mathcal{Z}_{1}/\mathcal{Z}_{0}$. We obtain a Bayes
factor of $\mathcal{B}^{1}_{0}=0.72$, which implies that the current data does
not strongly favor one hypothesis over the other (see Appendix .6 for more
discussions on the Bayesian parameter estimation).

Though the search yields no results, the posterior distribution of
$\Delta|\bar{\phi}|$ and $\bar{f}$ can be utilized to set constraints on the
phase transition parameters. In Fig. 3 we present the calculated
$\Delta|\bar{\phi}|$ and $\bar{f}$ for a GW170817-like system with a chirp
mass of $\mathcal{M}=1.186\,{M}_{\odot}$ and a mass ratio range of
$q=0.7$–$1$. The posterior distribution for $\Delta|\bar{\phi}|$ and $\bar{f}$
together with their 95% and 68% contours are also plotted. It can be seen that
for a certain transition density, the hybrid star models with small resonant
frequencies should be excluded. For example, those models with $p_{\rm
t}/\varepsilon_{\rm t}=0.11$ and $\bar{f}\leq 485\;\text{Hz}$ are ruled out at
95% confidence level for QMF+CSS. Converting these constraints on $\bar{f}$ to
phase transition properties, we conclusively exclude the possibility of a weak
phase transition occurring at low densities with a confidence level of 95%.
This exclusion is represented in Fig. 2 by the red region. Our conclusion here
is in agreement with our previous findings based on static tidal deformability
constraints from GW170817 (Miao et al., 2020) and the study of proton and
lambda production during relativistic heavy-ion collisions at several
GeV/nucleon incident beam energies (Li et al., 2023).

##  4 Discussion

The investigation of phase transitions in QCD matter is a central area of
interest in physics and has relevance for a range of open problems, such as
early Universe (Hogan, 1983; Cai et al., 2022), core-collapse supernovae
(Fischer et al., 2018), cosmological gamma-ray bursts (Cheng & Dai, 1996), as
well as the nonperturbative properties of QCD itself (Huber, 2020). We propose
to probe the phase transition properties possible in dense NS matter through
the tidal resonant excitation of g-modes.

Our analysis demonstrates that the tidal excitation of g-modes from
inspiralling neutron stars can induce significant phase shifts in the
gravitational waveform, on the order of $\sim\mathcal{O}(0.1)-\mathcal{O}(1)$.
These shifts are detectable by both aLIGO and 3rd generation detectors.
Notably, these low-frequency signatures complement previous efforts to
identify first-order phase transitions using high-frequency post-merger
gravitational waves (Bauswein et al., 2019), expanding the range of density
regimes that can be probed for potential phase transitions. We provide two
analytical formulas to clarify how the phase shift in the gravitational
waveform relates to the properties of the phase transition. These formulas
will serve as valuable tools for future attempts to directly estimate phase
transition parameters using GW data. Our analysis suggests that the estimation
errors are moderate, with a maximum of $\sim 15\%$ for the mode frequencies
and 25% for the tidal overlap integrals. This implies a potential
underestimation or overestimation of the phase shift by no more than a factor
of 2 to 3, as detailed in the Appendix .3-.4. Even in cases where no mode
resonance signal is detected, the search results can effectively rule out a
significant parameter space, as demonstrated earlier with GW170817.

To calculate the correction to the gravitational waveform induced by the
g-mode resonance, we use the form presented in Eq. (3). However, it is
important to note that this form is derived under the assumption that the
energy transfer of the resonance is instantaneous (Yu & Weinberg, 2017b).
Therefore, it is only valid when the duration of the resonance, $t_{\rm
res}\simeq 0.009s\,(\mathcal{M}/1.2\,{M}_{\odot})^{-5/6}(f/600\,{\rm
Hz})^{-11/6}$, is much shorter than the orbital decay time-scale, $t_{D}\simeq
0.08s\,(\mathcal{M}/1.2\,{M}_{\odot})^{-5/3}(f/600\;\text{Hz})^{-8/3}$ (Lai,
1994). For resonant frequencies that are very high (e.g., $\gtrsim
1000\;\text{Hz}$), the waveform correction $\Delta\Psi(f)$ should be modified
to account for the continuous process of energy transfer during the resonance.
It is also important to consider the relativistic corrections to the Newtonian
orbit as the binary approaches closer, especially for high resonant
frequencies. As demonstrated by Bildsten & Cutler (1992), tidal interaction is
not sufficient to spin up the binary NS system before the merger and an ir-
rotational merger is expected in most cases due to the spinning-down by
magnetic dipole radiation before the merger. At present, we do not take into
consideration the star rotation. Nevertheless, NS spin could be important in
cases such as dynamic captures, which can shift the mode frequency to lower
values and thus can contribute to larger GW phase shift (Ho & Lai, 1999).
These improvements will be the subject of future work.

In addition, it is crucial to identify the physical origin of the mode
resonance signature once it is detected. This is because there could be other
possible origins apart from the resonant excitation of hybrid star g-modes. As
suggested by Pan et al. (2020), other processes such as tidal-p-g coupling
(Weinberg et al., 2013) or those in modified gravity (Palenzuela et al., 2014;
Annulli et al., 2019; Krüger & Doneva, 2021; Mendes et al., 2021) predict
different corrections on the waveform. Therefore, one could perform a Bayesian
model selection to compare these different origins. On the other hand, crust-
core i-mode or core composition g-mode can produce a $\Delta\Psi(f)$ similar
to that studied in this work. For the former, its resonant frequency is
typically around $100\;\text{Hz}$, while for the latter, its frequency can
reach up to several hundred Hz according to recent studies using the Cowling
approximation (see e.g. Jaikumar et al. (2021)). Our Fisher analysis shows
that the resonant frequencies can be measured with a precision of
$\lesssim$10% using 3rd generation detectors (see Appendix .5). This precision
would enable us to discriminate between the discontinuity g-mode resonance
considered here and the crust-core i-mode or composition g-mode by comparing
their frequencies. It is worth mentioning that the existence and properties of
QCD phase transition are still very unclear in nuclear physics studies.
Various scenarios, including cross-over/first-order transitions and fast/slow
transitions, lack effective constraints. The analysis of this paper is based
on the slow first-order phase transition scenarios. We might expect future GW
observations of binary NS mergers could exclude or put constraints on such
scenarios, hence contributing to our understanding of the QCD phase diagram.

###  .1 Solving pulsation equations in general relativity

For the background star, the static and spherically symmetric line element is
given by

$ds^{2}=g_{\mu\nu}dx^{\mu}dx^{\nu}=-e^{2\Phi}dt^{2}+e^{2\lambda}dr^{2}+r^{2}(d\theta^{2}+\sin^{2}\theta
d\phi^{2})$ (1)

where the metric components $\Phi$ and $\lambda$ are functions of $r$ only. We
model the star as a perfect fluid, with the stress-energy-momentum tensor
given by

$T_{\mu\nu}=(p+\varepsilon)u_{\mu}u_{\nu}+pg_{\mu\nu},$ (2)

where $u^{\mu}$ is the four-velocity. We also assume the metric perturbation
as

$h_{\mu\nu}=-r^{l}H_{0}e^{i\omega t}Y_{lm}dt^{2}-r^{l}H_{0}e^{i\omega
t}Y_{lm}dr^{2}-r^{l}Ke^{i\omega t}Y_{lm}r^{2}(d\theta^{2}+\sin^{2}\theta
d\phi^{2})-2i\omega r^{l+1}H_{1}e^{i\omega t}Y_{lm}dtdr,$ (3)

and the Lagrangian displacement as

$\begin{split}\xi^{r}=&r^{l-1}e^{-\lambda}WY_{lm}e^{i\omega t},\\\
\xi^{\theta}=&-r^{l-2}V\partial_{\theta}Y_{lm}e^{i\omega t},\\\
\xi^{\phi}=&-r^{l-2}\sin^{-2}\theta V\partial_{\phi}Y_{lm}e^{i\omega
t},\end{split}$ (4)

where $Y_{lm}(\theta,\phi)$ is the spherical harmonics. Note the Eulerian
perturbations of energy density and pressure are given by

$\begin{array}[]{ll}\delta\varepsilon&=(p+\varepsilon)(\Delta
n/n)-\xi^{i}\frac{d\varepsilon}{dx^{i}},\\\ \delta p&=\gamma p(\Delta
n/n)-\xi^{i}\frac{dp}{dx^{i}}.\end{array}$ (5)

where $\gamma=[(p+\varepsilon)/p](\partial p/\partial\varepsilon)_{\rm ad}$ is
the adiabatic index. The Lagrangian perturbation of baryon number density
$\Delta n$ can be derived from the baryon number conservation
$\nabla_{\mu}(nu^{\mu})=0$, it reads

$\Delta n=-n\nabla_{i}\xi^{i}-\frac{1}{2}n\delta g^{(3)}/g^{(3)},$ (6)

with $g^{(3)}$ the determinant of the metric of the 3-geometry at constant
time. Plugging these perturbations (Eqs. 3-5) into linear Einstein equation
$\delta G_{\mu\nu}=8\pi\delta T_{\mu\nu}$, one finds (Lindblom & Detweiler,
1983; Detweiler & Lindblom, 1985)

$\displaystyle\frac{dH_{1}}{dr}=$
$\displaystyle-\bigg{[}l+1+\frac{2m}{r}e^{2\lambda}+4\pi
r^{2}(p-\varepsilon)e^{2\lambda}\bigg{]}\frac{H_{1}}{r}+\frac{e^{2\lambda}}{r}\bigg{[}H+K-16(\varepsilon+p)V\bigg{]},$
(7) $\displaystyle\frac{dK}{dr}=$
$\displaystyle\frac{l(l+1)}{2r}H_{1}+\frac{H_{0}}{r}-\left(l+1-r\Phi^{\prime}\right)\frac{K}{r}-\frac{8\pi}{r}(\varepsilon+p)e^{\lambda}W,$
(8) $\displaystyle\frac{dW}{dr}=$
$\displaystyle-\frac{l+1}{r}W+re^{\lambda}\bigg{[}\frac{1}{\gamma
p}e^{-\Phi}X-\frac{l(l+1)}{r^{2}}V+\frac{1}{2}H_{0}+K\bigg{]},$ (9)
$\displaystyle\frac{dX}{dr}=$
$\displaystyle-\frac{l}{r}X+(\varepsilon+p)\frac{e^{\Phi}}{r}\Bigg{\\{}\frac{1}{2}(1-r\Phi^{\prime})H_{0}+\frac{1}{2}\left(\omega^{2}r^{2}e^{-2\lambda}+\frac{l(l+1)}{2}\right)H_{1}-\left(\frac{1}{2}-\frac{3}{2}r\Phi^{\prime}\right)K-\frac{l(l+1)}{r}\Phi^{\prime}V$
$\displaystyle-\bigg{(}\omega^{2}e^{\lambda-2\Phi}+4\pi(\varepsilon+p)e^{\lambda}-r^{2}(\frac{1}{r^{2}}e^{-\lambda}\Phi^{\prime})^{\prime}\bigg{)}W\Bigg{\\}},$
(10)

with

$X\equiv\omega^{2}(p+\varepsilon)e^{-\Phi}V-r^{-1}e^{\Phi-\lambda}p^{\prime}W+\frac{1}{2}(p+\varepsilon)e^{\Phi}H_{0}.$
(11)

and

$\begin{split}(\frac{3m}{r}+\frac{(l+2)(l-1)}{2}+4\pi r^{2}p)H_{0}=&8\pi
r^{2}e^{-\Phi}X+r^{2}e^{-2\lambda}\bigg{[}\omega^{2}e^{-2\Phi}-\frac{l(l+1)}{2}\Phi^{\prime}\bigg{]}H_{1}\\\
&+\bigg{[}\frac{(l+2)(l-1)}{2}-\omega^{2}r^{2}e^{-2\Phi}+(r-3m-4\pi
r^{3}p)\Phi^{\prime}\bigg{]}K.\end{split}$ (12)

Here the prime denotes the derivative of $r$.

Figure 4: An illustration of our NS EOS with a first-order phase transition to
quark matter.

Figure 5: Left: Density profile for a $1.4\,{M}_{\odot}$ hybrid star, with
$\varepsilon_{\rm t}=300\,{\rm MeV/fm^{3}}$, $\Delta\varepsilon=40\,{\rm
MeV/fm^{3}}$ and $c_{s}=1$ within QMF+CSS; Right: Perturbation functions of
the discontinuity g-mode corresponding to the star shown in the left panel.
The unit of metric perturbation in the right panel is $153\,{\rm MeV/fm^{3}}$.

Quark matter is present from the first-order transition onset up to the
maximum central pressure of a star (see Fig. 4). Our description of the EOS is
written as follows:

$\varepsilon(p)=\left\\{\\!\begin{aligned} &\varepsilon_{\rm HM}(p),&p<p_{\rm
t}\\\ &\varepsilon_{\rm HM}(p_{\rm t})+\Delta\varepsilon+c_{s}^{-2}(p-p_{\rm
t}),&p>p_{\rm t}\end{aligned}\right.$ (13)

where ‘HM’ denotes hadronic matter. We then follow the method described in
Lindblom & Detweiler (1983) to solve these equations and find the eigenvalue
of $\omega$ for each mode. It’s particularly important to note that when
solving the aforementioned perturbations equations, we require the boundary
condition that $\Delta p=-r^{l}e^{-\Phi}XY_{lm}e^{i\omega t}$ to be continuous
at the hadron-quark interface. This necessitates that at the transition radius
$R_{\rm t}$,

$V_{+}=\frac{\varepsilon_{\rm t}+\Delta\varepsilon+p}{\varepsilon_{\rm
t}+p}V_{-}+\frac{\Delta\varepsilon}{\omega^{2}(\varepsilon_{\rm
t}+p)}\left(\frac{e^{2\Phi-\lambda}}{R_{\rm
t}}\Phi^{\prime}W+\frac{1}{2}e^{2\Phi}H_{0}\right)$ (14)

where the subscripts ‘$+$’ and ‘$-$’ denote the values at $R_{\rm t}+0^{+}$
and $R_{\rm t}-0^{+}$, respectively. This implies that at both ends of the
interface, although the radial displacement W is continuous, the tangential
displacement V is discontinuous due to the presence of a density
discontinuity. It is this particularity that gives rise to the discontinuity
g-mode studied in this work. Fig. 5 provides a solved example, allowing us to
observe the density jump and the corresponding variation of $V$ at the
interface.

###  .2 Relativistic tidal overlap integral

Table 1: Comparison for the numerical tidal overlap $Q$ and g-mode overlap
$\langle\xi_{g}|\xi_{\alpha}\rangle$ obtained by using different approaches.
The EOS model is QMF+CSS with $\varepsilon_{\rm t}=300\,{\rm MeV/fm^{3}}$ and
$\Delta\varepsilon=40\,{\rm MeV/fm^{3}}$.

$\alpha$ | Hybrid | GR | GR
---|---|---|---
(Newtonian inner product) | (relativistic inner product)
$\langle\xi_{g}|\xi_{\alpha}\rangle$ | $|Q_{\alpha}^{N}|$ | $\langle\xi_{g}|\xi_{\alpha}\rangle$ | $|Q_{\alpha}^{N}|$ | $\langle\xi_{g}|\xi_{\alpha}\rangle$ | $|Q_{\alpha}|$
g | 1 | 3.73e-2 | 1 | 4.17e-2 | 1 | 3.45e-2
f | 7.68e-2 | 5.92e-1 | 9.79e-3 | 5.86e-1 | 7.05e-5 | 7.08e-1
p1 | -2.71e-2 | 1.17e-1 | 3.95e-3 | 4.22e-2 | 1.65e-5 | 1.01e-2
p2 | 2.60e-2 | 5.25e-2 | -1.41e-3 | 8.88e-3 | -2.50e-5 | 2.69e-2
p3 | -1.82e-2 | 2.12e-3 | 2.21e-3 | 1.29e-2 | -2.95e-5 | 4.00e-2
$\sum_{\alpha}Q_{\alpha}^{2}$ | 0.3680 | 0.3073 | 0.5150
$\langle\nabla(r^{l}Y_{lm})|\nabla(r^{l}Y_{lm})\rangle$ | 0.3504 | 0.3504 | 0.5175

In some previous works, the hybrid (GR background + Newtonian perturbation)
approach was typically used to calculate oscillation modes and the Newtonian
formalism was employed to compute tidal overlap integral, i.e.,
$Q^{N}_{\alpha}\equiv\langle\xi_{\alpha}|\nabla(r^{l}Y_{lm})\rangle/\langle\xi_{\alpha}|\xi_{\alpha}\rangle$
with the inner product
$\langle\xi|\psi\rangle=\int\rho\boldsymbol{\xi}^{*}\cdot\boldsymbol{\psi}d^{3}x$.
This hybrid approach is inconsistent and will result in computed
eigenfunctions that are not orthogonal (Reisenegger, 1994), means
$\langle\xi_{\alpha}|\xi_{\beta}\rangle\neq 0$ when $\alpha\neq\beta$, see
Table 1. In this case, the obtained $Q^{N}$ for g-modes would contain some
f-mode contamination, leading to unreliable values of $Q^{N}$.

In this work, we calculate the stellar oscillation modes within the GR
framework. This naturally ensures that the eigenfunctions are orthogonal under
a relativistic extension of inner product 111Strictly speaking, the definition
of this inner product should also include the integration of metric
perturbations out to infinity (Detweiler & Ipser, 1973). However, this term is
small for g-modes considered in this work and we omit it here.

$\langle\xi^{\mu}|\psi^{\mu}\rangle=\int\sqrt{-g}e^{-2\Phi}(p+\varepsilon)\xi^{*\mu}\psi_{\mu}d^{3}x.$
(15)

Based on this inner product, we find that the relativistic tidal overlap
integral can still be defined as

$Q_{\alpha}=\frac{\langle\xi_{\alpha}|\nabla(r^{l}Y_{lm})\rangle}{\langle\xi_{\alpha}|\xi_{\alpha}\rangle}.$
(16)

In the following, we give a brief derivation and test the accuracy of our
numerical results.

We begin from the orthogonal part of the Euler equation
$\nabla_{\nu}T^{\mu\nu}=0$, which is expressed as

$(p+\varepsilon)u^{\nu}\nabla_{\nu}u^{\mu}+\perp^{\mu\nu}\nabla_{\nu}p=0,$
(17)

where $\perp^{\mu\nu}=g^{\mu\nu}+u^{\mu}u^{\nu}$ is the projection operator.
Plugging the the Lagrangian displacement of a certain mode $\alpha$,
$\xi_{\alpha}^{\mu}=(0,\boldsymbol{\xi}_{\alpha}(r)e^{i\omega_{\alpha}t})$,
into Eq. (17) one can derive the linear pulsation equation

$[\mathcal{L}-(p+\varepsilon)e^{-2\Phi}\omega_{\alpha}^{2}]\xi_{\alpha}^{\mu}=0,$
(18)

with

$\mathcal{L}\xi_{\alpha}^{\mu}=(\delta
p+\delta\varepsilon)u^{\nu}\nabla_{\nu}u^{\mu}+\perp^{\mu\nu}\nabla_{\nu}\delta
p+i\omega_{\alpha}e^{-\Phi}\left[(p+\varepsilon)\nabla_{\nu}u^{\mu}+u^{\mu}\nabla_{\nu}p\right]\xi_{\alpha}^{\nu}+i\omega_{\alpha}e^{-\Phi}\xi_{\alpha}^{\mu}u^{\nu}\nabla_{\nu}p.$
(19)

During the inspiral, the linear perturbation of the tidal potential is
governed by

$(p+\varepsilon)u^{\nu}\partial_{\nu}\delta
u^{\mu}+(p+\varepsilon)u^{\nu}\delta\Gamma^{\mu}_{\nu\rho}u^{\rho}+\mathcal{L}\xi^{\mu}=0.$
(20)

Plugging Eq. (18), $u^{\mu}=(e^{-\Phi},0,0,0)$, $\delta
u^{i}=e^{-\Phi}\partial_{t}\xi^{i}$ and
$\delta\Gamma^{\mu}_{00}=-1/2g^{\mu\rho}\partial_{\rho}h^{\rm tid}_{00}$ into
Eq. (20), one obtains

$(p+\varepsilon)e^{-2\Phi}\frac{\partial^{2}}{\partial
t^{2}}\xi^{\mu}+\mathcal{L}\xi^{\mu}=-(p+\varepsilon)e^{-2\Phi}\partial^{\mu}U,$
(21)

where $U=-1/2h^{\rm tid}_{00}=-M^{\prime}/|\boldsymbol{r}-\boldsymbol{D}|$ is
the tidal potential, $D$ is the distance of the secondary component
$M^{\prime}$. Eq. (21) is similar to its Newtonian counterpart, as seen, for
instance, in Eq.(2.5) of Lai (1994). Hence, similar to Lai (1994), we arrive
to the relativistic definition of tidal overlap integral as
$Q_{\alpha}=\langle\xi_{\alpha}|\nabla(r^{l}Y_{lm})\rangle/\langle\xi_{\alpha}|\xi_{\alpha}\rangle$.
In numerical calculation, we choose the normalization
$\langle\xi_{\alpha}|\xi_{\alpha}\rangle=MR^{2}$ for convenience.

To test the accuracy of our numerical results, we first calculate the inner
product $\langle\xi_{g}|\xi_{\alpha}\rangle$, the exemplary results are shown
in Table 1. It can be seen that for $\alpha\neq g$, the inner products we
obtained are very tiny, which verifies the orthogonality of the
eigenfunctions. We also check the sum rule of $Q_{\alpha}$, namely
$\sum_{\alpha}Q_{\alpha}^{2}=\langle\nabla(r^{l}Y_{lm})|\nabla(r^{l}Y_{lm})\rangle$,
which is a general property of the tidal overlap integral (Reisenegger, 1994).
From Table 1 we see that our numerical $Q$ obeys the sum rule.

###  .3 Approximated formulas of the mode frequency

Previously, Finn (1987) pointed out that when the discontinuity near the
surface (i.e., $R-R_{\rm t}\ll R_{\rm t}$), the frequency can be approximated
as

$\omega^{2}\approx l(l+1)\frac{\Delta\varepsilon}{\varepsilon_{\rm
t}+\Delta\varepsilon}\frac{R-R_{\rm t}}{R}\frac{M}{R^{3}}.$ (22)

And this approximation has been improved by Zhao & Lattimer (2022) with

$\omega^{2}\approx\frac{D\Delta\varepsilon/\varepsilon_{\rm
t}}{(1+\Delta\varepsilon/\varepsilon_{\rm t})/\tanh[D]+1/\tanh[D(R/R_{\rm
t}-1)]}\frac{M_{\rm t}}{R_{\rm t}^{3}}$ (23)

where $M_{\rm t}$ is the mass inside $R_{\rm t}$ and $D=1.21$ is a fitting
coefficient.

Figure 6: Estimated frequencies from Eq. (24) versus the frequencies
calculated numerically from relativistic pulsation equations.

However, both approximations of Finn (1987) and Zhao & Lattimer (2022) start
from the analytical solution of the frequency of surface gravity waves between
two incompressible fluids, obtained by Landau & Lifshitz (1959) assuming the
slab configuration. Thus it should be noted that their approximations can not
seamlessly cater to both the $R-R_{\rm t}\ll R$ and $R_{\rm t}\ll R$ cases,
see Fig. 13 in Zhao & Lattimer (2022). In this work, we present a more
accurate approximation that can accommodate both these two cases
simultaneously, i.e.,

$\omega^{2}\approx K\frac{(6/5)(\Delta\varepsilon/\varepsilon_{\rm
t})(1-x^{5})+1+\varsigma-\sqrt{(1+\varsigma)^{2}-(12/5)(\Delta\varepsilon/\varepsilon_{\rm
t})(\kappa-1)(1-x^{5})}}{(\varepsilon_{\rm
t}+\Delta\varepsilon)/\varepsilon_{\rm
t}-(2/5)(\Delta\varepsilon/\varepsilon_{\rm t})(1-x^{5})}\frac{M}{R^{3}},$
(24)

where $\varsigma=(1/5)(\Delta\varepsilon/\varepsilon_{\rm
t})[3(\kappa-1)+(3+2\kappa)x^{5}]$ with $x=R_{\rm t}/R$ and
$\kappa=(M_{t}/M)(R/R_{\rm t})^{3}$. $K=0.59$ is a fitting coefficient. Fig. 6
shows the approximation of Eq. (24). It can be seen that the estimation of Eq.
(24) performs well for both $R_{\rm t}<R/2$ and $R_{\rm t}>R/2$ regions. The
estimation errors of Eq. (24) are found to be not more than $\sim 15\%$.

We should emphasize here that, apart from a factor $K$, Eq. (24) is the exact
solution for a two-layer fluid in a spherical configuration within the
Newtonian limit. In the following we present a brief derivation.

We begin with the pulsation equation under the Cowling approximation (Cowling,
1941)

$\begin{split}\frac{d\xi_{r}}{dr}&=\left(\frac{gr}{c_{ad}^{2}}-2\right)\frac{\xi_{r}}{r}+\left[l(l+1)-\frac{\omega^{2}r^{2}}{c_{ad}^{2}}\right]\frac{\xi_{h}}{r},\\\
\frac{d\xi_{h}}{dr}&=\left(1-\frac{N^{2}}{\omega^{2}}\right)\frac{\xi_{r}}{r}+\left(\frac{N^{2}r}{g}-1\right)\frac{\xi_{h}}{r},\end{split}$
(25)

where $\xi_{r}$ and $\xi_{h}$ are the radial and tangential displacement,
respectively. $g$ is the gravitational acceleration.
$N^{2}=g^{2}(1/c_{e}^{2}-1/c_{ad}^{2})$ is the Brunt-Väisälä frequency with
$c_{e}$ and $c_{ad}$ the equilibrium and adiabatic sound speed, respectively.
For simplicity, we shall omit the gravity and composition gradient, i.e., we
assume $c_{ad}=c_{e}$, $gr/c_{ad}^{2}=0$ and $\omega^{2}r^{2}/c_{ad}^{2}=0$.
The density for the lower layer (i.e., quark core) and upper layer (i.e.,
hadronic envelope) are denoted as $\rho^{-}$ and $\rho^{+}$. Restrict to
$l=2$, one finds that Eq. (25) reduces to

$\begin{split}\frac{d\xi_{r}}{dr}&=-\frac{2\xi_{r}}{r}+\frac{6\xi_{h}}{r},\\\
\frac{d\xi_{h}}{dr}&=\frac{\xi_{r}}{r}-\frac{\xi_{h}}{r},\end{split}$ (26)

The general solutions can be obtained with

$\xi_{r}=\left\\{\\!\begin{aligned} &r,&0<r<R_{\rm t},\\\ &\frac{\alpha R_{\rm
t}^{4}}{r^{4}}+\frac{\beta r}{R_{\rm t}},&R_{\rm t}<r<R.\end{aligned}\right.$
(27)

and

$\xi_{h}=\left\\{\\!\begin{aligned} &\frac{1}{2}r,&0<r<R_{\rm t},\\\
&-\frac{1}{3}\frac{\alpha R_{\rm t}^{4}}{r^{4}}+\frac{1}{2}\frac{\beta
r}{R_{\rm t}},&R_{\rm t}<r<R.\end{aligned}\right.$ (28)

At the interface of two layers, $\xi_{r}$ and $\Delta p=\rho
g(\omega^{2}r\xi_{h}-\xi_{r})$ should be continuous, which requires

$\begin{split}\alpha&=\frac{\Delta\rho}{\rho^{+}}\left(\frac{6}{5}\frac{M_{t}}{R_{\rm
t}^{3}\omega^{2}}-\frac{3}{5}\right)R_{\rm t},\\\ \beta&=R_{\rm
t}-\alpha,\end{split}$ (29)

where $\Delta\rho\equiv\rho^{-}-\rho^{+}$. At the surface, the Lagrangian
perturbation of the pressure must vanish, i.e., $\Delta p=0$. Using Eq.
(27-29) and replacing $\rho^{+}$ ($\rho^{-}$) with $\varepsilon_{\rm t}$
($\varepsilon_{\rm t}+\Delta\varepsilon$), we finally obtain

$\omega^{2}=\frac{(6/5)(\Delta\varepsilon/\varepsilon_{\rm
t})(1-x^{5})+1+\varsigma-\sqrt{(1+\varsigma)^{2}-(12/5)(\Delta\varepsilon/\varepsilon_{\rm
t})(\kappa-1)(1-x^{5})}}{(\varepsilon_{\rm
t}+\Delta\varepsilon)/\varepsilon_{\rm
t}-(2/5)(\Delta\varepsilon/\varepsilon_{\rm t})(1-x^{5})}\frac{M}{R^{3}}.$
(30)

###  .4 Approximated formulas of the tidal overlap integral

We find a reliable approximation for the tidal overlap integral, as shown in
Fig. 7. We should mention that such an approximation is only valid for small
transition strength, i.e., $\Delta\varepsilon\lesssim 150\,{\rm MeV/fm^{3}}$.
However, higher transition strength ($\Delta\varepsilon\gtrsim 150\,{\rm
MeV/fm^{3}}$) corresponds to higher mode frequencies ($\gtrsim
1200\;\text{Hz}$) that we are not interested in, therefore we do not need to
worry about them when using this approximation. By defining $\tilde{Q}\equiv
Q[(\varepsilon_{\rm t}+\Delta\varepsilon)(\varepsilon_{\rm t}+p_{\rm
t})/\varepsilon_{\rm t}^{2}]^{-1}M_{1.4}^{-2}$, we find

$\tilde{Q}\simeq\frac{cx(1-dx)}{e^{a(1-x)}+b}.$ (31)

We fit the data points of the $1.4\,M_{\odot}$ stars and the resulting fitting
coefficients are presented in Table 2. It is found the estimation errors are
less than $\sim 25\%$.

Table 2: Fit parameters of the tidal overlap integral for $1.4\,M_{\odot}$
stars in the functional form of Eq. (31).

Model | a | b | c | d
---|---|---|---|---
NL3$\omega\rho$ (stiff) +CSS | 6.304 | 17.793 | 2.827 | 1.004
QMF (soft)+CSS | 5.726 | 12.014 | 3.183 | 1.064

Figure 7: Scaled tidal overlap integral $\tilde{Q}$ versus $x$, together with
the fitting lines for $1.4\,M_{\odot}$ stars. The associated fitting
coefficients can be found in Table 2.

###  .5 Detectability analysis

Following Cutler & Flanagan (1994); Yu & Weinberg (2017b), we use the Fisher
information matrix to estimate the detectability of the phase shift due to
resonance, assuming that the GW signal have a large signal-to-noise ratio
(SNR). The Fisher information matrix for a set of parameters
$\\{\theta_{i}\\}$ is defined by

$\Gamma_{ij}=\left(\frac{\partial h}{\partial\theta_{i}}|\frac{\partial
h}{\partial\theta_{i}}\right),$ (32)

where the inner product is defined by

$(h_{1}|h_{2})=2\int_{0}^{\infty}\frac{\tilde{h}_{1}^{*}(f)\tilde{h}_{2}(f)+\tilde{h}_{1}(f)\tilde{h}_{2}^{*}(f)}{S_{n}(f)}df,$
(33)

with $\tilde{h}$ and $S_{n}$ are the Fourier transform of GW strain data and
the spectral noise density of the detector. The SNR is given by ${\rm
SNR}=(h|h)^{1/2}$. The root-mean-square (rms) measurement error of
$\theta_{i}$ can be written as

$\Delta\theta_{i}=\sqrt{(\Gamma^{-1})_{ii}},$ (34)

where $\Gamma^{-1}$ is the inverse of the Fisher information matrix. If the
rms error of phase shift $\Delta|\delta\bar{\phi}|$is smaller than
$|\delta\bar{\phi}|$, then the phase shift is detectable.

We employ the IMRPhenomD_NTidal waveform (Dietrich et al., 2017), together
with the addition of phase correction given by Eq. (5), to calculate the GW
signal. The waveform parameters are

$\\{\theta_{i}\\}=\\{\mathcal{M},q,\Lambda_{1},\Lambda_{2},d_{L},t_{c},\phi_{c},|\delta\bar{\phi}|,\bar{f}\\}.$
(35)

Here for simplicity, we omit the angular dependence and the stellar rotation.
We also fix the tidal defromability as $\Lambda_{1.4}=400$ and
$\Lambda_{1.8}=300$, because it is found that $\Delta|\delta\bar{\phi}|$ is
insensitive to the choice of tidal deformability. In Table 3 we list some
quantitative results of rms measurement errors for a
$1.4\,M_{\odot}-1.4\,M_{\odot}$ system.

Table 3: The rms measurement errors of $|\delta\bar{\phi}|$ and $\bar{f}$ for
a $1.4\,{M}_{\odot}-1.4\,{M}_{\odot}$ system.

$\bar{f}\,({\rm Hz})$ | $|\delta\bar{\phi}|$ | $d_{L}\,({\rm Mpc})$ | Detector | SNR | $\Delta\bar{f}\,({\rm Hz})$ | $\Delta|\delta\bar{\phi}|$
---|---|---|---|---|---|---
500 | 1 | 100 | aLIGO | 49 | 140 | 1.497
ET | 502 | 18 | 0.202
CE | 1495 | 11 | 0.127
750 | 1 | 100 | aLIGO | 49 | 359 | 2.842
ET | 502 | 57 | 0.336
CE | 1495 | 36 | 0.217
1000 | 1 | 100 | aLIGO | 49 | 885 | 4.060
ET | 502 | 145 | 0.664
CE | 1495 | 93 | 0.429

###  .6 Bayesian parameter estimation for GW170817

We perform a Bayesian parameter estimation to analysis the GW170817 data by
incorporating the resonant parameter $|\delta\bar{\phi}|$ and $\bar{f}$. We
utilize the IMRPhenomD_NTidal waveform for the basic hypothesis
$\mathcal{H}_{0}$, whereas add the phase correction for the hypothesis
$\mathcal{H}_{1}$. The waveform parameters are given by

$\\{\theta_{i}\\}=\\{\mathcal{M},q,\Lambda_{1},\Lambda_{2},{\chi}_{1z},{\chi}_{2z},\theta_{\rm
jn},t_{c},\phi_{c},\Psi,|\delta\bar{\phi}|,\bar{f}\\}$ (36)

We fix the location of the source to the position determined by
electromagnetic observations (Abbott et al., 2017c; Levan et al., 2017) with
$\alpha(J2000)=197.45^{\circ}$, $\delta(J2000)=-23.38^{\circ}$ and $z=0.0099$.
The priors of the parameters are chosen by following those used in Abbott et
al. (2019), with the exception of the priors for $|\delta\bar{\phi}|$ and
$\bar{f}$. For the mode resonance parameters, we’ve chosen priors as
$|\delta\bar{\phi}|\sim{\rm LogU}[10^{-4},10^{2}]$ and $\bar{f}\sim{\rm
U}[50,1600]\;\text{Hz}$, where ‘U’ denotes the uniform distribution and ‘LogU’
represents the logarithmic uniform distribution.

In Fig. 8 we show the posterior distribution of waveform parameters. Since we
use the PyMultiNest to sample from the posterior distribution, we can obtain
the evidence ($\mathcal{Z}$) of each hypothesis directly. We then calculate
their Bayes Factor with $\mathcal{B}_{0}^{1}=\mathcal{Z}_{1}/\mathcal{Z}_{0}$.
We also calculate the Bayes Factor with the Savage-Dickey Density Ratio method
(Dickey & Lientz, 1970), by using the posterior distribution of hypothesis
$\mathcal{H}_{1}$ alone, as

$\frac{1}{\mathcal{B}_{0}^{1}}=\lim_{|\delta\bar{\phi}|\to
0}\frac{p(|\delta\bar{\phi}||data,\mathcal{H}_{1})}{p(|\delta\bar{\phi}||\mathcal{H}_{1})}.$
(37)

Both methods yield similar Bayes Factors. By evaluating the evidence, we
obtain $\mathcal{B}_{0}^{1}=0.72$, while using the Savage-Dickey Density Ratio
method, we find $\mathcal{B}_{0}^{1}=1.11$. According to the significance
quantification approach in Kass & Raftery (1995), when the Bayes Factor is
close to unity, it indicates that the current data does not strongly favor one
hypothesis over the other.

Figure 8: Posterior distributions of the waveform parameters inferred from
GW170817 strain data. The horizontal lines in last two diagonal panels denote
the prior distributions of corresponding parameters.

## Acknowledgments

We thank Dong Lai, Lap-Ming Lin, Yuxin Liu, Zhen Pan, Yong Gao and Zhenyu Zhu
for their helpful comments. The work is supported by the National SKA Program
of China (No. 2020SKA0120300) and the National Natural Science Foundation of
China (grant Nos. 12273028 and 12203017). This research has made use of data
and software obtained from the Gravitational Wave Open Science Center
$https://www.gw-openscience.org$, a service of LIGO Laboratory, the LIGO
Scientific Collaboration and the Virgo Collaboration. LIGO is funded by the
U.S. National Science Foundation. Virgo is funded by the French Centre
National de Recherche Scientifique (CNRS), the Italian Istituto Nazionale
della Fisica Nucleare (INFN) and the Dutch Nikhef, with contributions by
Polish and Hungarian institutes.

## References

  * Abbott et al. (2017a) Abbott, B. P., Abbott, R., Abbott, T. D., et al. 2017a, Phys. Rev. Lett., 119, 161101, doi: 10.1103/PhysRevLett.119.161101
  * Abbott et al. (2017b) —. 2017b, Classical and Quantum Gravity, 34, 044001, doi: 10.1088/1361-6382/aa51f4
  * Abbott et al. (2017c) —. 2017c, ApJ, 848, L12, doi: 10.3847/2041-8213/aa91c9
  * Abbott et al. (2018) —. 2018, Phys. Rev. Lett., 121, 161101, doi: 10.1103/PhysRevLett.121.161101
  * Abbott et al. (2019) —. 2019, Physical Review X, 9, 011001, doi: 10.1103/PhysRevX.9.011001
  * Alford et al. (2015) Alford, M. G., Burgio, G. F., Han, S., Taranto, G., & Zappalà, D. 2015, Phys. Rev. D, 92, 083002, doi: 10.1103/PhysRevD.92.083002
  * Alford et al. (2013) Alford, M. G., Han, S., & Prakash, M. 2013, Phys. Rev. D, 88, 083013, doi: 10.1103/PhysRevD.88.083013
  * Annulli et al. (2019) Annulli, L., Cardoso, V., & Gualtieri, L. 2019, Phys. Rev. D, 99, 044038, doi: 10.1103/PhysRevD.99.044038
  * Ashton et al. (2019) Ashton, G., Hübner, M., Lasky, P. D., et al. 2019, ApJS, 241, 27, doi: 10.3847/1538-4365/ab06fc
  * Bauswein et al. (2019) Bauswein, A., Bastian, N.-U. F., Blaschke, D. B., et al. 2019, Phys. Rev. Lett., 122, 061102, doi: 10.1103/PhysRevLett.122.061102
  * Bildsten & Cutler (1992) Bildsten, L., & Cutler, C. 1992, ApJ, 400, 175, doi: 10.1086/171983
  * Buchner et al. (2014) Buchner, J., Georgakakis, A., Nandra, K., et al. 2014, A&A, 564, A125, doi: 10.1051/0004-6361/201322971
  * Cai et al. (2022) Cai, R.-G., He, S., Li, L., & Wang, Y.-X. 2022, Phys. Rev. D, 106, L121902, doi: 10.1103/PhysRevD.106.L121902
  * Cheng & Dai (1996) Cheng, K. S., & Dai, Z. G. 1996, Phys. Rev. Lett., 77, 1210, doi: 10.1103/PhysRevLett.77.1210
  * Cowling (1941) Cowling, T. G. 1941, MNRAS, 101, 367, doi: 10.1093/mnras/101.8.367
  * Cutler & Flanagan (1994) Cutler, C., & Flanagan, É. E. 1994, Phys. Rev. D, 49, 2658, doi: 10.1103/PhysRevD.49.2658
  * Detweiler & Lindblom (1985) Detweiler, S., & Lindblom, L. 1985, ApJ, 292, 12, doi: 10.1086/163127
  * Detweiler & Ipser (1973) Detweiler, S. L., & Ipser, J. R. 1973, ApJ, 185, 685, doi: 10.1086/152447
  * Dickey & Lientz (1970) Dickey, J. M., & Lientz, B. P. 1970, The Annals of Mathematical Statistics, 41, 214.  http://www.jstor.org/stable/2239734
  * Dietrich et al. (2017) Dietrich, T., Bernuzzi, S., & Tichy, W. 2017, Phys. Rev. D, 96, 121501, doi: 10.1103/PhysRevD.96.121501
  * Finn (1987) Finn, L. S. 1987, MNRAS, 227, 265, doi: 10.1093/mnras/227.2.265
  * Fischer et al. (2018) Fischer, T., Bastian, N.-U. F., Wu, M.-R., et al. 2018, Nature Astronomy, 2, 980, doi: 10.1038/s41550-018-0583-0
  * Heiselberg & Pandharipande (2000) Heiselberg, H., & Pandharipande, V. 2000, Annual Review of Nuclear and Particle Science, 50, 481, doi: 10.1146/annurev.nucl.50.1.481
  * Ho & Lai (1999) Ho, W. C. G., & Lai, D. 1999, MNRAS, 308, 153, doi: 10.1046/j.1365-8711.1999.02703.x
  * Hogan (1983) Hogan, C. J. 1983, Phys. Rev. Lett., 51, 1488, doi: 10.1103/PhysRevLett.51.1488
  * Horowitz & Piekarewicz (2001) Horowitz, C. J., & Piekarewicz, J. 2001, Phys. Rev. Lett., 86, 5647, doi: 10.1103/PhysRevLett.86.5647
  * Huber (2020) Huber, M. Q. 2020, Phys. Rep., 879, 1, doi: 10.1016/j.physrep.2020.04.004
  * Huth et al. (2022) Huth, S., Pang, P. T. H., Tews, I., et al. 2022, Nature, 606, 276, doi: 10.1038/s41586-022-04750-w
  * Jaikumar et al. (2021) Jaikumar, P., Semposki, A., Prakash, M., & Constantinou, C. 2021, Phys. Rev. D, 103, 123009, doi: 10.1103/PhysRevD.103.123009
  * Kass & Raftery (1995) Kass, R. E., & Raftery, A. E. 1995, Journal of the American Statistical Association, 90, 773.  http://www.jstor.org/stable/2291091
  * Kokkotas & Schmidt (1999) Kokkotas, K. D., & Schmidt, B. G. 1999, Living Reviews in Relativity, 2, 2, doi: 10.12942/lrr-1999-2
  * Krüger & Doneva (2021) Krüger, C. J., & Doneva, D. D. 2021, Phys. Rev. D, 103, 124034, doi: 10.1103/PhysRevD.103.124034
  * Kurkela et al. (2010) Kurkela, A., Romatschke, P., & Vuorinen, A. 2010, Phys. Rev. D, 81, 105021, doi: 10.1103/PhysRevD.81.105021
  * Lai (1994) Lai, D. 1994, MNRAS, 270, 611, doi: 10.1093/mnras/270.3.611
  * Landau & Lifshitz (1959) Landau, L. D., & Lifshitz, E. M. 1959, Fluid mechanics (Oxford: Pergamon Press)
  * Lattimer & Prakash (2001) Lattimer, J. M., & Prakash, M. 2001, ApJ, 550, 426, doi: 10.1086/319702
  * Lau & Yagi (2021) Lau, S. Y., & Yagi, K. 2021, Phys. Rev. D, 103, 063015, doi: 10.1103/PhysRevD.103.063015
  * Levan et al. (2017) Levan, A. J., Lyman, J. D., Tanvir, N. R., et al. 2017, ApJ, 848, L28, doi: 10.3847/2041-8213/aa905f
  * Li et al. (2023) Li, A., Yong, G.-C., & Zhang, Y.-X. 2023, Phys. Rev. D, 107, 043005, doi: 10.1103/PhysRevD.107.043005
  * LIGO Scientific Collaboration et al. (2015) LIGO Scientific Collaboration, Aasi, J., Abbott, B. P., et al. 2015, Classical and Quantum Gravity, 32, 074001, doi: 10.1088/0264-9381/32/7/074001
  * Lindblom & Detweiler (1983) Lindblom, L., & Detweiler, S. L. 1983, ApJS, 53, 73, doi: 10.1086/190884
  * Maggiore et al. (2020) Maggiore, M., Van Den Broeck, C., Bartolo, N., et al. 2020, J. Cosmology Astropart. Phys, 2020, 050, doi: 10.1088/1475-7516/2020/03/050
  * Mendes et al. (2021) Mendes, R. F. P., Ortiz, N., & Stergioulas, N. 2021, Phys. Rev. D, 104, 104036, doi: 10.1103/PhysRevD.104.104036
  * Miao et al. (2020) Miao, Z., Li, A., Zhu, Z., & Han, S. 2020, ApJ, 904, 103, doi: 10.3847/1538-4357/abbd41
  * Palenzuela et al. (2014) Palenzuela, C., Barausse, E., Ponce, M., & Lehner, L. 2014, Phys. Rev. D, 89, 044024, doi: 10.1103/PhysRevD.89.044024
  * Pan et al. (2020) Pan, Z., Lyu, Z., Bonga, B., Ortiz, N., & Yang, H. 2020, Phys. Rev. Lett., 125, 201102, doi: 10.1103/PhysRevLett.125.201102
  * Reisenegger (1994) Reisenegger, A. 1994, ApJ, 432, 296, doi: 10.1086/174569
  * Reisenegger & Goldreich (1994) Reisenegger, A., & Goldreich, P. 1994, ApJ, 426, 688, doi: 10.1086/174105
  * Reitze et al. (2019) Reitze, D., Adhikari, R. X., Ballmer, S., et al. 2019, in Bulletin of the American Astronomical Society, Vol. 51, 35, doi: 10.48550/arXiv.1907.04833
  * Tonetto & Lugones (2020) Tonetto, L., & Lugones, G. 2020, Phys. Rev. D, 101, 123029, doi: 10.1103/PhysRevD.101.123029
  * Tsang et al. (2012) Tsang, D., Read, J. S., Hinderer, T., Piro, A. L., & Bondarescu, R. 2012, Phys. Rev. Lett., 108, 011102, doi: 10.1103/PhysRevLett.108.011102
  * Weber (1999) Weber, F. 1999, Pulsars as astrophysical laboratories for nuclear and particle physics (Routledge), doi: 10.1201/9780203741719
  * Weinberg et al. (2013) Weinberg, N. N., Arras, P., & Burkart, J. 2013, ApJ, 769, 121, doi: 10.1088/0004-637X/769/2/121
  * Xu & Lai (2017) Xu, W., & Lai, D. 2017, Phys. Rev. D, 96, 083005, doi: 10.1103/PhysRevD.96.083005
  * Yu & Weinberg (2017a) Yu, H., & Weinberg, N. N. 2017a, MNRAS, 464, 2622, doi: 10.1093/mnras/stw2552
  * Yu & Weinberg (2017b) —. 2017b, MNRAS, 470, 350, doi: 10.1093/mnras/stx1188
  * Zhao & Lattimer (2022) Zhao, T., & Lattimer, J. M. 2022, Phys. Rev. D, 106, 123002, doi: 10.1103/PhysRevD.106.123002
  * Zhu et al. (2023) Zhu, J., Wang, C., Xia, C., Zhou, E., & Ma, Y. 2023, Phys. Rev. D, 107, 083023, doi: 10.1103/PhysRevD.107.083023
  * Zhu et al. (2018) Zhu, Z.-Y., Zhou, E.-P., & Li, A. 2018, ApJ, 862, 98, doi: 10.3847/1538-4357/aacc28