# A note on $LP$-Kenmotsu manifolds admitting Ricci-Yamabe solitons

Mobin Ahmad, Gazala and Mohd Bilal

###### Abstract.

In the current note, we study Lorentzian para-Kenmotsu (in brief,
$LP$-Kenmotsu) manifolds admitting Ricci-Yamabe solitons (RYS) and gradient
Ricci-Yamabe soliton (gradient RYS). At last by constructing a 5-dimensional
non-trivial example we illustrate our result.

2010 Mathematics Subject Classification. 53C20, 53C21, 53C25, 53E20.

Keywords. Lorentzian para-Kenmotsu manifolds, Ricci-Yamabe solitons, Einstein
manifolds, $\nu$-Einstein manifolds.

##  1\. Introduction

In 2019, a scalar combination of Ricci and Yamabe flows was proposed by the
authors Güler and Crasmareanu [6], this advanced class of geometric flows
called Ricci-Yamabe (RY) flow of type $(\sigma,\rho)$ and is defined by

$\displaystyle\frac{\partial}{\partial t}g(t)+2\sigma S(g(t))+\rho
r(t)g(t)=0,\ \quad g(0)=g_{0}$

for some scalars $\sigma$ and $\rho$.

A solution to the RY flow is called RYS if it depends only on one parameter
group of diffeomorphism and scaling. A Riemannian (or semi-Riemannian)
manifold $M$ is said to have a RYS if

(1.1) $\displaystyle\pounds_{K}g+2\sigma S+(2\Lambda-\rho r)g=0,$

where $\sigma,\rho,\Lambda\in\mathbb{R}$ (the set of real numbers). If $K$ is
the gradient of a smooth function $v$ on $M$, then (1.1) is called the
gradient Ricci-Yamabe soliton (gradient RYS) and hence (1.1) turns to

(1.2) $\displaystyle\nabla^{2}v+\sigma S+(\Lambda-\frac{\rho r}{2})g=0,$

where $\nabla^{2}v$ is the Hessian of $v$. It is to be noted that a RYS of
types $(\sigma,0)$ and $(0,\rho)$ are known as $\sigma-$Ricci soliton and
$\rho-$Yamabe soliton, respectively. A RYS is said to be shrinking , steady or
expanding if $\Lambda<0,=0$ or $>0$, respectively. A RYS is said to be a
$\bullet$ Ricci soliton [7] if $\sigma=1,\rho=0$,
$\bullet$ Yamabe soliton [8] if $\sigma=0,\rho=1$,
$\bullet$ Einstein soliton [3] if $\sigma=1,\rho=-1$,

As a continuation of this study, we tried to study RYS in the frame-work of
$LP$-Kenmotsu manifolds of dimension $n$. We recommend the papers [1, 2, 5, 9,
10, 13, 15, 16, 17, 18, 19] and the references therein for more details about
the related studies.

##  2\. Preliminaries

An $n$-dimensional differentiable manifold $M$ with structure
$(\varphi,\zeta,\nu,g)$ is said to be a Lorentzian almost paracontact metric
manifold, if it admits a $(1,1)$-tensor field $\varphi$, a contravariant
vector field $\zeta$, a 1-form $\nu$ and a Lorentzian metric $g$ satisfying

(2.1) ${\nu{(\zeta)}+1=0},$ (2.2) ${\varphi}^{2}{E}={E}+{\nu}(E){\zeta},$
(2.3) $\varphi\zeta=0,\quad\nu(\varphi E)=0,$ (2.4)
$g(\varphi{E},\varphi{F})=g(E,F)+{\nu}(E){\nu}(F),$ (2.5)
$g(E,{\zeta})={\nu}(E),$ (2.6) $\varphi(E,F)=\varphi(F,E)=g(E,\varphi F)$

for any vector fields $E,F\in\chi(M)$, where $\chi(M)$ is the Lie algebra of
vector fields on $M$.
If $\zeta$ is a killing vector field, the (para) contact structure is called a
$K$-(para) contact. In such a case, we have

(2.7) $\nabla_{E}\zeta=\varphi E.$

Recently, the authors Haseeb and Prasad defined and studied the following
notion:

######  Definition 2.1.

A Lorentzian almost paracontact manifold $M$ is called Lorentzian para-
Kenmostu manifold if [11]

(2.8) $(\nabla_{E}\varphi)F=-g(\varphi E,F)\zeta-\nu(F)\varphi E$

for any $E,F$ on $M.$

In an $LP$-Kenmostu manifold, we have

(2.9) $\nabla_{E}\zeta=-E-\nu(E)\zeta,$ (2.10)
$(\nabla_{E}\nu)F=-g(E,F)-\nu(E)\nu(F),$

where $\nabla$ denotes the Levi-Civita connection respecting to the Lorentzian
metric $g$.

Furthermore, in an $LP$-Kenmotsu manifold, the following relations hold [11]:

(2.11) $g(R(E,F)G,\zeta)=\nu(R(E,F)G)=g(F,G)\nu(E)-g(E,G)\nu(F),$ (2.12)
$R(\zeta,E)F=-R(E,\zeta)F=g(E,F)\zeta-\nu(F)E,$ (2.13)
$R(E,F)\zeta=\nu(F)E-\nu(E)F,$ (2.14) $R(\zeta,E)\zeta=E+\nu(E)\zeta,$ (2.15)
$S(E,\zeta)=(n-1)\nu(E),~{}~{}S(\zeta,\zeta)=-(n-1),$ (2.16)
$Q\zeta=(n-1)\zeta$

for any $E,F,G\in\chi(M)$, where $R,S$ and $Q$ represent the curvature tensor,
the Ricci tensor and the $Q$ Ricci operator, respectively.

######  Definition 2.2.

[21] An $LP$-Kenmotsu manifold $M$ is said to be $\nu$-Einstein manifold if
its $S(\neq 0)$ is of the form

(2.17) $S(E,F)=ag(E,F)+b\nu(E)\nu(F),$

where $a$ and $b$ are smooth functions on $M$. In particular, if $b=0$, then
$M$ is termed as an Einstein manifold.

######  Remark 2.3.

[12] In an $LP$-Kenmotsu manifold of $n$-dimension, $S$ is of the form

(2.18) $S(E,F)=(\frac{r}{n-1}-1)g(E,F)+(\frac{r}{n-1}-n)\nu(E)\nu(F),$

where $r$ is the scalar curvature of the manifold.

######  Lemma 2.4.

In an $n$-dimensional $LP$-Kenmotsu manifold, we have

(2.19) $\zeta(r)=2(r-n(n-1)),$ (2.20) $(\nabla_{E}Q)\zeta=QE-(n-1)E,$ (2.21)
$(\nabla_{\zeta}Q)E=2QE-2(n-1)E$

for any $E$ on $M$.

###### Proof.

Equation (2.18) yields

(2.22) $QE=(\frac{r}{n-1}-1)E+(\frac{r}{n-1}-n)\nu(E)\zeta.$

Taking the covariant derivative of (2.22) with respect to $F$ and making use
of (2.9) and (2.10), we lead to

$\displaystyle(\nabla_{F}Q)E=\frac{F(r)}{n-1}(E+\nu(E)\zeta)-(\frac{r}{n-1}-n)(g(E,F)\zeta+\nu(E)F+2\nu(E)\nu(F)\zeta).$

By contracting $F$ in the foregoing equation and using trace
$\\{F\to(\nabla_{F}Q)E\\}=\frac{1}{2}E(r)$, we find

$\displaystyle\frac{n-3}{2(n-1)}E(r)=\big{\\{}\frac{\zeta(r)}{n-1}-(r-n(n-1))\big{\\}}\nu(E),$

which by replacing $E$ by $\zeta$ and using (2.1) gives (2.19). We refer the
readers to see [14] for the proof of (2.20) and (2.21). ∎

######  Remark 2.5.

From the equation $(\ref{2.20})$, it is noticed that if an $n$-dimensional
$LP$-Kenmotsu manifold possesses the constant scalar curvature, then
$r=n(n-1)$ and hence (2.18) reduces to $S(E,F)=(n-1)g(E,F)$. Thus, the
manifold under consideration is an Einstein manifold.

##  3\. Ricci-Yamabe solitons on $LP$-Kenmotsu manifolds

Let the metric of an $n$-dimensional $LP$-Kenmotsu manifold be a Ricci-Yamabe
soliton $(g,K,\Lambda,\sigma,\rho)$, then $(\ref{1.1})$ holds. By
differentiating $(\ref{1.1})$ covariantly with resprct to $G$, we have

(3.1) $\displaystyle(\nabla_{G}\pounds_{K}g)(E,F)$ $\displaystyle=$
$\displaystyle-2\sigma(\nabla_{G}S)(E,F)+\rho(Gr)g(E,F).$

Since $\nabla g=0$, then the following formula [20]

$(\pounds_{K}\nabla_{E}g-\nabla_{E}\pounds_{K}g-\nabla_{[K,E]}g)(F,G)=-g((\pounds_{K}\nabla)(E,F),G)-g((\pounds_{K}\nabla)(E,G),F)$

turns to

$(\nabla_{E}\pounds_{K}g)(F,G)=g((\pounds_{K}\nabla)(E,F),G)+g((\pounds_{K}\nabla)(E,G),F).$

Since the operator $\pounds_{K}\nabla$ is symmetric, therefore we have

$2g((\pounds_{K}\nabla)(E,F),G)=(\nabla_{E}\pounds_{K}g)(F,G)+(\nabla_{F}\pounds_{K}g)(E,G)-(\nabla_{G}\pounds_{K}g)(E,F),$

which by using (3.1) takes the form

(3.2) $\displaystyle 2g((\pounds_{K}\nabla)(E,F),G)$ $\displaystyle=$
$\displaystyle-2\sigma[(\nabla_{E}S)(F,G)+(\nabla_{F}S)(G,E)+(\nabla_{G}S)(E,F)]$
$\displaystyle+\rho[(Er)g(F,G)+(Fr)g(G,E)+(Gr)g(E,F)].$

Putting $F=\zeta$ in (3.2) and using (2.5), we find

(3.3) $\displaystyle 2g((\pounds_{K}\nabla)(E,\zeta),G)$ $\displaystyle=$
$\displaystyle-2\sigma[(\nabla_{E}S)(\zeta,G)+(\nabla_{\zeta}S)(G,E)-(\nabla_{G}S)(E,\zeta)]$
$\displaystyle+\rho[(Er)\nu(G)+2(r-n(n-1))g(E,G)-(Gr)\nu(E)]$

By virtue of (2.20) and (2.21), (3.3) leads to

$\displaystyle 2g((\pounds_{K}\nabla)(E,\zeta),G)$ $\displaystyle=$
$\displaystyle-4\sigma[S(E,G)-(n-1)g(E,G)]$
$\displaystyle+\rho[(Er)\nu(G)+2(r-n(n-1))g(E,G)-(Gr)\nu(E)].$

By eliminating $G$ from the foregoing equation, we have

$\displaystyle 2(\pounds_{K}\nabla)(F,\zeta)$ $\displaystyle=$
$\displaystyle\rho g(Dr,F)\zeta-\rho(Dr)\nu(F)-4\sigma QF$
$\displaystyle+[4\sigma(n-1)+2\rho(r-n(n-1))]F.$

If we take $r$ as constant, then from (2.19) we find $r=n(n-1)$, and hence (3)
reduces to

(3.5) $\displaystyle(\pounds_{K}\nabla)(F,\zeta)$ $\displaystyle=$
$\displaystyle-2\sigma QF+2\sigma(n-1)F.$

Taking covariant derivative of (3.5) with respect to $E$, we have

$\displaystyle(\nabla_{E}\pounds_{K}\nabla)(F,\zeta)$ $\displaystyle=$
$\displaystyle(\pounds_{K}\nabla)(F,E)-2\sigma\nu(E)[QF-(n-1)F]$
$\displaystyle-$ $\displaystyle 2\sigma(\nabla_{E}Q)F.$

Again from [20], we have

$\displaystyle(\pounds_{K}R)(E,F)G=(\nabla_{E}\pounds_{K}\nabla)(F,G)-(\nabla_{F}\pounds_{K}\nabla)(E,G),$

which by putting $G=\zeta$ and using (3) takes the form

$\displaystyle(\pounds_{K}R)(E,F)\zeta$ $\displaystyle=$ $\displaystyle
2\sigma\nu(F)(QE-(n-1)E)-2\sigma\nu(E)(QF-(n-1)F)$
$\displaystyle-2\sigma((\nabla_{E}Q)F-(\nabla_{F}Q)E).$

Putting $F=\zeta$ in (3) then using (2.1), (2.2), (2.20) and (2.21), we arrive
at

(3.8) $\displaystyle(\pounds_{K}R)(E,\zeta)\zeta=0.$

The Lie derivative of $R(E,\zeta)\zeta=-E-\nu(E)\zeta$ along $K$ leads to

(3.9)
$\displaystyle(\pounds_{K}R)(E,\zeta)\zeta-g(E,\pounds_{K}\zeta)\zeta+2\nu(\pounds_{K}\zeta)E=-(\pounds_{K}\nu)(E)\zeta.$

From (3.8) and (3.9), we have

(3.10)
$\displaystyle(\pounds_{K}\nu)(E)\zeta=-2\nu(\pounds_{K}\zeta)E+g(E,\pounds_{K}\zeta)\zeta.$

Taking the Lie derivative of $g(E,\zeta)=\nu(E)$, we find

(3.11)
$\displaystyle(\pounds_{K}\nu)(E)=g(E,\pounds_{K}\zeta)+(\pounds_{K}g)(E,\zeta).$

By putting $F=\zeta$ in (1.1) and using (2.15), we have

(3.12) $\displaystyle(\pounds_{K}g)(E,\zeta)=-\\{2\sigma(n-1)+2\Lambda-\rho
n(n-1)\\}\nu(E),$

where $r=n(n-1)$ being used.
The Lie derivative of $g(\zeta,\zeta)=-1$ along $K$ we lead to

(3.13) $\displaystyle(\pounds_{K}g)(\zeta,\zeta)=-2\nu(\pounds_{K}\zeta).$

From (3.12) and (3.16), we find

(3.14) $\displaystyle\nu(\pounds_{K}\zeta)=-\\{\sigma(n-1)+\Lambda-\frac{\rho
n(n-1)}{2}\\}.$

Now, combining the equations (3.10), (3.11), (3.12) and (3.17), we find

(3.15) $\displaystyle\Lambda=\frac{\rho n(n-1)}{2}-\sigma(n-1).$

Thus, we have

######  Theorem 3.1.

Let $(M,g)$ be an $n$-dimensional $LP$-Kenmotsu manifold admitting Ricci-
Yamabe soliton $(g,K,\Lambda,\sigma,\rho)$ with constant scalar curvature
tensor, then $\Lambda=\frac{\rho n(n-1)}{2}-\sigma(n-1).$

For $\sigma=1$ and $\rho=0$, from (3.15) we have $\Lambda=-(n-1)$. Thus, we
have the following:

######  Corollary 3.2.

If an $n$-dimensional $LP$-Kenmotsu manifold admits a Ricci soliton with
constant scalar curvature, then the soliton is shrinking.

For $\sigma=0$ and $\rho=1$, from (3.15) we have $\Lambda=\frac{n(n-1)}{2}$.
Thus, we have the following:

######  Corollary 3.3.

If an $n$-dimensional $LP$-Kenmotsu manifold admits a Yamabe soliton with
constant scalar curvature, then the soliton is shrinking.

For $\sigma=1$ and $\rho=-1$, from (3.15) we have
$\Lambda=-\frac{(n^{2}-1)}{2}$. Thus, we have the following:

######  Corollary 3.4.

If an $n$-dimensional $LP$-Kenmotsu manifold admits an Einstein soliton with
constant scalar curvature, then the soliton is shrinking.

Now, we consider the metric of an $n$-dimensional $LP$-Kenmotsu manifold as a
Ricci-Yamabe soliton $(g,\zeta,\Lambda,\sigma,\rho)$, then from $(\ref{1.1})$
and (2.9) we have

(3.16) $\displaystyle\ \ \ \ S(E,F)=-\frac{1}{\sigma}(\Lambda-1-\frac{\rho
r}{2})g(E,F)+\frac{1}{\sigma}\nu(E)\nu(F),\ \ where\ \ \sigma\neq 0.$

By putting $F=\zeta$ in (3.16) and using (2.15), we find

(3.17) $\displaystyle\Lambda=\frac{\rho r}{2}-\sigma(n-1).$

Now, comparing (2.18) and (3.17), we have $r=\frac{n-1}{\sigma}+n(n-1)$, which
by using in (3.17) it follows that
$\Lambda=-\sigma(n-1)+\frac{\rho(n-1)(1+n\sigma)}{2\sigma}.$ Thus, we have the
following theorem:

######  Theorem 3.5.

An $n$-dimensional $LP$-Kenmotsu manifold with constant scalar curvature
admitting Ricci-Yamabe soliton $(g,\zeta,\Lambda,\sigma,\rho)$ is an
$\nu$-Einstein manifold. Moreover, the soliton is expanding, steady or
shrinking according to $\frac{\rho}{\sigma}>2\sigma-\rho n$,
$\frac{\rho}{\sigma}=2\sigma-\rho n$, or $\frac{\rho}{\sigma}<2\sigma-\rho n$.

##  4\. Gradient Ricci-Yamabe solitons on $LP$-Kenmotsu manifolds

######  Definition 4.1.

A Riemannian $($or semi-Riemannian$)$ metric $g$ on $M$ is called a gradient
RYS, if

(4.1) $\displaystyle Hessv+\sigma S+(\Lambda-\frac{\rho r}{2})g=0,$

where $Hessv$ denotes the Hessian of a smooth function $v$ on $M$ and defined
by $Hessv=\nabla\nabla v$.

Let $M$ be an $n$-dimensional $LP$-Kenmotsu manifold with $g$ as a gradient
RYS. Then equation (4.1) can be written as

(4.2) $\displaystyle\nabla_{E}Dv+\sigma QE+(\Lambda-\frac{\rho r}{2})E=0,$

for all vector fields $E$ on $M$, where $D$ denotes the gradient operator of
$g$. Taking the covariant derivative of (4.2) with respect to $F$, we have

(4.3) $\displaystyle\ \ \
\nabla_{F}\nabla_{E}Dv=-\sigma\\{(\nabla_{F}Q)E+Q(\nabla_{F}E)\\}+\rho\frac{F(r)}{2}E-(\Lambda-\frac{\rho
r}{2})\nabla_{F}E.$

Interchanging $E$ and $F$ in (4.3), we lead to

(4.4) $\displaystyle\ \ \
\nabla_{E}\nabla_{F}Dv=-\sigma\\{(\nabla_{E}Q)F+Q(\nabla_{E}F)\\}+\rho\frac{E(r)}{2}F-(\Lambda-\frac{\rho
r}{2})\nabla_{E}F.$

By making use of (4.2)-(4.4), we find

(4.5) $\displaystyle
R(E,F)Dv=\sigma\\{(\nabla_{F}Q)E-(\nabla_{E}Q)F\\}+\frac{\rho}{2}\\{E(r)F-F(r)E\\}.$

Now, from (2.18), we find

$\displaystyle QE=(\frac{r}{n-1}-1)E+(\frac{r}{n-1}-n)\nu(E)\zeta,$

which on taking covariant derivative with repect to $F$ leads to

$\displaystyle(\nabla_{F}Q)E$ $\displaystyle=$
$\displaystyle\frac{F(r)}{n-1}(E+\nu(E)\zeta)-(\frac{r}{n-1}-n)(g(E,F)\zeta$
$\displaystyle+2\nu(E)\nu(F)\zeta+\nu(E)F).$

By using (4) in (4.5), we have

(4.7) $\displaystyle R(E,F)Dv$ $\displaystyle=$
$\displaystyle\frac{(n-1)\rho-2\sigma}{2(n-1)}\\{E(r)F-F(r)E\\}+\frac{\sigma}{n-1}\\{F(r)\nu(E)\zeta-E(r)\nu(F)\zeta\\}$
$\displaystyle-\sigma(\frac{r}{n-1}-n)(\nu(E)F-\nu(F)E).$

Contracting forgoing equation along $E$ gives

$\displaystyle S(F,Dv)$ $\displaystyle=$
$\displaystyle\big{\\{}\frac{(n-1)^{2}\rho-2\sigma(n-2)}{n-1}\big{\\}}F(r)$
$\displaystyle+\frac{\sigma(n-3)(r-n(n-1))}{n-1}\nu(F).$

From the equation (2.18), we can write

(4.9) $\displaystyle
S(F,Dv)=(\frac{r}{n-1}-1)F(v)+(\frac{r}{n-1}-n)\nu(F)\zeta(v).$

Now, by equating (4) and (4.9), then putting $F=\zeta$ and using (2.1),
(2.19), we find

(4.10)
$\displaystyle\zeta(v)=\frac{r-n(n-1)}{n-1}\\{2(n-1)\rho-\frac{\sigma(5n-13)}{n-1}\\}.$

Taking the inner product of (4.7) with $\zeta$, we get

$\displaystyle
F(v)\nu(E)-E(v)\nu(F)=\frac{\rho}{2}\\{E(r)\nu(F)-F(r)\nu(E)\\},$

which by replacing $E$ by $\zeta$ and using (2.19), (4.10), we infer

(4.11) $\displaystyle
F(v)=-(r-n(n-1))\\{3\rho-\frac{\sigma(5n-13)}{(n-1)^{2}}\\}\nu(F)-\frac{\rho}{2}F(r).$

If we take $r$ as constant, then from Remark 2.5, we get $r=n(n-1)$. Thus,
(4.11) leads to $F(v)=0$. This implies that $v$ is constant. Thus, the soliton
under the consideration is trivial. Hence we state:

######  Theorem 4.2.

If the metric of an $LP$-Kenmotsu manifold of constant scalar curvature tensor
admitting a special type of vector field is gradient RYS, then the soliton is
trivial.

For $v$ constant, (1.2) turns to

$\displaystyle\sigma QE=-(\Lambda-\frac{\rho r}{2})E,$

which leads to

(4.12) $\displaystyle S(E,F)=-\frac{1}{\sigma}(\Lambda-\frac{\rho
n(n-1)}{2})g(E,F),\ \ \ \sigma\neq 0.$

By putting $E=F=\zeta$ in (4.12) and using (2.15), we obtain

(4.13) $\displaystyle\Lambda=\frac{\rho n(n-1)}{2}-\sigma(n-1).$

######  Corollary 4.3.

If an $n$-dimensional $LP$-Kenmotsu manifold admits a gradient Ricci soliton
with the constant scalar curvature, then the manifold under the consideration
is an Einstein manifold and $\Lambda=\frac{\rho n(n-1)}{2}-\sigma(n-1).$

For $\sigma=1$ and $\rho=0$, from (4.13) we find $\Lambda=-(n-1)$. Thu, we
have the following:

######  Corollary 4.4.

If an $n$-dimensional $LP$-Kenmotsu manifold admits a gradient Ricci soliton
with the constant scalar curvature, then the soliton is shrinking.

For $\sigma=1$ and $\rho=-1$, from (4.13) we have
$\Lambda=-\frac{(n-1)(n+2)}{2}$. Thus, we have the following:

######  Corollary 4.5.

If an $n$-dimensional $LP$-Kenmotsu manifold admits an gradient Einstein
soliton with constant scalar curvature, then the soliton is shrinking.

Example. We consider the 5-dimensional manifold
$M^{5}=\left\\{(x_{1},x_{2},x_{3},x_{4},x_{5})\in\mathbb{R}^{5}:x_{5}>0\right\\}$,
where $(x_{1},x_{2},x_{3},x_{4},x_{5})$ are the standard coordinates in
$\mathbb{R}^{5}$. Let $\varrho_{1},$ $\varrho_{2}$, $\varrho_{3}$,
$\varrho_{4}$ and $\varrho_{5}$ be the vector fields on $M^{5}$ given by

$\varrho_{1}=e^{x_{5}}\frac{\partial}{\partial
x_{1}},~{}~{}\varrho_{2}=e^{x_{5}}\frac{\partial}{\partial{x_{2}}},~{}~{}\varrho_{3}=e^{x_{5}}\frac{\partial}{\partial{x_{3}}},~{}~{}\varrho_{4}=e^{x_{5}}\frac{\partial}{\partial
x_{4}},~{}~{}\varrho_{5}=\frac{\partial}{\partial x_{5}}=\zeta,$

which are linearly independent at each point of $M^{5}$. Let $g$ be the
Lorentzian metric defined by

$\ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ g(\varrho_{i},\varrho_{i})=1,\
\quad{\text{for}}\ \quad 1\leq i\leq 4\ \quad{\text{and}}\ \quad
g(\varrho_{5},\varrho_{5})=-1,$ $g(\varrho_{i},\varrho_{j})=0,\
\quad{\text{for}}\ \quad i\neq j,\ \quad 1\leq i,j\leq 5.$

Let $\nu$ be the 1-form defined by $\nu(E)=g(E,\varrho_{5})=g(\varrho,\zeta)$
for all $E\in\chi{(M^{5})}$, and let $\varphi$ be the $(1,1)$-tensor field
defined by

$\varphi\varrho_{1}=-\varrho_{2},~{}~{}\varphi\varrho_{2}=-\varrho_{1},~{}~{}\varphi\varrho_{3}=-\varrho_{4},~{}~{}\varphi\varrho_{4}=-\varrho_{3},~{}~{}\varphi\varrho_{5}=0.$

By applying linearity of $\varphi$ and $g$, we have

$\nu(\zeta)=g(\zeta,\zeta)=-1,~{}~{}\varphi^{2}E=E+\nu(E)\zeta~{}~{}{\text{and}}~{}~{}g(\varphi
E,\varphi F)=g(E,F)+\nu(E)\nu(F)$

for all $E,F\in\chi{(M^{5})}$. Thus for $\varrho_{5}=\zeta$, the structure
$(\varphi,\zeta,\nu,g)$ defines a Lorentzian almost paracontact metric
structure on $M^{5}$. Then we have

$[\varrho_{i},\varrho_{j}]=-\varrho_{i},\ \quad{\text{for}}\ \quad 1\leq i\leq
4,j=5,\ \ \ $ $[\varrho_{i},\varrho_{j}]=0,\ \quad{\text{otherwise}}.\ \ \ \ \
\ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ $

By using Koszul’s formula, we can easily find we obtain

$\displaystyle\nabla_{\varrho_{i}}{\varrho}_{j}=\begin{cases}-\varrho_{5},\
\quad 1\leq i=j\leq 4,\\\ -\varrho_{i},\ \quad 1\leq i\leq 4,j=5,\\\ 0,\ \quad
otherwise.\end{cases}$

Also one can easily verify that

$\nabla_{E}\zeta=-E-\eta(E)\zeta\ \quad{\text{and}}\
\quad(\nabla_{E}\varphi)F=-g(\varphi E,F)\zeta-\nu(F)\varphi E.$

Therefore, the manifold is an $LP$-Kenmotsu manifold.
From the above results, we can easily obtain the non-vanishing components of
$R$ as follows:

$R(\varrho_{1},\varrho_{2})\varrho_{1}=-\varrho_{2},~{}~{}~{}R(\varrho_{1},\varrho_{2})\varrho_{2}=\varrho_{1},~{}~{}~{}R(\varrho_{1},\varrho_{3})\varrho_{1}=-\varrho_{3},~{}~{}~{}R(\varrho_{1},\varrho_{3})\varrho_{3}=\varrho_{1},$
$R(\varrho_{1},\varrho_{4})\varrho_{1}=-v_{4},~{}~{}~{}R(\varrho_{1},\varrho_{4})\varrho_{4}=\varrho_{1},~{}~{}~{}R(\varrho_{1},\varrho_{5})\varrho_{1}=-\varrho_{5},~{}~{}~{}R(\varrho_{1},\varrho_{5})\varrho_{5}=-\varrho_{1},$
$R(\varrho_{2},\varrho_{3})\varrho_{2}=-\varrho_{3},~{}~{}~{}R(\varrho_{2},\varrho_{3})\varrho_{3}=\varrho_{2},~{}~{}~{}R(\varrho_{2},\varrho_{4})\varrho_{2}=-\varrho_{4},~{}~{}~{}R(\varrho_{2},\varrho_{4})\varrho_{4}=\varrho_{2},$
$R(\varrho_{2},\varrho_{5})\varrho_{2}=-\varrho_{5},~{}~{}~{}R(\varrho_{2},\varrho_{5})\varrho_{5}=-\varrho_{2},~{}~{}~{}R(\varrho_{3},\varrho_{4})\varrho_{3}=-\varrho_{4},~{}~{}~{}R(\varrho_{3},\varrho_{4})\varrho_{4}=\varrho_{3},$
$R(\varrho_{3},\varrho_{5})\varrho_{3}=-\varrho_{5},~{}~{}~{}R(\varrho_{3},\varrho_{5})\varrho_{5}=-\varrho_{3},~{}~{}~{}R(\varrho_{4},\varrho_{5})\varrho_{4}=-\varrho_{5},~{}~{}~{}R(\varrho_{4},\varrho_{5})\varrho_{5}=-\varrho_{4}.$

Also, we calculate the Ricci tensors as follows:

$S(\varrho_{1},\varrho_{1})=S(\varrho_{2},\varrho_{2})=S(\varrho_{3},\varrho_{3})=S(\varrho_{4},\varrho_{4})=4,\
\ \ \ S(\varrho_{5},\varrho_{5})=-4.$

Therefore, we have

$r=S(\varrho_{1},\varrho_{1})+S(\varrho_{2},\varrho_{2})+S(\varrho_{3},\varrho_{3})+S(\varrho_{4},\varrho_{4})-S(\varrho_{5},\varrho_{5})=20.$

Now by taking
$Dv=(\varrho_{1}v)\varrho_{1}+(\varrho_{2}v)\varrho_{2}+(\varrho_{3}v)\varrho_{3}+(\varrho_{4}v)\varrho_{4}+(\varrho_{5}v)\varrho_{5}$,
we have

$\displaystyle\nabla_{\varrho_{1}}Dv=(\varrho_{1}(\varrho_{1}v)-(\varrho_{5}v)){\varrho_{1}}+{(\varrho_{1}}(\varrho_{2}v)){\varrho_{2}}+{(\varrho_{1}}(\varrho_{3}v)){\varrho_{3}}+{(\varrho_{1}}(\varrho_{4}v)){\varrho_{4}}+(\varrho_{1}(\varrho_{5}v)-(\varrho_{1}v)){\varrho_{5}},$
$\displaystyle\nabla_{\varrho_{2}}Dv=(\varrho_{2}(\varrho_{1}v))\varrho_{1}+(\varrho_{2}(\varrho_{2}v)-(\varrho_{5}v)){\varrho_{2}}+{(\varrho_{2}}(\varrho_{3}v)){\varrho_{3}}+{(\varrho_{2}}(\varrho_{4}v)){\varrho_{4}}+(\varrho_{2}(\varrho_{5}v)-(\varrho_{2}v)){\varrho_{5}},$
$\displaystyle\nabla_{\varrho_{3}}Dv=(\varrho_{3}(\varrho_{1}v))\varrho_{1}+(\varrho_{3}(\varrho_{2}v))\varrho_{2}+(\varrho_{3}(\varrho_{3}v)-(\varrho_{5}v)){\varrho_{3}}+(\varrho_{3}(\varrho_{4}v)){\varrho_{4}}+(\varrho_{3}(\varrho_{5}v)-(\varrho_{3}v)){\varrho_{5}},$
$\displaystyle\nabla_{\varrho_{4}}Dv=(\varrho_{4}(\varrho_{1}v))\varrho_{1}+(\varrho_{4}(\varrho_{2}v))\varrho_{2}+(\varrho_{4}(\varrho_{3}v))\varrho_{3}+(\varrho_{4}(\varrho_{4}v)-(\varrho_{5}v)){\varrho_{4}}+(\varrho_{4}(\varrho_{5}v)-(\varrho_{4}v)){\varrho_{5}},$
$\displaystyle\nabla_{\varrho_{5}}Dv=(\varrho_{5}(\varrho_{1}v))\varrho_{1}+(\varrho_{5}(\varrho_{2}v))\varrho_{2}+(\varrho_{5}(\varrho_{3}v))\varrho_{3}+(\varrho_{5}(\varrho_{4}v))\varrho_{4}+(\varrho_{5}(\varrho_{5}v)){\varrho_{5}}.$

Thus, by virtue of (4.2), we obtain

(4.14)
$\displaystyle\begin{cases}\varrho_{1}(\varrho_{1}v)-\varrho_{5}v=-(\Lambda+4\sigma-10\rho),\\\
\varrho_{2}(\varrho_{2}v)-\varrho_{5}v=-(\Lambda+4\sigma-10\rho),\\\
\varrho_{3}(\varrho_{3}v)-\varrho_{5}v=-(\Lambda+4\sigma-10\rho),\\\
\varrho_{4}(\varrho_{4}v)-\varrho_{5}v=-(\Lambda+4\sigma-10\rho),\\\
\varrho_{5}(\varrho_{5}v)=-(\Lambda+4\sigma-10\rho),\\\
\varrho_{1}(\varrho_{2}v)=\varrho_{1}(\varrho_{3}v)=\varrho_{1}(\varrho_{4}v)=0,\\\
\varrho_{2}(\varrho_{1}v)=\varrho_{2}(\varrho_{3}v)=\varrho_{2}(\varrho_{4}v)=0,\\\
\varrho_{3}(\varrho_{1}v)=\varrho_{3}(\varrho_{2}v)=\varrho_{3}(\varrho_{4}v)=0,\\\
\varrho_{4}(\varrho_{1}v)=\varrho_{4}(\varrho_{2}v)=\varrho_{4}(\varrho_{3}v)=0,\\\
\varrho_{1}(\varrho_{5}v)-(\varrho_{1}v)=\varrho_{2}(\varrho_{5}v)-(\varrho_{2}v)=0,\\\
\varrho_{3}(\varrho_{5}v)-(\varrho_{3}v)=\varrho_{4}(\varrho_{5}v)-(\varrho_{4}v)=0.\end{cases}$

Thus, the equations in (4.14) are respectively amounting to

$\displaystyle e^{2x_{5}}\frac{\partial^{2}v}{\partial
x_{1}^{2}}-\frac{\partial v}{\partial x_{5}}=-(\Lambda+4\sigma-10\rho),$
$\displaystyle e^{2x_{5}}\frac{\partial^{2}v}{\partial
x_{2}^{2}}-\frac{\partial v}{\partial x_{5}}=-(\Lambda+4\sigma-10\rho),$
$\displaystyle e^{2x_{5}}\frac{\partial^{2}v}{\partial
x_{3}^{2}}-\frac{\partial v}{\partial x_{5}}=-(\Lambda+4\sigma-10\rho),$
$\displaystyle e^{2x_{5}}\frac{\partial^{2}v}{\partial
x_{4}^{2}}-\frac{\partial v}{\partial x_{5}}=-(\Lambda+4\sigma-10\rho),$
$\displaystyle\frac{\partial^{2}v}{\partial
x_{5}^{2}}=-(\Lambda+4\sigma-10\rho),$
$\displaystyle\frac{\partial^{2}v}{\partial x_{1}\partial
x_{2}}=\frac{\partial^{2}v}{\partial x_{1}\partial
x_{3}}=\frac{\partial^{2}v}{\partial x_{1}\partial
x_{4}}=\frac{\partial^{2}v}{\partial x_{2}\partial
x_{3}}=\frac{\partial^{2}v}{\partial x_{2}\partial
x_{4}}=\frac{\partial^{2}v}{\partial x_{3}\partial x_{4}}=0,$ $\displaystyle
e^{x_{5}}\frac{\partial^{2}v}{\partial x_{5}\partial x_{1}}+\frac{\partial
v}{\partial x_{1}}=e^{x_{5}}\frac{\partial^{2}v}{\partial x_{5}\partial
x_{2}}+\frac{\partial v}{\partial
x_{2}}=e^{x_{5}}\frac{\partial^{2}v}{\partial x_{5}\partial
x_{3}}+\frac{\partial v}{\partial
x_{3}}=e^{x_{5}}\frac{\partial^{2}v}{\partial x_{5}\partial
x_{4}}+\frac{\partial v}{\partial x_{4}}=0.$

From the above equations it is observed that $v$ is constant for
$\Lambda=-4\sigma+10\rho$. Hence, equation (4.2) is satisfied. Thus, $g$ is a
gradient RYS with the soliton vector field $K=Dv$, where $v$ is constant and
$\Lambda=-4\sigma+10\rho$. Hence, Theorem 4.2 is verified.

## References

  * [1] Blaga, A. M., Solitons and geometrical structure in a perfect fluid spacetime, Rocky Mt. J. Math. (2020).
  * [2] Blaga, A. M., Some geometrical aspects of Einstein, Ricci and Yamabe solitons, J. Geom. Symmetry Phys., 52 (2019), 17-26.
  * [3] Catino, G. and Mazzieri, L., Gradient Einstein solitons, Nonlinear Anal., 132(2016), 66-94.
  * [4] Catino, G., Cremaschi, L., Djadli, Z., Mantegazza, C. and Mazzieri, L., The Ricci Bourguignon flow, Pacific J. Math., 28(2017), 337-370.
  * [5] Chidananda, S., and Venkatesha, V., Yamabe soliton and Riemann soliton on Lorentzian para-Sasakian manifold, Commun. Korean Math. Soc., https://doi.org/10.4134/CKMS.c200365.
  * [6] Güler, S. and Crasmareanu, M., Ricci-Yamabe maps for Riemannian flows and their volume variation and volume entropy, Turk. J. Math., 43 (2019), 2631-2641.
  * [7] Hamilton, R. S., Lectures on Geometric Flows (Unpublished manuscript, 1989).
  * [8] Hamilton, R. S., The Ricci Flow on Surfaces, Mathematics and General Relativity (Santa Cruz, CA, 1986), Contemp. Math., A.M.S., 71 (1988), 237-262.
  * [9] Haseeb, A. and De, U. C., $\eta$-Ricci solitons in $\epsilon$-Kenmotsu manifolds, J. Geom. 110, 34 (2019).
  * [10] Haseeb, A. and Almusawa, H., Some results on Lorentzian para-Kenmotsu manifolds admitting $\eta$-Ricci solitons, Palestine Journal of Mathematics, 11(2)(2022), 205-213.
  * [11] Haseeb, A. and Prasad, R., Certain results on Lorentzian para-Kenmotsu manifolds, Bol. Soc. Parana. Mat., 39(3) (2021), 201-220.
  * [12] Haseeb, A. and Prasad, R., Some results on Lorentzian para-Kenmotsu manifolds, Bull. Transilvania Univ. of Brasov, 13(62) (2020), no. 1, 185-198.
  * [13] Haseeb, A., Prasad, R. and Mofarreh, F., Sasakian manifolds admitting $*$-$\eta$-Ricci-Yamabe solitons, Advances in Mathematical Physics, Vol. 2022, Article ID 5718736, 7 pages. doi: https:// doi.org/10.1155/2022/5718736
  * [14] Li, Y., Haseeb, A. and Ali, M., $LP$-Kenmotsu manifolds admitting $\eta$-Ricci solitons and spacetime, Journal of Mathematics, 2022, Article ID 6605127, 10 pages.
  * [15] Lone, M. A. and Harry, I. F., Ricci Solitons on Lorentz-Sasakian space forms, Journal of Geometry and Physics, 104547, doi: https://doi.org/10.1016/j.geomphys.2022.104547.
  * [16] Pankaj, Chaubey, S. K and Prasad, R., Three dimensional Lorentzian para-Kenmotsu manifolds and Yamabe soliton, Honam Mathematical J., 43(4) (2021), 613-626.
  * [17] Singh, J. P. and Khatri, M., On Ricci-Yamabe soliton and geometrical structure in a perfect fluid spacetime, Afr. Mat., 32(2021), 1645-1656.
  * [18] Venkatesha, Kumara, H. A., Ricci soliton and geometrical structure in a perfect fluid spacetime with torse-forming vector field, Afr. Mat. 30 (2019), 725-736
  * [19] Yoldas, H. I., On Kenmotsu manifolds admitting $\eta$-Ricci-Yamabe solitons,Int. J. Geom. Met. Mod. Phy., 18(12) (2021), 2150189.
  * [20] Yano, K., Integral Formulas in Riemannian geometry, Pure and Applied Mathematics, Vol. I, Marcel Dekker, New York, 1970.
  * [21] Yano, K. and Kon, M., Structures on manifolds, World Scientific, (1984).

Mobin Ahmad
Department of Mathematics,
Integral University, Kursi Road,
Lucknow-226026.
Email<EMAIL_ADDRESS>

Gazala
Department of Mathematics,
Integral University, Kursi Road,
Lucknow-226026.
Email<EMAIL_ADDRESS>

Mohd. Bilal
Department of Mathematical Sciences,
Umm Ul Qura University,
Makkah, Saudi Arabia.
Email<EMAIL_ADDRESS>


# Constraining primordial non-Gaussianity from the large scale structure two-
point and three-point correlation functions

Z. Brown1, R. Demina1, A. G. Adame2, S. Avila3, E. Chaussidon4, S. Yuan4, V.
Gonzalez-Perez2,5, J. García-Bellido2, J. Aguilar4, S. Ahlen6, R. Blum7, D.
Brooks8, T. Claybaugh4, S. Cole9, A. de la Macorra10, B. Dey11, P. Doel8, K.
Fanning12,13, J. E. Forero-Romero14,15, E. Gaztañaga16,17,18, S. Gontcho A
Gontcho4, K. Honscheid19,20,21, C. Howlett22, S. Juneau7, R. Kehoe23, T.
Kisner4, M. Landriau4, L. Le Guillou24, M. Manera3,25, R. Miquel3,26, E.
Mueller27, A. Muñoz-Gutièrrez10, A. D. Myers28, J. Nie29, G. Niz30,31, N.
Palanque-Delabrouille4,32, C. Poppett 4,37,38, M. Rezaie35, G. Rossi36, E.
Sanchez37, E. Schlafly38, D. Schlegel4, M. Schubnell39,40, J. H. Silber4, D.
Sprayberry7, G. Tarlé40, M. Vargas-Magaña10, B. A. Weaver7, Z. Zhou29, and H.
Zou29

The authors’ affiliations are listed in Appendix A E-mail:
<EMAIL_ADDRESS>
(Accepted XXX. Received YYY; in original form ZZZ)

###### Abstract

Surveys of cosmological large-scale structure (LSS) are sensitive to the
presence of local primordial non-Gaussianity (PNG), and may be used to
constrain models of inflation. Local PNG, characterized by $f_{\mathrm{NL}}$,
the amplitude of the quadratic correction to the potential of a Gaussian
random field, is traditionally measured from LSS two-point and three-point
clustering via the power spectrum and bi-spectrum. We propose a framework to
measure $f_{\mathrm{NL}}$ using the configuration space two-point correlation
function (2pcf) monopole and three-point correlation function (3pcf) monopole
of survey tracers. Our model estimates the effect of the scale-dependent bias
induced by the presence of PNG on the 2pcf and 3pcf from the clustering of
simulated dark matter halos. We describe how this effect may be scaled to an
arbitrary tracer of the cosmological matter density. The 2pcf and 3pcf of this
tracer are measured to constrain the value of $f_{\mathrm{NL}}$. Using
simulations of luminous red galaxies observed by the Dark Energy Spectroscopic
Instrument (DESI), we demonstrate the accuracy and constraining power of our
model, and forecast the ability to constrain $f_{\mathrm{NL}}$ to a precision
of $\sigma_{f_{\mathrm{NL}}}\approx 22$ with one year of DESI survey data.

###### keywords:

cosmology: inflation – cosmology: large-scale structure of Universe –
cosmology: early Universe

††pubyear: 2023††pagerange: Constraining primordial non-Gaussianity from the
large scale structure two-point and three-point correlation functions–C

##  1 Introduction

The fields driving (and/or present during) inflation leave imprints on
primordial density fluctuations. In the case of a single scalar field,
primordial deviations from average density are expected to be distributed
according to a Gaussian of mean zero (Maldacena, 2003; Bartolo et al., 2004),
while more complex multi-field models (Creminelli et al., 2011) are expected
to deviate from Gaussianity. Thus, measurements of primordial non-Gaussianity
(PNG) can distinguish between different inflationary models.

For local non-Gaussianity, which depends only on the local value of the
potential, the primordial potential is written as

$\Phi=\phi+f_{\mathrm{NL}}(\phi^{2}-\langle\phi\rangle^{2})\ ,$ (1)

where $\phi$ denotes the potential of a Gaussian random field, and
$f_{\mathrm{NL}}$ is the amplitude of the quadratic correction to this
potential (Gangui et al., 1993; Mueller et al., 2019).

A variety of methods have been employed to measure $f_{\mathrm{NL}}$ from both
the early and late-time clustering properties of cosmological matter tracers.
Initially, $f_{\mathrm{NL}}$ was measured from the angular clustering of the
cosmic microwave background (CMB; Komatsu & Spergel, 2001, 2002). At present,
the bi-spectrum of CMB anisotropies observed by the Planck satellite
(Collaboration, 2020) offers the tightest constraints, where the most probable
value and 68% confidence limits (CL) are $f_{\mathrm{NL}}=-0.9\pm 5.1$.
Studies of cosmological large-scale structure (LSS) are also promising tests
of PNG. Specifically, Cunha et al. (2010) and Hotchkiss (2011) have shown that
the number density of galaxy clusters is sensitive to $f_{\mathrm{NL}}$.
Beyond the number density, or one-point statistics, additional sensitivity to
PNG has been observed in two-point (McDonald, 2008; Dalal et al., 2008; Slosar
et al., 2008; Alvarez et al., 2014; Biagetti, 2019; Mueller et al., 2021) and
three-point statistics (D’Amico et al., 2022; Cabass et al., 2022), as
characterized by the power spectrum and bi-spectrum, which are the Fourier
transforms of the configuration space two-point and three-point correlation
functions (2pcf and 3pcf, respectively). It was demonstrated that a
simultaneous study of the power spectrum and bi-spectrum of galaxy
distributions increases the sensitivity to PNG (Gualdi et al., 2021).

Analyses of LSS clustering, e.g. involving the detection of the baryon
acoustic oscillations (BAO) signal, are typically carried out in both
configuration and Fourier space (Alam et al., 2017). Here, we propose an
analogous strategy for the measurement of local PNG. In the following study,
we develop a method to measure $f_{\mathrm{NL}}$ by simultaneously using the
configuration-space 2pcf and 3pcf.

The effect of PNG on the clustering of cosmological matter tracers is a scale-
dependent bias, where the largest deviations in clustering from the Gaussian
case occur at large scales. Hence, greater sensitivity is expected from
current and future large-volume spectroscopic surveys, e.g. the Dark Energy
Spectroscopic Instrument (DESI). DESI is a spectroscopic surveyor attached to
the Mayall 4-meter telescope at Kitt Peak National Observatory (DESI
Collaboration et al., 2022; Miller et al., 2023; Silber et al., 2023; DESI
Collaboration et al., 2016b). Over an estimated five year period, DESI will
measure spectra of 40 million galaxies and quasars (DESI Collaboration et al.,
2016a). Although the primary goal of DESI is to investigate the nature of dark
energy through the Universe’s expansion history (Levi et al., 2013), its large
volume (about a third of the sky) also makes it ideal for measuring PNG from
LSS clustering. DESI will record spectra for a variety of cosmological matter
tracers at different redshift ranges (Allende Prieto et al., 2020; Ruiz-Macias
et al., 2020; Zhou et al., 2020; Raichoor et al., 2020; Yèche et al., 2020;
Lan et al., 2023; Alexander et al., 2023; Cooper et al., 2022; Hahn et al.,
2022; Zhou et al., 2023; Raichoor et al., 2023b; Chaussidon et al., 2023). In
this analysis, however, we study the sensitivity achievable with DESI using
luminous red galaxies (LRGs) as an example tracer.

We begin by briefly introducing the algorithm used to evaluate the correlation
functions, ConKer, in Section 2. We describe the method to extract the
constraints on $f_{\mathrm{NL}}$ using the 2pcf and 3pcf in configuration-
space in Section 3. In Section 4 we discuss the application of this method to
simulated DESI LRGs and present forecasts of DESI sensitivity to PNG. We
discuss the implications and challenges of this method in Section 5 and
conclude in Section 6.

##  2 Evaluation of 2pcf and 3pcf using ConKer

We compute the 2pcf and 3pcf using the ConKer algorithm described in Brown et
al. (2022). Employing techniques developed in Brown et al. (2021), this
algorithm estimates correlation functions by convolving the matter density
field with spherical shell kernels. Here, we provide an overview of the method
and introduce the cases of relevance for this study.

For an arbitrary matter tracer, let $\rho(\vec{r})$ be the tracer density at
position $\vec{r}$. The expected density of a random distribution of that
tracer is $\bar{\rho}(\vec{r})$. We define the deviation in density from
average as

$\Delta(\vec{r})=\rho(\vec{r})-\bar{\rho}(\vec{r})\ .$ (2)

For some correlation order $n$, let us consider all possible configurations of
$n$ points, referred to as $n$-plets, with one vertex at point $0$,
characterized by a vector $\vec{r}$, and the other vertices are defined by
vectors $\vec{r_{i}}$, where $i=1,...n-1$ (see Fig. 1). Each point is
connected to point $0$ by the vector $\vec{s_{i}}=\vec{r_{i}}-\vec{r}$. The
$n$-point correlation functions ($n$pcfs) characterize the excess of a
particular $n$-plet over that of a random distribution of tracers.

Figure 1: An illustration of the two ($0-1$) and three ($0-1-2$) point
correlation function evaluation using ConKer. The blue spherical kernels of
radii $s$ for the 2pcf and $s_{1}$ and $s_{2}$ for the 3pcf are convolved with
the matter distributions. The integration over $\vec{r}$ is performed by
scanning the position of the kernel center $0$ over the surveyed volume.

In this approach, the $n$pcfs are given by

$\begin{split}\xi_{n}(\vec{s}_{1},..\vec{s}_{(n-1)})=\frac{1}{R_{n}^{0}}\int\Delta(\vec{r})\Delta(\vec{r}_{1})..\Delta(\vec{r}_{(n-1)})d\vec{r}\
,\end{split}$ (3)

where the integration over $\vec{r}$ implies all possible positions of point
$0$ in the surveyed volume. The normalization $R_{n}^{0}$ is evaluated using
the same integration performed over the random distribution of tracers.

The angle-averaged $n$pcfs (or “monopole” in terms of the Legendre expansion)
are defined as

$\begin{split}\xi_{n}^{0}(s_{1},..s_{(n-1)})=\frac{1}{R_{n}^{0}}\int\Delta(\vec{r})\Delta(\vec{r}_{1})..\Delta(\vec{r}_{(n-1)})d\hat{\vec{s}}_{1}..d\hat{\vec{s}}_{(n-1)}d\vec{r}\
,\end{split}$ (4)

where the integration over the unit vectors $\hat{\vec{s}}_{i}$ implies
averaging over the orientations of $\vec{s_{i}}$. With ConKer, these integrals
are evaluated by convolving spherical shell kernels of radii $s_{i}$ with the
matter density field. After averaging over angles, $n$pcfs depend only on the
absolute distances $s_{i}$, and not on the orientation of the corresponding
vectors.

While ConKer can be used to evaluate the Legendre multipoles of $n$pcfs of
arbitrary order, in this study we consider only the monopoles of the 2pcf and
3pcf. Thus, for brevity we drop the superscript $0$ in the notation for the
$n$pcfs.

The distances $s$ are binned in $N_{\mathrm{bin}}$ bins. Then, the 2pcf can be
presented as an $N_{\mathrm{bin}}$-dimensional vector, while the 3pcf is a
symmetric $N_{\mathrm{bin}}\times N_{\mathrm{bin}}$-dimensional matrix. We
refer to the special case where $s_{1}=s_{2}$ as the diagonal 3pcf, because it
corresponds to the elements that sit along the main diagonal of this matrix.
We convert the part of the matrix below the diagonal into a vector of length
$(N_{\mathrm{bin}}^{2}-N_{\mathrm{bin}})/2$. We define the observable vector
$O_{i}$ as the concatenated vectors of the 2pcf, diagonal 3pcf, and non-
diagonal 3pcf,

$O_{i}=[\xi_{2}(s),\xi_{3}(s_{1}=s_{2}),\xi_{3}(s_{1}>s_{2})]\ .$ (5)

##  3 Description of the method

###  3.1 Theoretical Considerations

The relative deviation in density from the average cosmological matter density
is defined as

$\delta(\vec{r})=(\rho(\vec{r})-\bar{\rho})/\bar{\rho}\ .$ (6)

The matter overdensity $\delta_{\phi}$ associated with the primordial
gravitational potential $\Phi$ is

$\delta_{\phi}=\alpha\Phi\ .$ (7)

The expression for $\alpha$, a function of the wave-number $k$ associated with
density perturbations on scales $s\sim k^{-1}$ and the redshift $z$ is

$\alpha(k,z)=\frac{2D(z)}{3\Omega_{m,0}}T(k)\frac{g(z_{\mathrm{rad}})}{g(0)}(kd_{H})^{2}\
,$ (8)

where $D(z)$ is the growth parameter normalized to 1 at $z=0$. The factor
$g(z_{\mathrm{rad}})/g(0)$ takes into account the difference of this
normalization with respect to the early time matter-dominated epoch. The
Hubble distance, $d_{H}=c/H_{0}\approx 3000\ h^{-1}\mathrm{Mpc}$, is evaluated
using the speed of light $c$ and the present day Hubble’s constant $H_{0}$.
The reduced Hubble constant is defined as $h=H_{0}/(100\ \mathrm{km/s/Mpc})$,
and $\Omega_{m,0}$ is the present day relative fraction of the matter density.

The transfer function $T(k)$ is close to unity for scales
$s>L_{0}=16(\Omega_{m,0}h^{2})^{-1}\approx 3\ h^{-1}\mathrm{Mpc}$, which are
typically considered in LSS clustering analyses. All terms preceding the
dimensionless scale-dependent parameter $(kd_{H})^{2}$ are also at the order
of unity; hence, for typically relevant scales $s\ll d_{H}$, the factor
$\alpha^{-1}$ can be considered small.

Cosmological surveys do not directly observe the distribution of matter,
mostly comprised of dark matter (DM). The degree to which an arbitrary tracer
$t$ adheres to the underlying matter distribution is characterized by the bias
$b_{t}$ (Desjacques et al., 2018). At large scales, this is described by a
linear relation between the matter ($\delta$) and tracer ($\delta_{t}$)
overdensities:

$\delta_{t}(z)=b_{t}(z)\delta(z)\ .$ (9)

For Gaussian initial conditions, and at linear scales, the bias can be
considered scale-independent and equal to

$b_{t}=\tilde{b}_{1t}(z)=b_{1t}+f\mu^{2}\ ,$ (10)

where in the first term the redshift dependence of the linear bias of a given
tracer is parametrized as $b_{1t}=b_{0t}/D(z)$ and the second term is the so-
called Kaiser term (Kaiser, 1987), which accounts for redshift space
distortions (RSD) caused primarily by the peculiar motion of the tracers.
Here, $f=-d\log D(z)/d\log(1+z)$ is the logarithmic growth rate, and $\mu$ is
the cosine of the angle with respect to the line of sight. Defining the linear
bias $\tilde{b}_{1t}$ to include RSD via the Kaiser factor is similar to the
formalism used in Castorina et al. (2019). The effects of RSD on two point and
three point clustering have been extensively studied (Slepian & Eisenstein,
2017). In the monopoles of the 2pcf and 3pcf, RSD affects an overall
renormalization of the linear bias. Importantly, this effect is scale
independent.

In the presence of local PNG, the bias acquires scale dependence encoded in
the factor $\alpha^{-1}(z,k)$ of the form

$b_{t}(z,k)=\tilde{b}_{1t}(z)+b_{\phi t}(z)f_{\mathrm{NL}}\alpha^{-1}(z,k)\ ,$
(11)

where $b_{\phi t}$ is a bias related to the gravitational potential, a.k.a.
PNG bias, which we note is degenerate with $f_{\mathrm{NL}}$. A typical
parameterization for $b_{\phi t}$ is given by the relation

$b_{\phi t}=2\delta_{c}(b_{1t}-p_{t})\ ,$ (12)

where $\delta_{c}=1.69$ is the critical density of spherical collapse, and
$p_{t}$ is a potentially tracer-dependent parameter describing the formation
of said tracer in primordial overdensities (Barreira, 2020). In LSS analyses
the common assumption for galaxies is that $p_{t}=1$. While we do not fix the
value of $p$, we do impose a corresponding prior with 1 as the central value.

Recently, Barreira (2022) and Lazeyras et al. (2023) have found variation in
$b_{\phi t}$ due to varying galaxy or secondary halo properties like
concentration. Additionally, Ross et al. (2017) and Rezaie et al. (2021)
showed how observational systematics in LSS surveys (the Sloan Digital Sky
Survey in this case) induce scale-dependent effects on the two-point
clustering of cosmological tracers, regardless of PNG. It is therefor
necessary to validate and confirm any potential detection of PNG in the most
robust way possible to rule out any effects that mimic PNG.

The $n$pcfs of any tracer $\xi_{n}$ depend on the density variation to the
$n^{\mathrm{th}}$ power and are thus approximately proportional to
$b_{t}^{n}$. The dependence on $f_{\mathrm{NL}}$ can be expanded to the first
order of the small parameter $\alpha^{-1}$, and as a result to the first order
in $f_{\mathrm{NL}}$:

$\xi_{n}(s,f_{\mathrm{NL}})=\xi^{b_{1}}_{n}(s)+\xi^{\mathrm{PNG}}_{n}(s,f_{\mathrm{NL}})\
,$ (13)

where $\xi^{b_{1}}_{n}(s)$ are the $n$pcfs of the underlying matter
distribution in the absence of PNG ($f_{\mathrm{NL}}=0$) and
$\xi^{\mathrm{PNG}}_{n}(s,f_{\mathrm{NL}})$ is the contribution due to PNG
proportional to the first power of $f_{\mathrm{NL}}$.

The difference between $n$pcfs corresponding to non-zero and zero
$f_{\mathrm{NL}}$ is linear in $f_{\mathrm{NL}}$ with coefficients

$A_{n}(s,b_{1t})=(\xi_{n}(s,b_{1t},f_{\mathrm{NL}})-\xi_{n}(s,b_{1t},0))/f_{\mathrm{NL}}\
.$ (14)

We use this linearity to interpolate the value of the $n$pcfs at a given scale
$s$ for arbitrary values of $f_{\mathrm{NL}}$. In Appendix B we discuss tests
of linearity and show the scales at which this assumption breaks down. The
coefficients $A_{n}$ are determined in each bin in $s$ using DM halo
simulations ($t=h$) with $f_{\mathrm{NL}}=0,100$. This describes the size of
the change in the scale dependent bias for a given $f_{\mathrm{NL}}$.

###  3.2 Statistical model

To construct a model of $n$pcfs of an arbitrary tracer, we use simulations of
the tracer without PNG, referred to as the fiducial simulations, and halos
from DM simulations with varied $f_{\mathrm{NL}}$. For the simulated dark
matter halos, we use an approximated N-body, or “FastPM” technique (Feng et
al., 2016) to generate halos in $(3\ h^{-1}\mathrm{Gpc})^{3}$ boxes. We refer
to these simulations as FastPM-L3 (where L3 indicates the length of the
simulation box in $h^{-1}\mathrm{Gpc}$). We use independent samples of the
tracer of interest, referred to as the test simulations, to evaluate the
performance of the model. At the final stage of the analysis the test
simulation will be replaced by data. RSD are included for all tracer
simulations, but not for the FastPM-L3 halos, which is reflected in the
definition of the corresponding bias parameters.

Using the fiducial simulations, we evaluate the fiducial $n$pcfs,
$\xi^{\mathrm{fid}}_{n}(s)$. These simulations are of the same tracer type as
the test sample, and approximately match their spatial distribution. It is
possible that the linear bias in the fiducial simulations,
$\tilde{b}_{1t}^{\mathrm{fid}}$, does not match the value observed in the test
samples, $\tilde{b}_{1t}$. To allow for this difference, we scale
$\xi^{b_{1}}_{n}(s)$ by:

$\xi^{b_{1}}_{n}(s)=\left(\frac{\tilde{b}_{1t}}{\tilde{b}_{1t}^{\mathrm{fid}}}\right)^{n}\xi^{\mathrm{fid}}_{n}(s)\
.$ (15)

To construct $\xi^{\mathrm{PNG}}_{n}(s,f_{\mathrm{NL}})$ we use coefficients
$A_{n}(s,b_{1h})$ derived from the FastPM-L3 simulations. Since halos, $h$,
and tracers, $t$, have different values of the linear bias — $b_{1h}$ and
$\tilde{b}_{1t}$ respectively, we scale
$\xi^{\mathrm{PNG}}_{n}(s,f_{\mathrm{NL}})$ according to the redshift
dependence of $\alpha$ given by Eq. 8. The PNG bias for the FastPM-L3 halos
and tracers could also differ, which is reflected in the different values of
parameter $p$: $p_{h}$ and $p_{t}$, respectively. Thus,
$\xi^{\mathrm{PNG}}_{n}(s,f_{\mathrm{NL}})$ has the following form:

$\xi^{\mathrm{PNG}}_{n}(s,f_{\mathrm{NL}})=\left(\frac{\tilde{b}_{1t}}{b_{1h}}\right)^{n}A_{n}(s)r_{\phi}f_{\mathrm{NL}}\
.$ (16)

The difference between the ratio of the PNG to the linear bias for the
FastPM-L3 halos compared to the tracers is encoded in the factor $r_{\phi}$,
where

$r_{\phi}=\frac{b_{1h}}{\tilde{b}_{1t}}\frac{b_{\phi t}}{b_{\phi
h}}=\frac{b_{1h}}{\tilde{b}_{1t}}\frac{b_{1t}-p_{t}}{b_{1h}-p_{h}}\ .$ (17)

We note that the linear tracer bias within the definition of $b_{\phi t}$
appears without the Kaiser term. The relationship between $b_{1t}$ and
$\tilde{b}_{1t}$ is given by the Eq. 10. The Kaiser factor was decomposed into
a Legendre basis in Eq. 2.19 of Castorina et al. (2019). Following this, we
relate the two linear biases for the monopole term as $\tilde{b}_{1t}\approx
b_{1t}+f/3$. The logarithmic growth rate may be estimated from the redshift
dependent matter density given by $f\approx\Omega_{m}(z)^{0.55}$. Here, the
$z$-dependent relative fraction of the matter density is

$\Omega_{m}(z)=\frac{\Omega_{m,0}(z+1)^{3}}{\Omega_{m,0}(z+1)^{3}+\Omega_{k,0}(z+1)^{2}+\Omega_{\Lambda,0}}\
,$ (18)

which depends on the present day matter, curvature, and cosmological constant
relative densities $\Omega_{m,0}$, $\Omega_{k,0}$, and $\Omega_{\Lambda,0}$
respectively.

The final expressions for the expected $n$pcfs have the following form:

$\xi_{n}(s,f_{\mathrm{NL}})=\left(\frac{\tilde{b}_{1t}}{\tilde{b}_{1t}^{\mathrm{fid}}}\right)^{n}\xi^{\mathrm{fid}}_{n}(s)+\left(\frac{\tilde{b}_{1t}}{b_{1h}}\right)^{n}A_{n}(s)r_{\phi}f_{\mathrm{NL}}\
.$ (19)

This equation provides the expected values of the observable vector $O_{i}$,
defined in Eq. 5, which we denote as $\tilde{O_{i}}$. We use the observed and
expected values of $O_{i}$ to construct a test statistic $\chi^{2}$ that
depends on the parameters of interest (POIs): $f_{\mathrm{NL}}$ and the linear
bias $b_{0t}$, and on a set of nuisance parameters $\theta$:

$\theta=[b_{0h},b_{0t}^{\mathrm{fid}},p_{h},p_{t}]\ .$ (20)

The nuisance parameters are constrained using Gaussian priors with central
values given by a vector $\Theta$ with Gaussian widths $\sigma$. The test
statistic is then defined

$\chi^{2}=\sum_{ij}\left[O_{i}-\tilde{O}_{i}\right]^{T}C^{-1}_{ij}\left[O_{j}-\tilde{O}_{j}\right]+\sum_{f}\frac{(\theta_{f}-\Theta_{f})^{2}}{(\sigma_{f})^{2}}\
,$ (21)

where $C_{ij}$ denotes the covariance matrix. The first summation over indices
$i$ and $j$ is performed over bins in $s$ and the second over nuisance
parameters $\theta_{f}$.

$\chi^{2}$ is minimized to find the optimal values of the POIs in the test
sample. In the following, we discuss the application of this method to the
specific case of DESI luminous red galaxies, LRGs.

##  4 Application to simulated sample of DESI LRGs

###  4.1 Simulations

All simulations employed in this study are based on a flat $\Lambda$CDM
cosmology. Simulations of LRGs are generated with $\Omega_{m,0}=0.315$ and all
FastPM-L3 halo simulations have $\Omega_{m,0}=0.3089$. Due to the difference
in the present day value of the matter density, the second term of Eq. 19 is
scaled by their ratio. For the fiducial simulations, we use an ensemble of
1000 mocks based on the effective Zel’dovich (EZ) approximation (Chuang et
al., 2015). For each realization, northern galactic cap (NGC) and southern
galactic cap (SGC) catalogs are generated according to the angular footprint
of the complete DESI Year-5 LRG survey. The redshift distribution $n(z)$ is
subsampled to match that of the LRGs in the survey validation catalog (DESI
Collaboration et al., 2023a, b). We refer to these simulations as subsampled
Year-5 DESI LRGs (SY5). Each realization consists of approximately 2.70M
objects in the combined NGC and SGC. The cosmic variance observed across this
ensemble of simulations is expected to be larger than the full DESI Year-5 LRG
survey.

DESI will, however, be releasing a partial LRG catalog as a part of the Year-1
(Y1) data release. We also use an ensemble of 1000 EZ mocks corresponding to
the Y1 angular and redshift coverage, allowing us to forecast the sensitivity
to PNG achievable using only the Y1 LRG survey. Each realization of these
simulations consists of approximately 3.27M objects, matching the DESI Y1
LRGs.

Figure 2:  The angular coverage (top left) and the redshift distribution
(bottom right) of the DESI SY5 (grey) and Y1 (cyan) LRG simulations.

As our test samples, we consider an ensemble of higher resolution N-body
simulations. They consist of 25 realizations of AbacusSummit (Maksimova et
al., 2021; Garrison et al., 2021) LRGs, distributed about halos using a halo
occupation distribution (HOD) model defined in Bose et al. (2022); Yuan et al.
(2023). Like the SY5 EZ mocks, the AbacusSummit LRGs are convolved with the
NGC and SGC DESI Year-5 angular window and subsampled to match the survey
validation redshift distribution.

For both the AbacusSummit and SY5 EZ mock catalogs, we include only tracers
within the redshift range $0.4<z_{g}<1.4$, where $g$ denotes galaxies. The
angular coverage and redshift distributions of these simulations are shown in
Fig. 2.

In the FastPM-L3 halo simulations, $N=2560^{3}$ DM particles of mass
$m_{p}=1.4\times 10^{11}\ h^{-1}M_{\odot}$ in a $(3\ h^{-1}\mathrm{Gpc})^{3}$
volume were evolved from $z=99$ to $z=1$ using 50 linearly spaced time steps.
100 realizations with $f_{\mathrm{NL}}=0,\ 100$ were generated with matched
phases, such that individual realizations at both $f_{\mathrm{NL}}=0$ and
$f_{\mathrm{NL}}=100$ share the same initial conditions (Avila & Adame, 2023).
This results in two ensembles of DM halo simulations, one with and one without
the presence of PNG. Catalogs were created by selecting halos with
$M>1.15\times 10^{13}\ h^{-1}M_{\odot}$. Due to their size, these simulations
are ideal for detecting the effects of the linear change in $n$pcfs due to
PNG. Although they are likely too small to include larger modes leading to
non-linear effects, our model assumes linearity in $f_{\mathrm{NL}}$ as
discussed in Appendix B.

In addition to the FastPM-L3 halo simulations, we also use a set of four
larger FastPM simulations generated in $(5.52\ h^{-1}\mathrm{Gpc})^{3}$ boxes.
We refer to these as the FastPM-L5.52 halo simulations. Their use in this
study and a description of their properties are found in Appendices B & C. A
complete summary of all simulations used may be found in Table 1.

Name | $f_{\mathrm{NL}}$ | Tracer | Geometry | Volume | $N_{\mathrm{real}}$ | Purpose
---|---|---|---|---|---|---
|  |  |  | [($h^{-1}\mathrm{Gpc})^{3}$] |  |
EZ Mocks SY5 | 0 | LRGs | Subsampled DESI Y5 cut sky | 35.7 | 1000 | Fiducial model
EZ Mocks Y1 | 0 | LRGs | DESI Y1 cut sky | 26.3 | 1000 | Fiducial model
AbacusSummit | 0 | LRGs | Subsampled DESI Y5 cut sky | 35.7 | 25 | Validation tests
FastPM-L3 | 0, 100 | DM Halos |  $L=3\ h^{-1}\mathrm{Gpc}$ box | 27 | 100 | Estimate scale-dependant bias
FastPM-L5.52 | -25, 0, 12, 25 | DM Halos |  $L=5.52\ h^{-1}\mathrm{Gpc}$ box | 168.2 | 1 | Validation tests and linearity checks
Table 1: Descriptions of all simulations used, including the values of
$f_{\mathrm{NL}}$, the type of tracer, the geometry of the catalogs, their
volume, and their purpose in this study. $N_{\mathrm{real}}$ is the number of
realizations of the simulation for each value of $f_{\mathrm{NL}}$.

###  4.2 Fiducial $n$pcfs and Covariance Matrix

The fiducial $n$pcfs are evaluated using the average over the ensemble of EZ
mocks and are shown for the SY5 simulations in Fig. 3. The top panel shows the
2pcf and the diagonal part of the 3pcf. In the bottom panel of Fig. 3, we show
the 3pcf as a function of the triangle index. For comparison we also present
the results of the AbacusSummit LRGs, used as test samples. The $n$pcf
uncertainties, evaluated using the covariance matrices, are also shown for the
SY5 LRGs and the Y1 LRGs.

Figure 3: Top: The average diagonal $n$pcfs $\xi_{n}^{\mathrm{di}}$ when
$n=2,3$ for the SY5 EZ mocks (grey) and AbacusSummit simulation (green), as a
function of the distance scale, $s$. The shaded region shows uncertainties
estimated from the covariance matrix for the SY5 EZ mocks (grey) and Y1 EZ
mocks (cyan). The thin green lines show the $n$pcfs for each realization of
the AbacusSummit simulations independently. The diagonal $n$pcfs are
multiplied by $\tilde{s}^{n}$, where $\tilde{s}=(s/100\ h^{-1}\mathrm{Mpc})$,
to emphasize features at large scales. Bottom: The 3pcf $\xi_{3}$ for the SY5
EZ mocks (grey) and AbacusSummit simulation (green), as a function of the
triangle index where $s_{1}>s_{2}$. The 3pcf is multiplied by
$\tilde{s}_{1}^{2}\tilde{s}_{2}^{2}$. In all panels, the dashed black
horizontal line denotes 0.

We evaluate the covariance matrix from both the ensemble of SY5 EZ mocks and
from the Y1 EZ mocks. The normalized correlation matrix for SY5 is shown in
Fig. 4. Bins corresponding to the 2pcf and diagonal 3pcf are located in the
bottom left of the matrix, and bins corresponding to the off-diagonal 3pcf are
located in the top right.

Figure 4: The correlation matrix $C_{ij}=C_{ij}/\sqrt{C_{ii}C_{jj}}$, derived
from the ensemble of SY5 EZ mocks.

###  4.3 Sensitivity to PNG

We use the FastPM-L3 halo simulations with different values of
$f_{\mathrm{NL}}$ to evaluate the sensitivity of the $n$pcfs to the presence
of PNG as presented in Fig. 5. The lower panels show the coefficients $A_{n}$,
which are computed from the deviation in the $n$pcfs for $f_{\mathrm{NL}}=100$
from that of $f_{\mathrm{NL}}=0$. We observe a particular sensitivity in the
2pcf where the effect on halo clustering from PNG is clearly visible.

Figure 5:  The $n$pcfs $\xi_{n}^{\mathrm{di}}$ for $n=2,3$ (upper panels) for
the ensemble of the FastPM-L3 simulated halos with $f_{\mathrm{NL}}=0$ (blue),
$100$ (red), and the estimate of $A_{n}$ (lower panels). The shaded regions
show uncertainties estimated from the covariance matrix of each ensemble. Top:
The diagonal $n$pcfs and $A_{n}$ multiplied by $\tilde{s}^{n}$ are shown as
functions of the distance scale, $s$. Bottom: 3pcf and $A_{3}$ multiplied by
$\tilde{s}_{1}^{2}\tilde{s}_{2}^{2}$ are shown as a function of the triangle
index. In all panels, the dashed black horizontal line denotes 0.

In all cases of the 3pcf, we find considerably less sensitivity to the
presence of PNG, and large uncertainties for these monopole configurations. We
do note, however, that the FastPM-L3 halos show a degree of PNG sensitivity at
small scales (see the left-most portion of the $A_{3}$ panel). At large
scales, the uncertainties are considerably larger than any deviation due to
the presence of PNG. The increased statistics of the full DESI sample will
reduce these uncertainties, thus enhancing the sensitivity to PNG at large
scales.

###  4.4 Linear and PNG biases

The priors on the linear biases $b_{0h}$ and $b_{0g}^{\mathrm{fid}}$ are
derived using the FastPM-L3 halo and fiducial simulations respectively. To do
this, we measure the ratio of the 2pcf for each tracer to the 2pcf predicted
by linear theory. The theoretical 2pcf was computed using the open source
cosmology toolkit nbodykit (Hand, 2018). For both cases, we estimate the
linear bias by performing a fit according to:

$b_{1}=(\xi_{2}/\xi_{2}^{\mathrm{lin}})^{1/2}\ .$ (22)

We show the results of this fit and the measured 2pcfs again in Fig. 6 for the
FastPM-L3 halos and the SY5 EZ mocks. The fit is stable at the scales
considered in this analysis. For the SY5 LRGs, we measure $\tilde{b}_{1}$ with
the fit, remove the Kaiser term, and normalize to $b_{0g}$ to report the value
of $b_{0}$ in the legend of Fig. 6.

Here, $\xi_{2}^{\mathrm{lin}}$ is evaluated at the same redshift, or effective
redshift of the sample. For the FastPM-L3 halo simulations, the entire volume
has been evolved to $z_{h}=1$ and the effective redshift for the SY5 EZ mock
LRG simulations is $z_{g}=0.82$. These redshifts are used to compute $D(z)$ to
normalize each linear bias.

The most probable values of $b_{0h}$ and $b_{0g}^{\mathrm{fid}}$ arising from
the 2pcf fits are chosen as the central values of the Gaussian priors, and we
choose a width of $\sigma=0.04$ for the FastPM-L3 halos, and $\sigma=0.055$
for the LRGs, reflecting the uncertainty in our estimates. In this study we
choose $1$ as the central value of the Gaussian prior on $p_{h}$ and $p_{g}$
with a width of $0.1$ for both.

Figure 6: Left: (Upper panel) The 2pcf $\xi_{2}$ multiplied by $\tilde{s}^{2}$
for the ensemble of FastPM-L3 halo simulations with $f_{\mathrm{NL}}=0$ (solid
line), the uncertainties (shaded region), and the theoretical 2pcf computed
from linear theory (dashed line) scaled by the fit value of the linear bias,
$[b_{0h}D(z_{h})]^{2}$, as a function of the distance, $s$. (Bottom panel) the
residual between the observed 2pcf and the scaled theoretical 2pcf. Right: The
same for the ensemble of SY5 EZ mocks and the fiducial linear galaxy bias,
$b_{0g}^{\mathrm{fid}}$.

###  4.5 Test samples

To explore the projected sensitivity of the SY5 DESI LRGs, we perform two
tests. First, we validate our model by measuring the parameters using the
average $n$pcfs of the AbacusSummit LRGs as our observation. We treat this as
an observation with the covariance matrix derived from the SY5 EZ mocks.

Secondly, we generate toy data randomly distributed around the values of the
expected $n$pcfs given in Eq. 19, with uncertainties predicted by the
covariance matrix of the SY5 EZ mocks. While the AbacusSummit simulations were
generated with $f_{\mathrm{NL}}=0$, toy experiments can be produced with
arbitrary values of $f_{\mathrm{NL}}^{\mathcal{T}}$ and
$b_{0g}^{\mathcal{T}}$, thus allowing us to calibrate the method on sample
observations with PNG. We may also vary the values of the nuisance parameters,
$\theta$, in the generation of toy model data.

###  4.6 Results of the fit

We compute the 2pcf using ConKer in $N_{\mathrm{bin}}=34$ bins of width
$\Delta s=10$ $h^{-1}$Mpc from $s=50\ h^{-1}\mathrm{Mpc}\ \mbox{---}\ 380$
$h^{-1}$Mpc, unless stated otherwise, for all simulation samples. The 3pcf is
similarly computed in in $N_{\mathrm{bin}}=16$ bins of width $\Delta s=10$
$h^{-1}$Mpc from $s=50\ h^{-1}\mathrm{Mpc}\ \mbox{---}\ 200$ $h^{-1}$Mpc. The
combination of these two measurements results in a vector with
$N_{\mathrm{obs}}^{\mathrm{2pcf}}+N_{\mathrm{obs}}^{\mathrm{3pcf}}=170$ total
bins. Our choices of scales are motivated by the limitation of the linearity
assumption discussed in Appendix B.

Figure 7: Marginalized distributions of the POIs and nuisance parameters
$\theta$, for SY5 toy model $n$pcfs (grey), and for the mean of the
AbacusSummit simulations (green). The dark and light regions represent
$1\sigma$ and $2\sigma$ contours, respectively. The most probable values and
$1\sigma$ CL are labeled for each parameter above the respective panel. Priors
are shown by dashed lines in the one-dimensional histograms. The stars in the
$(f_{\mathrm{NL}}$, $b_{0g})$ panel show the most probable values for each
realization of the AbacusSummit mocks. Figure 8: A distribution of the most
probable values of $f_{\mathrm{NL}}$ for each realization of the AbacusSummit
mocks (green). The histogram bins are fit to a Gaussian (dashed line) where
the width and central value of the fit are given in the legend. Figure 9:
Marginalized distributions of the POIs and nuisance parameters $\theta$, for
SY5 toy model $n$pcfs with $f_{\mathrm{NL}}^{\mathcal{T}}=25$,
$b_{0g}^{\mathcal{T}}=1.34$ (purple) and $f_{\mathrm{NL}}^{\mathcal{T}}=0$,
$b_{0g}^{\mathcal{T}}=1.21$ (orange). The dark and light regions represent
$1\sigma$ and $2\sigma$ contours, respectively. The most probable values and
$1\sigma$ CL are labeled for each parameter above the respective panel. Priors
are given by dashed lines in the one-dimensional histograms. The truth values
indicated by the vertical and horizontal dotted lines correspond to the values
of the parameters used to generate these toy model data.

To sample the test statistic described in Eq. 21, we use a Markov Chain Monte-
Carlo (MCMC) algorithm with 50 walkers across 10K steps, and remove the
initial 500 steps of burn-in. We show the marginalized distributions over the
POIs and nuisance parameters in Fig. 7. For the SY5 toy data, generated with
$f_{\mathrm{NL}}^{\mathcal{T}}=0$ and $b_{0g}^{\mathcal{T}}=1.34$, the
expected uncertainty on $f_{\mathrm{NL}}$ is approximately 16. For the
validation case (the AbacusSummit simulations) we observe that
$f_{\mathrm{NL}}$ is consistent with the expected value (namely
$f_{\mathrm{NL}}=0$) with a slightly larger uncertainty of approximately 18.
We note that for these two cases, the value of $b_{0g}$ differs, which can be
expected for different simulations. The ability to correctly measure
$f_{\mathrm{NL}}$ is retained even when the tracer bias differs from the
fiducial value.

To understand the degree to which this measurement is limited by cosmic
variance, we repeat the procedure (measuring the most probable values of the
parameters) for each realization of the AbacusSummit mocks. We present these
measurements of the POIs as stars in the $(f_{\mathrm{NL}}$, $b_{0g})$ panel
of Fig. 7. We observe marginally larger variation in the measurements
performed across each realization of the simulation, as compared to the
estimated uncertainty. We show a distribution over the measured values of
$f_{\mathrm{NL}}$ for each realization in Fig. 8. The central value of the
Gaussian fit to the histogram is within the $\sigma/2$ of the true value. The
variation between realizations of the AbacusSummit mocks is slightly larger
than the projected $f_{\mathrm{NL}}$ uncertainty (32 vs expected 19), which is
not surprising given the small number of AbacusSummit mocks.

Figure 10: Upper left: The measured values of $f_{\mathrm{NL}}$, with
$1\sigma$ uncertainties, for SY5 toy model data where
$f_{\mathrm{NL}}^{\mathcal{T}}$ is varied from $-50$ to $50$. The dashed line
has a slope of 1 and intercept of 0. Upper right: The size of the
uncertainties on the $f_{\mathrm{NL}}$ measurement,
$\sigma_{f_{\mathrm{NL}}}$, for SY5 toy model data where
$b_{0g}^{\mathcal{T}}$ is varied from $1.01$ to $1.68$. The uncertainties are
given as a function of the ratio between the toy model linear bias, and the
toy model fiducial linear bias,
$[b_{0g}/b_{0g}^{\mathrm{fid}}]^{\mathcal{T}}$. Lower left: The measured
values of $f_{\mathrm{NL}}$, with $1\sigma$ uncertainties, for SY5 toy model
data where $f_{\mathrm{NL}}^{\mathcal{T}}=25$ and the values of the $p$
parameters have been modified. Both $p_{h}^{\mathcal{T}}$ (cyan) and
$p_{g}^{\mathcal{T}}$ (pink) were varied from $0.5$ to $2.0$. The priors on
$p_{h}$ and $p_{g}$ have not been modified. The dashed line denotes
$f_{\mathrm{NL}}=25$. Lower right: The size of the uncertainties on the
$f_{\mathrm{NL}}$ measurement, $\sigma_{f_{\mathrm{NL}}}$, for SY5 toy model
data where $p$ is varied. This panel shows the magnitude of the uncertainties
in the lower left panel.

We use an SY5 toy model to produce 2pcfs and 3pcfs for two different cases -
($f_{\mathrm{NL}}^{\mathcal{T}}=25$ and $b_{0g}^{\mathcal{T}}=1.34$) and
($f_{\mathrm{NL}}^{\mathcal{T}}=0$, $b_{0g}^{\mathcal{T}}=1.21$). We show the
marginalized distributions over the POIs and nuisance parameters for these two
additional cases in Fig. 9. The POIs were correctly evaluated in both cases.
We note that the uncertainty on $f_{\mathrm{NL}}$ increases with the decrease
of $b_{0g}$.

We perform additional systematic tests to understand the effects of the POIs
and nuisance parameters on the uncertainties on $f_{\mathrm{NL}}$. To ensure
we retain the ability to measure $f_{\mathrm{NL}}$ in the extreme cases, we
generate SY5 toy model data varying $f_{\mathrm{NL}}^{\mathcal{T}}$ from $-50$
to $50$. We show the measured $f_{\mathrm{NL}}$ vs. input
$f_{\mathrm{NL}}^{\mathcal{T}}$ in the upper left panel of Fig 10. The
magnitude of the uncertainties are consistent across this range, and the
measured $f_{\mathrm{NL}}$ is in agreement with the toy model value for all
cases. While this test is performed using only toy model data, we also verify
our model’s ability to accurately constrain $f_{\mathrm{NL}}$ for a variety of
non-zero PNG cases. This test is described in detail in Appendix C.

As shown in Fig. 9, a change in the linear galaxy bias relative to the
fiducial simulation results in a change in the size of the $f_{\mathrm{NL}}$
uncertainties. To further explore this dependence, we generate a set of SY5
toy model data with $f_{\mathrm{NL}}^{\mathcal{T}}=0$, and varied linear bias,
$b_{0g}^{\mathcal{T}}$. In the upper right panel of Fig 10, we show the
magnitude of the $f_{\mathrm{NL}}$ uncertainty as a function of the ratio
between the toy model linear bias and the toy model fiducial linear bias. It
is clear that tracers with a larger linear bias with respect to the fiducial
model offer improved constraining power on $f_{\mathrm{NL}}$.

While we have chosen to center our priors on $p_{h}$ and $p_{g}$ on 1, this
may not correspond to the true value for either tracer. To explore the
sensitivity to the values of $p$, we generated toy model data with varied
$p_{h}^{\mathcal{T}}$ and $p_{g}^{\mathcal{T}}$. Since this most dramatically
affects the sensitivity when $f_{\mathrm{NL}}$ is nonzero, these toy model
data sets were generated with $f_{\mathrm{NL}}^{\mathcal{T}}=25$. We apply the
same fitting procedure, with the priors centered on 1. We show the measured
$f_{\mathrm{NL}}$ and its uncertainty in the lower left and right panels of
Fig 10. We retain the ability to measure $f_{\mathrm{NL}}$ to within $1\sigma$
of the true value when either $p_{g}$ or $p_{h}$ is varied by about 0.5. For a
larger variation, our measurement of $f_{\mathrm{NL}}$ begins to deviate
significantly from the true value. We also note that while this mismatch in
$p$ between the true value and the prior affects the size of the uncertainties
on $f_{\mathrm{NL}}$, this effect is small compared to the effect of varying
the linear galaxy bias with respect to the fiducial simulation. This effect is
also demonstrated by the shape of the $(f_{\mathrm{NL}}$, $p_{h})$ and
$(f_{\mathrm{NL}}$, $p_{g})$ contours in Fig. 9. If the priors on both $p$
parameters are accurate, we can expect larger uncertainties on
$f_{\mathrm{NL}}$ in the case of higher $p_{g}$ and lower $p_{h}$.

To summarize, the parameter that is most significantly correlated with
$f_{\mathrm{NL}}$ is $b_{0g}$, specifically the measured value relative to the
value in the fiducial simulation. The dependence on the values of $p_{h}$ and
$p_{g}$ is considerably less pronounced.

###  4.7 Advantages of large scales and higher correlation orders

By combining the information contained in the 2pcf and 3pcf, we achieve
greater sensitivity to the POIs than for either case alone. We demonstrate
this by constructing the test statistic for SY5 toy model data using the 2pcf
and 3pcf combined, as well as each separately. For all of these cases, we
generate toy observations with $f_{\mathrm{NL}}^{\mathcal{T}}=0$ and
$b_{0g}^{\mathcal{T}}=1.34$. We show the marginalized distributions of the
POIs for this test in Fig. 11. By itself, the 3pcf monopole lacks the ability
to meaningfully constrain the value of $f_{\mathrm{NL}}$. However, in
combination with the information contained in the 2pcf, we improve the
sensitivity.

Figure 11: Marginalized distributions of the POIs for SY5 toy model $n$pcfs, using test statistics constructed from the 2pcf and 3pcf combined (grey), 2pcf only (red) and 3pcf only (blue). The dark and light regions represent $1\sigma$ and $2\sigma$ contours, respectively. The most probable values and $1\sigma$ CL are labeled for each parameter above the respective panel. 2pcf scales | 3pcf scales | $N_{\mathrm{obs}}^{\mathrm{tot}}$ | $\sigma_{f_{\mathrm{NL}}}$ | $\sigma_{b_{0g}}$
---|---|---|---|---
[$h^{-1}$Mpc] | [$h^{-1}$Mpc] |  |  |
50 — 200 | 50 — 200 | 152 | 22 | 0.029
50 — 260 | 50 — 200 | 158 | 19 | 0.030
50 — 320 | 50 — 200 | 164 | 18 | 0.029
50 — 380 | 50 — 200 | 170 | 16 | 0.029
Table 2: $\sigma_{f_{\mathrm{NL}}}$ and $\sigma_{b_{0g}}$, the uncertainties
on $f_{\mathrm{NL}}$ and $b_{0g}$ respectively, extracted using the SY5 toy
model with $f_{\mathrm{NL}}^{\mathcal{T}}=0$ and $b_{0g}^{\mathcal{T}}=1.34$
by varying $s_{\mathrm{max}}$ of the 2pcf. $N_{\mathrm{obs}}^{\mathrm{tot}}$
is the total length of the combined 2pcf and 3pcf observable. Figure 12:
$\sigma_{f_{\mathrm{NL}}}$ and $\sigma_{b_{1g}}$, the uncertainties on
$f_{\mathrm{NL}}$ and $b_{0g}$ respectively, extracted from SY5 toy model data
with $f_{\mathrm{NL}}^{\mathcal{T}}=0$ and $b_{0g}^{\mathcal{T}}=1.34$ as a
function of a varying $s_{\mathrm{max}}$ in the vector of 2pcf observations.

Given that most of the sensitivity to $f_{\mathrm{NL}}$ is extracted from the
2pcf, we may also wonder how much of the constraining power of the 2pcf comes
from large scales. As cosmological surveys increase in volume, measurements of
$n$pcfs at these scales will improve in their precision. We can demonstrate
how galaxy clustering at large scales improves our own sensitivity to the POIs
by constructing the test statistic using a limited range of scales. We
generate SY5 toy model $n$pcfs, again corresponding to
$f_{\mathrm{NL}}^{\mathcal{T}}=0$ and $b_{0g}^{\mathcal{T}}=1.34$, but vary
the maximum clustering scale $s_{\mathrm{max}}$ of the 2pcf. We keep the
separation between bins, $\Delta s=10$ $h^{-1}$Mpc, and generate 2pcfs from
$s=50\ h^{-1}\mathrm{Mpc}\ \mbox{---}\ s_{\mathrm{max}}$. The scales used to
generate the 3pcf for these test are held constant ($s=50\ h^{-1}\mathrm{Mpc}\
\mbox{---}\ 200\ h^{-1}\mathrm{Mpc}$). The procedure is repeated for
$s_{\mathrm{max}}=200,260,320,380$ $h^{-1}$Mpc.

We show the resulting $1\sigma$ uncertainties on the POIs in Table 2 and Fig.
12. In the table, we round the uncertainties on $b_{0g}$ to an additional
decimal to show the minimal effect. Increasing the scales of the 2pcf does
little to further constrain the linear bias. The sensitivity is nearly scale
independent. The ability to measure $f_{\mathrm{NL}}$, however, depends
strongly on $s_{\mathrm{max}}$. As we include larger clustering scales in the
vector of 2pcf observations, we reduce the uncertainties on $f_{\mathrm{NL}}$.
The diminishing nature of this effect is due to increasing cosmic variance in
the $n$pcfs for surveys of the size used in this study. Thus, a larger survey
(or a sample of LRGs combined with other tracers) occupying larger volume
would offer even greater constraining power. Beyond the scales chosen for this
analysis, however, we would likely need to include higher order
$f_{\mathrm{NL}}$ terms.

Figure 13: Marginalized distributions of the POIs for toy model $n$pcfs
corresponding to DESI SY5 LRGs (grey) and DESI Y1 LRGs (cyan). The dark and
light regions represent $1\sigma$ and $2\sigma$ contours, respectively. The
most probable values and $1\sigma$ CL are labeled for each parameter above the
respective panel.

###  4.8 DESI Y1 LRG forecast

The DESI Y1 LRGs will cover a smaller volume than the full sample, as the
survey is not yet complete. Thus, $n$pcfs corresponding to the ensemble of Y1
EZ mocks, whose coverage is shown in Fig. 2, have larger variance. To forecast
the sensitivity of our model to $f_{\mathrm{NL}}$ in the Y1 LRG survey, we
repeat the procedure with the Y1 EZ mocks defining the fiducial model. We
generate 2pcf and 3pcf Y1 toy model data using the covariance of this
ensemble. In Fig. 13, we show the marginalized distributions of the POIs for
the Y1 and SY5 toy model data with $f_{\mathrm{NL}}^{\mathcal{T}}=0$. The
effective redshift of the Y1 sample differs from that of the SY5 simulations.
In the Y1 sample, we derive the Gaussian prior on $b_{0t}$ using $z_{g}=0.73$.
For the normalized linear bias, we find good agreement between the priors for
the SY5 and Y1 samples, which informs the toy model value,
$b_{0g}^{\mathcal{T}}=1.34$ for both. We compare the projected precision for
each, finding that for the Y1 LRGs, we can expect
$\sigma_{f_{\mathrm{NL}}}\approx 22$ compared to
$\sigma_{f_{\mathrm{NL}}}\approx 16$ for the SY5 toy model data. It is
reasonable to assume that this sensitivity will continue to improve with
larger volumes covered by future galaxy surveys.

##  5 Discussion

Any measurement of local PNG from LSS would benefit from a robust
confirmation, using two-point and three-point statistics in configuration and
Fourier space. We have provided a method to achieve that.

In this study, we limited the cases considered to the 2pcf and 3pcf monopole
term for comparison to existing Fourier space studies. ConKer could be
employed to measure $n$pcfs for higher correlation orders, which may be even
more sensitive to PNG. However, this approach presents several difficulties.
First, higher orders of correlation introduce additional PNG biases that may
have non-trivial values. This would require a careful reconsideration of the
second term in Eq. 19. The parameter $r_{\phi}$, which characterizes the ratio
of PNG biases between the DM halos and the tracers, may not be sufficient.
Furthermore, correlations of order $n>3$ may need to be parsed into their
connected and disconnected terms (Philcox et al., 2021), which may introduce
additional systematics in the case of non-zero PNG.

To apply the proposed method to observed LRGs, we also need to account for
survey systematics, which can induce a scale-dependent deviation in the galaxy
clustering signal. Imaging systematics from the DESI Legacy Imaging Surveys
(Zou et al., 2017; Dey et al., 2019; Schlegel et al., 2023) and survey
completeness (Schlafly et al., 2023; Kirkby et al., 2023; Myers et al., 2023)
are typically mitigated by weights, which should be introduced into the model
as additional nuisance parameters and marginalized. Spurious effects of the
spectroscopic pipeline (Guy et al., 2023; Bailey et al., 2023; Brodzeller et
al., 2023) will also be mitigated with redshift failure weights. If the
analysis considers scales below $s=50\ h^{-1}\mathrm{Mpc}$, we may also
introduce nuisance parameters corresponding to variations in the Halo
Occupancy Distribution (HOD). In addition, smaller scales are subject to
variation arising from the process of assigning DESI fibers to targets
(Raichoor et al., 2023a). These extensions of our model are the focus of
ongoing research, targeted for a separate publication.

Finally, we note that while this study has demonstrated the application of our
method to DESI LRGs, the model is applicable to an arbitrary matter tracer.
All that is required is an ensemble of fiducial simulations for that tracer to
evaluate the covariance matrix, and the expected clustering signal in the
absence of PNG (like the EZ mocks), and a choice of the parameter, $p$. These
simulations are generally less computationally expensive than full N-body
simulations, and often accompany LSS data catalog releases. Therefore, this
strategy is applicable to other DESI matter tracers, or those observed in
future surveys like Euclid.

##  6 Conclusions

In summary, we developed a method by which local PNG, parameterized by
$f_{\mathrm{NL}}$, is extracted from configuration space correlation
functions. Using the 2pcf and 3pcf of FastPM-L3 halo simulations to
characterize the scale dependent bias, we have estimated the sensitivity of
DESI LRGs to this signal. We find competitive constraints on local PNG using
toy model data, and confirm the robustness of this result using independent,
high fidelity simulations. With an appropriate treatment of observational
systematics, we believe this framework is applicable to current and future LSS
surveys like DESI or Euclid, and will provide additional insight into the
mechanisms of inflation.

## Acknowledgements

The authors would like to thank Z. Slepian, J. Hou, H. Gil-Marín, M. Wang, and
L. Samushia for their interest and insightful comments, and S. BenZvi, K.
Douglass, and J. Bermejo Climent for useful discussions.

Z. Brown and R. Demina acknowledge support from the U.S. Department of Energy
under the grant DE-SC0008475.0.

This material is based upon work supported by the U.S. Department of Energy
(DOE), Office of Science, Office of High-Energy Physics, under Contract No.
DE–AC02–05CH11231, and by the National Energy Research Scientific Computing
Center, a DOE Office of Science User Facility under the same contract.

Additional support for DESI was provided by the U.S. National Science
Foundation (NSF), Division of Astronomical Sciences under Contract No.
AST-0950945 to the NSF’s National Optical-Infrared Astronomy Research
Laboratory; the Science and Technology Facilities Council of the United
Kingdom; the Gordon and Betty Moore Foundation; the Heising-Simons Foundation;
the French Alternative Energies and Atomic Energy Commission (CEA); the
National Council of Science and Technology of Mexico (CONACYT); the Ministry
of Science and Innovation of Spain (MICINN), and by the DESI Member
Institutions: https://www.desi.lbl.gov/collaborating-institutions.

Any opinions, findings, and conclusions or recommendations expressed in this
material are those of the author(s) and do not necessarily reflect the views
of the U. S. National Science Foundation, the U. S. Department of Energy, or
any of the listed funding agencies.

The authors are honored to be permitted to conduct scientific research on
Iolkam Du’ag (Kitt Peak), a mountain with particular significance to the
Tohono O’odham Nation.

## Data Availability

All simulations of LRGs included in this study (EZ mocks and AbacusSummit
simulations) will be made available in upcoming data releases from the Dark
Energy Spectroscopic Instrument. Please contact the authors for access to the
collection of Fast-PM simulations. The data shown in this paper’s figures are
available in machine-readable format at https://zenodo.org/records/10794474.

## References

  * Alam et al. (2017) Alam S., et al., 2017, Monthly Notices of the Royal Astronomical Society, 470, 2617
  * Alexander et al. (2023) Alexander D. M., et al., 2023, AJ, 165, 124
  * Allende Prieto et al. (2020) Allende Prieto C., et al., 2020, Research Notes of the American Astronomical Society, 4, 188
  * Alvarez et al. (2014) Alvarez M., et al., 2014, Technical report, Testing Inflation with Large Scale Structure: Connecting Hopes with Reality. Brookhaven National Laboratory (BNL)
  * Avila & Adame (2023) Avila S., Adame A. G., 2023, Monthly Notices of the Royal Astronomical Society, 519, 3706
  * Bailey et al. (2023) Bailey et al. 2023
  * Barreira (2020) Barreira A., 2020, Journal of Cosmology and Astroparticle Physics, 2020, 031
  * Barreira (2022) Barreira A., 2022, Journal of Cosmology and Astroparticle Physics, 2022, 033
  * Bartolo et al. (2004) Bartolo N., Komatsu E., Matarrese S., Riotto A., 2004, Physics Reports, 402, 103
  * Biagetti (2019) Biagetti M., 2019, Galaxies, 7, 71
  * Biagetti et al. (2017) Biagetti M., Lazeyras T., Baldauf T., Desjacques V., Schmidt F., 2017, Monthly Notices of the Royal Astronomical Society, 468, 3277
  * Bose et al. (2022) Bose S., Eisenstein D. J., Hadzhiyska B., Garrison L. H., Yuan S., 2022, Monthly Notices of the Royal Astronomical Society, 512, 837
  * Brodzeller et al. (2023) Brodzeller A., et al., 2023, arXiv e-prints, p. arXiv:2305.10426
  * Brown et al. (2021) Brown Z., Mishtaku G., Demina R., Liu Y., Popik C., 2021, Astronomy & Astrophysics, 647, A196
  * Brown et al. (2022) Brown Z., Mishtaku G., Demina R., 2022, Astronomy & Astrophysics, 667, A129
  * Cabass et al. (2022) Cabass G., Ivanov M. M., Philcox O. H., Simonović M., Zaldarriaga M., 2022, Physical Review D, 106, 043506
  * Castorina et al. (2019) Castorina E., et al., 2019, Journal of Cosmology and Astroparticle Physics, 2019, 010
  * Chaussidon (2024) Chaussidon E., 2024, In prep.
  * Chaussidon et al. (2023) Chaussidon E., et al., 2023, ApJ, 944, 107
  * Chuang et al. (2015) Chuang C.-H., Kitaura F.-S., Prada F., Zhao C., Yepes G., 2015, Monthly Notices of the Royal Astronomical Society, 446, 2621
  * Collaboration (2020) Collaboration P., 2020, Astron. Astrophys, 641, A9
  * Cooper et al. (2022) Cooper A. P., et al., 2022, arXiv e-prints, p. arXiv:2208.08514
  * Creminelli et al. (2011) Creminelli P., D’Amico G., Musso M., Norena J., 2011, Journal of Cosmology and Astroparticle Physics, 2011, 038
  * Cunha et al. (2010) Cunha C., Huterer D., Doré O., 2010, Physical Review D, 82, 023004
  * D’Amico et al. (2022) D’Amico G., Lewandowski M., Senatore L., Zhang P., 2022, arXiv preprint arXiv:2201.11518
  * DESI Collaboration et al. (2016a) DESI Collaboration et al., 2016a, arXiv e-prints, p. arXiv:1611.00036
  * DESI Collaboration et al. (2016b) DESI Collaboration et al., 2016b, arXiv e-prints, p. arXiv:1611.00037
  * DESI Collaboration et al. (2022) DESI Collaboration et al., 2022, AJ, 164, 207
  * DESI Collaboration et al. (2023a) DESI Collaboration et al., 2023a, arXiv e-prints, p. arXiv:2306.06307
  * DESI Collaboration et al. (2023b) DESI Collaboration et al., 2023b, arXiv e-prints, p. arXiv:2306.06308
  * Dalal et al. (2008) Dalal N., Dore O., Huterer D., Shirokov A., 2008, Physical Review D, 77, 123514
  * Desjacques et al. (2018) Desjacques V., Jeong D., Schmidt F., 2018, Physics reports, 733, 1
  * Dey et al. (2019) Dey A., et al., 2019, AJ, 157, 168
  * Feng et al. (2016) Feng Y., Chu M.-Y., Seljak U., McDonald P., 2016, Monthly Notices of the Royal Astronomical Society, 463, 2273
  * Gangui et al. (1993) Gangui A., Lucchin F., Matarrese S., Mollerach S., 1993, arXiv preprint astro-ph/9312033
  * Garrison et al. (2021) Garrison L. H., Eisenstein D. J., Ferrer D., Maksimova N. A., Pinto P. A., 2021, Monthly Notices of the Royal Astronomical Society, 508, 575
  * Gualdi et al. (2021) Gualdi D., Gil-Marín H., Verde L., 2021, Journal of Cosmology and Astroparticle Physics, 2021, 008
  * Guy et al. (2023) Guy J., et al., 2023, AJ, 165, 144
  * Hahn et al. (2022) Hahn C., et al., 2022, arXiv e-prints, p. arXiv:2208.08512
  * Hand (2018) Hand N. e. a., 2018, Astron.J., 4, 156
  * Hotchkiss (2011) Hotchkiss S., 2011, Journal of Cosmology and Astroparticle Physics, 2011, 004
  * Kaiser (1987) Kaiser N., 1987, Monthly Notices of the Royal Astronomical Society, 227, 1
  * Kirkby et al. (2023) Kirkby et al. 2023
  * Komatsu & Spergel (2001) Komatsu E., Spergel D. N., 2001, Physical Review D, 63, 063002
  * Komatsu & Spergel (2002) Komatsu E., Spergel D. N., 2002, in The Ninth Marcel Grossmann Meeting: On Recent Developments in Theoretical and Experimental General Relativity, Gravitation and Relativistic Field Theories (In 3 Volumes). pp 2009–2010
  * Lan et al. (2023) Lan T.-W., et al., 2023, ApJ, 943, 68
  * Lazeyras et al. (2023) Lazeyras T., Barreira A., Schmidt F., Desjacques V., 2023, Journal of Cosmology and Astroparticle Physics, 2023, 023
  * Levi et al. (2013) Levi M., et al., 2013, arXiv e-prints, p. arXiv:1308.0847
  * Maksimova et al. (2021) Maksimova N. A., Garrison L. H., Eisenstein D. J., Hadzhiyska B., Bose S., Satterthwaite T. P., 2021, Monthly Notices of the Royal Astronomical Society, 508, 4017
  * Maldacena (2003) Maldacena J., 2003, Journal of High Energy Physics, 2003, 013
  * McDonald (2008) McDonald P., 2008, Physical Review D, 78, 123519
  * Miller et al. (2023) Miller T. N., et al., 2023, arXiv e-prints, p. arXiv:2306.06310
  * Mueller et al. (2019) Mueller E.-M., Percival W. J., Ruggeri R., 2019, Monthly Notices of the Royal Astronomical Society, 485, 4160
  * Mueller et al. (2021) Mueller E.-M., et al., 2021, arXiv preprint arXiv:2106.13725
  * Myers et al. (2023) Myers A. D., et al., 2023, AJ, 165, 50
  * Philcox et al. (2021) Philcox O. H. E., Slepian Z., Hou J., Warner C., Cahn R. N., Eisenstein D. J., 2021, ENCORE: Estimating Galaxy $N$-point Correlation Functions in $\mathcal{O}(N_{\rm g}^{2})$ Time (arXiv:2105.08722)
  * Raichoor et al. (2020) Raichoor A., et al., 2020, Research Notes of the American Astronomical Society, 4, 180
  * Raichoor et al. (2023a) Raichoor et al. 2023a
  * Raichoor et al. (2023b) Raichoor A., et al., 2023b, AJ, 165, 126
  * Rezaie et al. (2021) Rezaie M., et al., 2021, Monthly Notices of the Royal Astronomical Society, 506, 3439
  * Ross et al. (2017) Ross A. J., et al., 2017, Monthly Notices of the Royal Astronomical Society, 464, 1168
  * Roy et al. (2014) Roy F., Bouillot V. R., Rasera Y., 2014, Astronomy & Astrophysics, 564, A13
  * Ruiz-Macias et al. (2020) Ruiz-Macias O., et al., 2020, Research Notes of the American Astronomical Society, 4, 187
  * Schlafly et al. (2023) Schlafly E. F., et al., 2023, arXiv e-prints, p. arXiv:2306.06309
  * Schlegel et al. (2023) Schlegel et al. 2023
  * Silber et al. (2023) Silber J. H., et al., 2023, AJ, 165, 9
  * Slepian & Eisenstein (2017) Slepian Z., Eisenstein D. J., 2017, Monthly Notices of the Royal Astronomical Society, 469, 2059
  * Slosar et al. (2008) Slosar A., Hirata C., Seljak U., Ho S., Padmanabhan N., 2008, Journal of Cosmology and Astroparticle Physics, 2008, 031
  * Yèche et al. (2020) Yèche C., et al., 2020, Research Notes of the American Astronomical Society, 4, 179
  * Yuan et al. (2023) Yuan S., Hadzhiyska B., Abel T., 2023, Monthly Notices of the Royal Astronomical Society, 520, 6283
  * Zhou et al. (2020) Zhou R., et al., 2020, Research Notes of the American Astronomical Society, 4, 181
  * Zhou et al. (2023) Zhou R., et al., 2023, AJ, 165, 58
  * Zou et al. (2017) Zou H., et al., 2017, PASP, 129, 064101

##  Appendix A Author Affiliations

1 Department of Physics and Astronomy, University of Rochester, 500 Joseph C.
Wilson Boulevard, Rochester, NY 14627, USA
2 Instituto de Física Teórica (IFT) UAM/CSIC, Universidad Autónoma de Madrid,
Cantoblanco, E-28049, Madrid, Spain
3 Institut de Física d’Altes Energies (IFAE), The Barcelona Institute of
Science and Technology, Campus UAB, 08193 Bellaterra Barcelona, Spain
4 Lawrence Berkeley National Laboratory, 1 Cyclotron Road, Berkeley, CA 94720,
USA
5 Centro de Investigación Avanzada en Física Fundamental (CIAFF), Facultad de
Ciencias, Universidad Autónoma de Madrid, ES-28049 Madrid, Spain
6 Physics Dept., Boston University, 590 Commonwealth Avenue, Boston, MA 02215,
USA
7 NSF’s NOIRLab, 950 N. Cherry Ave., Tucson, AZ 85719, USA
8 Department of Physics & Astronomy, University College London, Gower Street,
London, WC1E 6BT, UK
9 Institute for Computational Cosmology, Department of Physics, Durham
University, South Road, Durham DH1 3LE, UK
10 Instituto de Física, Universidad Nacional Autónoma de México, Cd. de México
C.P. 04510, México
11 Department of Physics & Astronomy and Pittsburgh Particle Physics,
Astrophysics, and Cosmology Center (PITT PACC), University of Pittsburgh, 3941
O’Hara Street, Pittsburgh, PA 15260, USA
12 Kavli Institute for Particle Astrophysics and Cosmology, Stanford
University, Menlo Park, CA 94305, USA
13 SLAC National Accelerator Laboratory, Menlo Park, CA 94305, USA
14 Departamento de Física, Universidad de los Andes, Cra. 1 No. 18A-10,
Edificio Ip, CP 111711, Bogotá, Colombia
15 Observatorio Astronómico, Universidad de los Andes, Cra. 1 No. 18A-10,
Edificio H, CP 111711 Bogotá, Colombia
16 Institut d’Estudis Espacials de Catalunya (IEEC), 08034 Barcelona, Spain
17 Institute of Cosmology & Gravitation, University of Portsmouth, Dennis
Sciama Building, Portsmouth, PO1 3FX, UK
18 Institute of Space Sciences, ICE-CSIC, Campus UAB, Carrer de Can Magrans
s/n, 08913 Bellaterra, Barcelona, Spain
19 Center for Cosmology and AstroParticle Physics, The Ohio State University,
191 West Woodruff Avenue, Columbus, OH 43210, USA
20 Department of Physics, The Ohio State University, 191 West Woodruff Avenue,
Columbus, OH 43210, USA
21 The Ohio State University, Columbus, 43210 OH, USA
22 School of Mathematics and Physics, University of Queensland, 4072,
Australia
23 Department of Physics, Southern Methodist University, 3215 Daniel Avenue,
Dallas, TX 75275, USA
24 Sorbonne Université, CNRS/IN2P3, Laboratoire de Physique Nucléaire et de
Hautes Energies (LPNHE), FR-75005 Paris, France
25 Departament de Física, Serra Húnter, Universitat Autònoma de Barcelona,
08193 Bellaterra (Barcelona), Spain
26 Institució Catalana de Recerca i Estudis Avançats, Passeig de Lluís
Companys, 23, 08010 Barcelona, Spain
27 Department of Physics and Astronomy, University of Sussex, Brighton BN1
9QH, U.K
28 Department of Physics & Astronomy, University of Wyoming, 1000 E.
University, Dept. 3905, Laramie, WY 82071, USA
29 National Astronomical Observatories, Chinese Academy of Sciences, A20 Datun
Rd., Chaoyang District, Beijing, 100012, P.R. China
30 Departamento de Física, Universidad de Guanajuato - DCI, C.P. 37150, Leon,
Guanajuato, México
Instituto Avanzado de Cosmología A. C., San Marcos 11 - Atenas 202.
31 Magdalena Contreras, 10720. Ciudad de México, México 32 IRFU, CEA,
Université Paris-Saclay, F-91191 Gif-sur-Yvette, France
33 Space Sciences Laboratory, University of California, Berkeley, 7 Gauss Way,
Berkeley, CA 94720, USA
34 University of California, Berkeley, 110 Sproul Hall #5800 Berkeley, CA
94720, USA
35 Department of Physics, Kansas State University, 116 Cardwell Hall,
Manhattan, KS 66506, USA
36 Department of Physics and Astronomy, Sejong University, Seoul, 143-747,
Korea
37 CIEMAT, Avenida Complutense 40, E-28040 Madrid, Spain
38 Space Telescope Science Institute, 3700 San Martin Drive, Baltimore, MD
21218, USA
39 Department of Physics, University of Michigan, Ann Arbor, MI 48109, USA
40 University of Michigan, Ann Arbor, MI 48109, USA

##  Appendix B Assumption of Linearity

The primary assumption of the statistical model described in Section 3 is that
the deviations in the 2pcf and 3pcf from the Gaussian case are linear in
$f_{\mathrm{NL}}$. However, the theoretical 2pcf in the presence of PNG
contains a term which scales as $f_{\mathrm{NL}}^{2}$. Similarly, the 3pcf
with PNG contains a term which scales as $f_{\mathrm{NL}}^{3}$. However, the
smallness of $\alpha^{-1}$ allows us to claim that for the scales considered
in this analysis, the linear term is dominant.

We verify this assumption using an additional set of larger volume FastPM
simulations with PNG (Chaussidon, 2024). Using the same FastPM method as the
simulations described in Section 4.1, we generate four simulations in $(5.52\
h^{-1}\mathrm{Gpc})^{3}$ boxes evolved to $z=1.5$. These larger simulations
are run with $f_{\mathrm{NL}}=-25,0,12,25$. We again create catalogs by
selecting halos with $M>1.15\times 10^{13}\ h^{-1}M_{\odot}$.

For each one, we measure the 2pcf and diagonal 3pcf using ConKer at the same
scales as the ones used in our model. Fig. 14 shows the magnitude of the 2pcf
in selected bins plotted as a function of $f_{\mathrm{NL}}$. For each selected
bin, we perform both a linear and quadratic fit. At the lowest scales, the
linear assumption is extremely well motivated, and we observe almost no
difference between the linear and quadratic fits. As the scale increases we
begin to see deviations from linearity. The decision to truncate our model at
$s=380\ h^{-1}\mathrm{Mpc}$ is motivated by this breakdown of linearity. A
more conservative upper limit at a lower scale would better preserve the
linearity of the 2pcf, but increase the uncertainties on $f_{\mathrm{NL}}$. In
order to extend the model beyond these scales, we would likely need to model
both the linear and quadratic $f_{\mathrm{NL}}$ terms, which would
significantly complicate the procedure.

We perform a similar test for the diagonal 3pcf, which is shown in Fig. 15.
Here, we also include a cubic fit as the 3pcf contains a term which scales as
$f_{\mathrm{NL}}^{3}$. For the 3pcf, we truncate the model at $s=200\
h^{-1}\mathrm{Mpc}$. Fig. 15 shows that below this scale, the assumption of
linearity in $f_{\mathrm{NL}}$ is justified. For the 3pcf when $s_{1}$ is not
equal to $s_{2}$, we exclude any bins where either $s_{1}$ or $s_{2}$ is
larger than $s=200\ h^{-1}\mathrm{Mpc}$.

These tests, performed using the FastPM-L5.52 simulations, motivate the choice
of scales in our statistical model. For the 2pcf when $50\
h^{-1}\mathrm{Mpc}<s<380\ h^{-1}\mathrm{Mpc}$ and for the 3pcf when $50\
h^{-1}\mathrm{Mpc}<s_{1},s_{2}<200\ h^{-1}\mathrm{Mpc}$, we are confident that
changes in the clustering signal due to PNG are linear in $f_{\mathrm{NL}}$.

Figure 14: The 2pcf of the FastPM-L5.52 simulations in selected bins, shown as
a function of $f_{\mathrm{NL}}$ (blue markers). The grey dashed line denotes a
linear fit while the red dot-dash line denotes a quadratic fit. Figure 15: The
diagonal 3pcf of the FastPM-L5.52 simulations in selected bins, shown as a
function of $f_{\mathrm{NL}}$ (blue markers). The grey dashed line denotes a
linear fit, the red dot-dash line denotes a quadratic fit, and the green
dotted line denotes a cubic fit.

##  Appendix C Non-zero PNG validation

In Section 4.6, we demonstrate the ability of our method to correctly measure
$f_{\mathrm{NL}}$ in independent simulations without PNG
($f_{\mathrm{NL}}=0$). To validate our method, we must also verify that it can
correctly measure $f_{\mathrm{NL}}$ in the presence of PNG. The AbacusSummit
simulations are unsuitable for this test, as they have only been generated
with Gaussian initial conditions. However, we can use the FastPM-L5.52
simulations described in Appendix B for such a test.

We begin by constructing the statistical model using the FastPM-L3 simulations
to interpolate the values of $A_{n}$. For the test case, we use each of the
FastPM-L5.52 simulations, where $f_{\mathrm{NL}}=-25,0,12,25$. The fiducial
$n$pcfs, or $\xi_{n}^{\mathrm{fid}}$, are derived from the $f_{\mathrm{NL}}=0$
case of these larger simulations. The covariance matrix is estimated by
scaling the covariance of the FastPM-L3 simulations with $f_{\mathrm{NL}}=0$
according to the difference in volume. This is required because the
computational burden of generating sufficiently many FastPM-L5.52 simulations
to evaluate the covariance is too great.

In Fig. 16, we show the marginalized distributions of the POIs for each case.
For the zero and positive $f_{\mathrm{NL}}$ cases, the true value falls within
the $1\sigma$ contours. For the case where $f_{\mathrm{NL}}=-25$, we slightly
underestimate the value at between $1\sigma$ and $2\sigma$. We do, however,
expect a slight deviation between the true and measured values, due in part to
the halo finding procedure employed in the simulations. The halos in these
FastPM simulations are determined from a Friends-of-Friends algorithm (Roy et
al., 2014), and Chaussidon (2024) showed that it reduces the amplitude of the
$f_{\mathrm{NL}}$ signal up to 30%, which was also noted in Biagetti et al.
(2017). The values of $f_{\mathrm{NL}}$ that we measure are compatible with
the ones measured in Chaussidon (2024). We may also expect deviation between
the true and measures values due to our linear parameterization. From Fig. 14,
we see that points corresponding to positive values of $f_{\mathrm{NL}}$ are
still in the linear regime, while points corresponding to negative
$f_{\mathrm{NL}}$ deviate from linearity. This is because the positive 0th and
negative 1st order terms cancel, making the 2nd order term more significant.
Since the parameterization was derived from mocks with positive
$f_{\mathrm{NL}}$, it is not surprising that the results of the fit show
better agreement for positive $f_{\mathrm{NL}}$, while deviating more
significantly from the true value for the negative $f_{\mathrm{NL}}$ cases.

This test demonstrates the ability of our model to accurately constrain
$f_{\mathrm{NL}}$ for an independent set of simulations at a different
effective redshift than those used to generate the values of $A_{n}$. From
this, we remain confident in the ability to apply these techniques to the
upcoming DESI data releases.

Figure 16: Marginalized distributions of the POIs for each of the FastPM-L5.52
simulations (green, blue, orange, and red corresponding to
$f_{\mathrm{NL}}=-25,0,12,25$ respectively). The dark and light regions
represent $1\sigma$ and $2\sigma$ contours, respectively. The most probable
values and $1\sigma$ CL are labeled for each parameter above the respective
panel. The dashed lines indicate the true values, which are also given as
$\tilde{f}_{\mathrm{NL}}$ in the legend, where the tilde indicates truth.