#  Parity-violating scalar trispectrum from a rolling axion during inflation

Tomohiro Fujita  Tomoaki Murata  Ippei Obata  Maresuke Shiraishi

###### Abstract

We study a mechanism of generating the trispectrum (4-point correlation) of
curvature perturbation through the dynamics of a spectator axion field and
U(1) gauge field during inflation. Owing to the Chern-Simons coupling, only
one helicity mode of gauge field experiences a tachyonic instability and
sources scalar perturbations. Sourced curvature perturbation exhibits parity-
violating nature which can be tested through its trispectrum. We numerically
compute parity-even and parity-odd component of the sourced trispectrum. It is
found that the ratio of parity-odd to parity-even mode can reach
$\mathcal{O}(10\%)$ in an exact equilateral momentum configuration. We also
investigate a quasi-equilateral shape where only one of the momenta is
slightly longer than the other three, and find that the parity-odd mode can
reach, and more interestingly, surpass the parity-even one. This may help us
to interpret a large parity-odd trispectrum signal extracted from BOSS galaxy-
clustering data.

RUP-23-20

##  1 Introduction

Is the parity symmetry of our universe broken? While the standard model of
elementary particles (SM) violates the parity symmetry, its effect would be
negligible for physics on large scales in our universe. Also, general
relativity (GR), which governs the time evolution of our universe, respects
the parity symmetry. Hence, exploring the cosmological parity-violating
phenomena may help us to search for new physics beyond SM or GR, such as
identifying the nature of the dark sector described in the standard $\Lambda$
cold dark matter ($\Lambda$CDM) scenario.

Temperature and E/B-mode polarization of the cosmic microwave background (CMB)
are powerful observables of cosmological parity violation. A cross correlation
between the temperature/E-mode polarization field and the B-mode polarization
one, i.e., a $TB$/$EB$ correlation, that is parity-odd quantity, is a
frequently-used diagonostic tool [1]. Refs. [2, 3, 4, 5, 6] have found
$>3\sigma$ signals in the $EB$ correlation from the Planck and WMAP data.
These are likely to be produced by rotating intrinsic E-mode polarization at
late stages of our universe (a.k.a the cosmic birefringence effect); hence, a
variety of late-time phenomena relevant to new particles beyond SM, such as
axions, have been examined [7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19,
20, 21, 22, 23, 24]. If the parity symmetry is broken by primordial
gravitational waves, it would be also imprinted in the $TB$ and $EB$
correlations. Although such a primordial tensor signal has not been detected
with $>2\sigma$ level so far in these cross-correlations [25, 26, 27, 28] as
well as 3-point correlations (bispectra) [29, 30, 31], many theoretical works
have been done [32, 33, 34, 35, 36, 37, 38, 39, 40, 41, 42, 43, 44, 45, 46,
47, 48, 49, 50, 51, 52, 53, 54, 55, 56, 57, 58, 59, 60, 61, 62, 63, 64, 65,
66, 67, 68, 69, 70, 71]. This is because the production of the parity-
violating gravitational waves is well motivated by axions. Axion-gauge
dynamics could trigger the particle production of one helicity mode of gauge
field and source the cosmological perturbations (also, see [72, 73] for recent
reviews).

Differently from the tensor sector, the parity-odd information in the
primordial scalar sector is inaccessible in the 2-point and 3-point statistics
due to the statistical isotropy of our universe. Hence, a 4-point correlation
(trispectrum) becomes the lowest-order diagonostic tool for the scalar-sector
parity violation [74]. Motivated by this fact, Refs. [75, 76, 77, 78, 79, 80]
have measured the trispectrum of galaxy number density fields with the BOSS
galaxy-clustering data and found $>3\sigma$ parity-odd signals. Measuring the
CMB temperature and E-mode polarization trispecta using the Planck data has
indicated weaker significance level [81, 82], while it is still worth
discussing the relation to new physics, and several generation mechanisms have
already been proposed [74, 83, 84, 79, 85, 80, 86, 87].

In anticipation of a physical interpretation of the measured large parity-odd
signal in the scalar trispectrum, this paper focuses on the axion-gauge
dynamics during inflation, and, for the first time, examines the trispectrum
generation in a model of spectator axion-$U(1)$ gauge field, called the
rolling axion model [41]. In this model, an axion is not an inflaton but a
spectator field and the axion velocity dynamically evolves in time. At one-
loop level, the sourced gauge field provides the scalar-mode correlation with
a unique scale dependence, and generates a sizable amount of non-Gaussianity
at around equilateral wavenumber configurations. Via massive analytical and
numerical analyses, we confirm that the induced scalar trispectrum has not
only parity-even but also parity-odd signals, which corresponds to real and
imaginary numbers, respectively. We also find that the ratio of parity-odd to
parity-even mode is ${\cal O}(10\%)$ at the exact equilateral limit
($k_{1}=k_{2}=k_{3}=k_{4}$), and more interestingly, can reach/exceed unity
for quasi-equilateral configurations ($k_{1}=k_{2}=k_{3}\sim k_{4}$) depending
on orientations of wavevectors.

Previous studies [74, 85, 87] have also computed the scalar trispectrum in
other axion models where an inflaton is identified with an axion, and the
axion velocity is constant in time.111For a tensor trispectrum analysis, see
Ref. [88]. The model of Ref. [74] (originally [40]), predicting a sizable
signal rather for collapsed configurations (e.g., $|\bm{k}_{1}+\bm{k}_{2}|\ll
k_{1},k_{2},k_{3},k_{4}$) due to a $f(\phi)(F^{2}+F\tilde{F})$ coupling, has
already been tested with the Planck and BOSS data, and no significant evidence
(with a maximal deviation of $2.0\sigma$) has been found [78, 81, 82]. The
model of Ref. [87] includes a nonvanishing mass term of the gauge field,
inducing different trispectrum shapes. Regarding the model of Ref. [85], like
our case, the trispectrum peaks for the equilateral configurations. However,
due to the difference of the generation mechanism, the parity-odd signal in
our model can become larger than that in Ref. [85]; thus, our model may be
more useful on a physical interpretation of the large parity-odd trispectrum
signal observed in the BOSS data.

This paper is organized as follows. In section 2, we briefly describe a model
of rolling axion-$U(1)$ model. In section 3, we present a perturbation
dynamics of sourced scalar mode by the enhanced gauge field. In section 4, we
analytically formulate the scalar-mode trispectrum. In section 5, we
numerically evaluate the parity-odd and parity-even component in the
trispectrum. We will summarize our result in section 6. Throughout this paper,
we set a natural unit $\hbar=c=1$.

##  2 Model

We consider the following axion-$U(1)$ spectator model (known as a rolling
axion model [41]) where a spectator axionic field $\sigma$ couples to $U(1)$
gauge field $A_{\mu}$ during inflation:

$S=\int d^{4}x\sqrt{-g}\left[\dfrac{1}{2}M_{\rm
Pl}^{2}R-\dfrac{1}{2}(\partial\varphi)^{2}-U(\varphi)-\dfrac{1}{2}(\partial\sigma)^{2}-V(\sigma)-\dfrac{1}{4}F_{\mu\nu}F^{\mu\nu}+\dfrac{\lambda}{4}\dfrac{\sigma}{f}F_{\mu\nu}\tilde{F}^{\mu\nu}\right]\
,$ (2.1)

where $R$ is the Ricci scalar, $\varphi$ is an inflaton, $U(\varphi)$ is its
potential, $V(\sigma)$ is a potential of the axion,
$F_{\mu\nu}=\partial_{\mu}A_{\nu}-\partial_{\nu}A_{\mu}$ is the field strength
of $U(1)$ gauge field,
$\tilde{F}^{\mu\nu}=\epsilon^{\mu\nu\rho\sigma}F_{\rho\sigma}/(2\sqrt{-g})$ is
its Hodge dual, $\epsilon^{\mu\nu\rho\lambda}$ is the Levi-Cività symbol with
$\epsilon^{0123}=1$, $M_{\rm Pl}\equiv 1/\sqrt{8\pi G}$ is the reduced Planck
mass, and $\lambda,f$ are constant parameters. We consider the flat FLRW
Universe, ${\rm d}s^{2}=-{\rm d}t^{2}+a^{2}(t){\rm d}{\bm{x}}^{2}$, with
$a(t)$ being the scale factor. In what follows, we assume that the unspecified
inflaton field dominates the inflationary universe (i.e. $V(\sigma)\ll
U(\varphi)$) and the Hubble parameter $H(t)\equiv\dot{a}/a$ is well
approximated by a constant (i.e. $\dot{\varphi}^{2}\ll U(\varphi)$), where dot
denotes the time derivative, $\dot{X}\equiv\partial_{t}X$. Our analysis in
this and the next section mostly follows that of Ref. [41].

Regarding the spectator potential, we adopt the cosine potential form:

$V(\sigma)=\Lambda^{4}\left[1-\cos\left(\dfrac{\sigma}{f}\right)\right]\ .$
(2.2)

The equation of motion (EoM) for the background spectator field $\sigma(t)$ is
given by

$\ddot{\sigma}+3H\dot{\sigma}+\partial_{\sigma}V=0\ ,$ (2.3)

where we neglect the backreaction term from the gauge field, $\lambda\langle
E_{i}B_{i}\rangle/f$. Further neglecting the $\ddot{\sigma}$ term and solving
the slow-roll equation $3H\dot{\sigma}=-\partial_{\sigma}V$, we obtain the
solution as

$\sigma(t)=2f\text{Arccot}\left(e^{\delta H(t-t_{*})}\right)\
,\qquad\delta\equiv\dfrac{\Lambda^{4}}{3H^{2}f^{2}}\ ,$ (2.4)

where a dimensionless parameter $\delta$ characterizes the width of cosine
potential and $t_{*}$ is a time when the axion passes through the inflection
point of the potential. Then, the slow-roll condition

$\frac{\ddot{\sigma}}{3H\dot{\sigma}}=-\frac{\delta}{3}\tanh{\left(\delta
H(t-t_{*})\right)}\ll 1,$ (2.5)

is valid when $\delta\ll 3$ is satisfied.

For later convenience, we define the dimensionless parameter

$\xi(t)\equiv-\lambda\dfrac{\dot{\sigma}(t)}{2fH}>0\quad(\dot{\sigma}<0)$
(2.6)

which describes the speed of spectator field and we assume it is positive in
this paper. The slow-roll solution for $\xi$ is given by

$\xi(t)=\dfrac{\xi_{*}}{\cosh(\delta H(t-t_{*}))}\
,\qquad\xi_{*}\equiv\dfrac{\lambda\delta}{2}\ ,$ (2.7)

where $\xi_{*}=\xi(t_{*})$ is the maximum value of $\xi$.

In analyzing the gauge field perturbation, we switch the time variable from
the cosmic time $t$ to the conformal time $\tau$. We fix the gauge by choosing
the Coulomb gauge $\partial_{i}A_{i}=0$ and the temporal gauge $A_{0}=0$. To
describe the production of quanta of the gauge field, we promote the classical
field $A_{i}$ to an operator

$\displaystyle\hat{A}_{i}(\tau,\bm{x})$
$\displaystyle=\sum_{\lambda=\pm}\int\dfrac{{\rm
d}^{3}k}{(2\pi)^{3}}\hat{A}^{\lambda}_{\bm{k}}(\tau)e^{\lambda}_{i}(\hat{\bm{k}})e^{i\bm{k}\cdot\bm{x}}\
,$ (2.8) $\displaystyle\hat{A}^{\lambda}_{\bm{k}}(\tau)$
$\displaystyle=\hat{a}^{\lambda}_{\bm{k}}A^{\lambda}_{k}(\tau)+\hat{a}^{\lambda\dagger}_{-\bm{k}}A^{\lambda*}_{k}(\tau),\quad(\lambda=\pm)\
,$ (2.9)

where $\hat{a}^{\lambda}_{\bm{k}},\hat{a}^{\lambda\dagger}_{\bm{k}}$ are the
creation and annihilation operators obeying the commutation relation,
$[\hat{a}^{\lambda}_{\bm{k}},\
\hat{a}^{\lambda^{\prime}\dagger}_{-\bm{k}^{\prime}}]=\delta^{\lambda\lambda^{\prime}}(2\pi)^{3}\delta^{3}(\bm{k}+\bm{k}^{\prime})$.
We employ the circular polarization vectors

$\bm{e}^{\pm}(\hat{\bm{k}})=\dfrac{1}{\sqrt{2}}(\cos\theta_{\hat{\bm{k}}}\cos\phi_{\hat{\bm{k}}}\mp
i\sin\phi_{\hat{\bm{k}}},\ \cos\theta_{\hat{\bm{k}}}\sin\phi_{\hat{\bm{k}}}\pm
i\cos\phi_{\hat{\bm{k}}},-\sin\theta_{\hat{\bm{k}}})\ ,$ (2.10)

in terms of
$\hat{\bm{k}}=(\sin\theta_{\hat{\bm{k}}}\cos\phi_{\hat{\bm{k}}},\sin\theta_{\hat{\bm{k}}}\sin\phi_{\hat{\bm{k}}},\cos\theta_{\hat{\bm{k}}})$,
and it satisfies the orthonormal condition
$e^{\lambda}_{i}(\hat{\bm{k}})e^{\lambda^{\prime}*}_{i}(\hat{\bm{k}})=\delta^{\lambda\lambda^{\prime}}$,
the transverse property $k_{i}e^{\pm}_{i}(\hat{\bm{k}})=0$ and the
eigenvectors of the curl $i\epsilon_{ijk}k_{j}e^{\pm}_{k}(\hat{\bm{k}})=\pm
ke^{\pm}_{i}(\hat{\bm{k}})$, where $k=|\bm{k}|$, $\hat{\bm{k}}=\bm{k}/k$ and
$\epsilon_{ijk}\equiv\epsilon_{0ijk}$.

The EoM for the mode function of the gauge field reads

$\left[\partial_{\tau}^{2}+k^{2}\pm\dfrac{2k\xi(\tau)}{\tau}\right]A^{\pm}_{k}(\tau)=0\
,\qquad\xi(\tau)=\dfrac{2\xi_{*}}{\left(\tfrac{\tau_{*}}{\tau}\right)^{\delta}+\left(\tfrac{\tau}{\tau_{*}}\right)^{\delta}}\
,$ (2.11)

where we use $a=-1/(H\tau)$ and $a_{*}=-1/(H\tau_{*})$. While the minus mode
$A^{-}_{k}$ is not amplified, we can see that the plus mode $A^{+}_{k}$
experiences a tachyonic instability for $-k\tau<2\xi(\tau)$. We are interested
in the parameter region $\xi(\tau)\gtrsim 1$, so that $A^{+}_{k}$ starts to
grow slightly before the horizon-crossing, $-k\tau=1$.

Although the solution of (2.11) cannot be written in a closed form, we can use
the WKB approximation. Then, for $\tau>\bar{\tau}$ obeying
$-k\bar{\tau}=2\xi(\bar{\tau})$, we obtain the following approximate solution:

$\displaystyle A^{+}_{k}(\tau)$
$\displaystyle\simeq\dfrac{N(\xi_{*},-k\tau_{*},\delta)}{\sqrt{2k}}\left(\dfrac{-k\tau}{2\xi(\tau)}\right)^{1/4}$
$\displaystyle\times\exp\left[-\dfrac{4\xi_{*}^{1/2}}{1+\delta}\left(\dfrac{\tau}{\tau_{*}}\right)^{\delta/2}(-k\tau)^{1/2}{{}_{2}}F_{1}\left(\tfrac{1}{2},\
\tfrac{1+\delta}{4\delta};\ \tfrac{5\delta+1}{4\delta};\
-\left(\tfrac{\tau}{\tau_{*}}\right)^{2\delta}\right)\right]\,,$
$\displaystyle\equiv\dfrac{1}{\sqrt{2k}}\mathcal{A}(\xi,-k\tau)\ ,$ (2.12)
$\displaystyle\partial_{\tau}A^{+}_{k}(\tau)$
$\displaystyle\simeq\sqrt{\dfrac{2k\xi(\tau)}{-\tau}}A^{+}_{k}(\tau)\ ,$
(2.13)

where the normalization factor $N$ is determined by matching $A^{+}_{k}$ at
late time to the full numerical solution of (2.11) and ${{}_{2}}F_{1}$ is a
hypergeometric function. In the deep super-horizon regime $-k\tau\rightarrow
0$, the energy density of gauge fields dilutes due to the expansion of the
universe. Thus, the gauge field makes the most contribution to the physical
enhancement of fluctuations at around the horizon-crossing. A similar
approximated expression was used in the previous work [41]. Here we find a
more accurate expression with the hypergeometric function, which is derived in
Appendix A.

##  3 Scalar mode production

We will estimate the four-point correlation function of the curvature
perturbation on a flat hypersurfaces $\zeta=-H\delta\varphi/\dot{\varphi}$
induced by the gauge field perturbation, which is depicted by a one-loop
diagram in Figure 1. To this end, we first calculate the fluctuation of the
spectator axion $\delta\sigma$ in this section. $\delta\sigma$ is produced by
the intrinsic inhomogeneity of the Chern-Simons gauge interaction
$\delta_{EB}\equiv\left(E_{i}B_{i}-\langle
E_{i}B_{i}\rangle\right)|_{\delta\sigma=0}$, where
$E_{i}=-\partial_{\tau}A_{i}/a^{2},\
B_{i}=\epsilon_{ijk}\partial_{j}A_{k}/a^{2}$. The EoM for $\delta\sigma$ is
given by

$\displaystyle\ddot{\delta\sigma}+3H\dot{\delta\sigma}-a^{-2}\nabla^{2}\delta\sigma\simeq-\dfrac{\lambda}{f}\delta_{EB}\
.$ (3.1)

Note that we do not include the mass term of $\delta\sigma$ and other higher
order interactions for simplicity. The effective mass squared of
$\delta\sigma$ is $V^{\prime\prime}(\sigma)=3\delta H^{2}\cos{(\sigma/f)}$, it
exactly vanishes at the inflection point, and it remains small in other
regions under the the slow-roll condition $\delta\ll 3$.

Since $\langle E_{i}B_{i}\rangle(t)$ contributes only to zero-mode, the
Fourier transform of the above EoM reads

$\left[\dfrac{\partial^{2}}{\partial\tau^{2}}+k^{2}-\dfrac{2}{\tau^{2}}\right](a\delta\hat{\sigma}_{\bm{k}})=\int\dfrac{{\rm
d}^{3}p}{(2\pi)^{3}}\hat{\mathcal{S}}(\tau,\bm{p},\ \bm{k}-\bm{p})\ ,$ (3.2)

where

$\hat{\mathcal{S}}(\tau,\ \bm{p},\
\bm{k}-\bm{p})\equiv-\lambda\dfrac{a^{3}}{2f}e^{+}_{i}(\hat{\bm{p}})e^{+}_{i}(\savestack{\tmpbox}{\stretchto{\scaleto{\scalerel*[\widthof{\bm{k}-\bm{p}}]{\kern-0.6pt\bigwedge\kern-0.6pt}{\rule[-642.0pt]{4.30554pt}{642.0pt}}}{}}{0.5ex}}\stackon[1pt]{\bm{k}-\bm{p}}{\tmpbox})\left(\hat{E}^{+}_{\bm{p}}\hat{B}^{+}_{\bm{k}-\bm{p}}+\hat{E}^{+}_{\bm{k}-\bm{p}}\hat{B}^{+}_{\bm{p}}\right),$
(3.3)

where we define
$\hat{E}^{+}_{\bm{p}}=-\partial_{\tau}\hat{A}^{+}_{\bm{p}}/a^{2},$ and
$\hat{B}^{+}_{\bm{p}}=|\bm{p}|\hat{A}^{+}_{\bm{p}}/a^{2}$, and we have
symmetrized the operators of the source term with respect to the momenta
$\bm{p}$ and $\bm{k}-\bm{p}$ to make the correlation functions independent
from the order of operators. $\delta\sigma_{k}$ has two solutions
$\delta\sigma_{k}=\delta\sigma^{\rm(V)}_{k}+\delta\sigma^{\rm(S)}_{k}$. The
homogeneous solution $\delta\sigma^{\rm(V)}_{k}$ is a normal vacuum
fluctuation and its expression for the Bunch-Davies initial condition is given
by the first kind of Hankel function $H_{\nu}^{(1)}=J_{\nu}+iY_{\nu}$:

$a\delta\sigma_{k}^{\rm(V)}=-\dfrac{1}{\sqrt{2k}}\sqrt{\dfrac{-\pi
k\tau}{2}}H^{(1)}_{3/2}(-k\tau)\equiv u(-k\tau)\ .$ (3.4)

In the super-horizon limit, this solution yields the well-known amplitude
$\langle(\delta\sigma_{k}^{\rm(V)})^{2}\rangle^{1/2}=H/(2\pi)$. In contrast,
the inhomogeneous solution $\delta\sigma_{k}^{\rm(S)}$ sourced by the gauge
field can be found as

$\displaystyle\delta\hat{\sigma}_{k}^{\rm(S)}$
$\displaystyle\simeq\dfrac{1}{a}\int_{-\infty}^{\infty}d\tau^{\prime}G_{k}(\tau,\tau^{\prime})\int\dfrac{{\rm
d}^{3}p}{(2\pi)^{3}}\hat{\mathcal{S}}(\tau^{\prime},\ \bm{p},\ \bm{k}-\bm{p})\
,$ (3.5)

where

$\displaystyle G_{k}(\tau,\tau^{\prime})$
$\displaystyle=i\Theta(\tau-\tau^{\prime})\left[u(-k\tau)u^{*}(-k\tau^{\prime})-u^{*}(-k\tau)u(-k\tau^{\prime})\right],$
$\displaystyle=\Theta(\tau-\tau^{\prime})\dfrac{\pi}{2}\sqrt{\tau\tau^{\prime}}\left[J_{3/2}(-k\tau)Y_{3/2}(-k\tau^{\prime})-Y_{3/2}(-k\tau)J_{3/2}(-k\tau^{\prime})\right],$
(3.6)

is a retarded Green function and in the super-horizon limit of
$\tau\rightarrow 0$ it reads

$G_{k}(\tau\rightarrow
0,\tau^{\prime})=\Theta(\tau-\tau^{\prime})\sqrt{\dfrac{\pi}{2}}\dfrac{\sqrt{-k\tau^{\prime}}}{-k^{2}\tau}J_{3/2}(-k\tau^{\prime})\
.$ (3.7)

$\zeta$$\zeta$$\zeta$$\zeta$$\sigma$$\sigma$$\sigma$$\sigma$$A^{+}$ Figure 1:
Diagram of the four-point function that we compute. The external solid line
represents the curvature perturbation $\zeta$. The internal dashed line shows
the axion $\sigma$. The intermediate wiggly line represents the one helicity
of gauge field $A^{+}$.

With the aid of slow-roll parameters of the background inflaton and the axion,

$\epsilon_{\varphi}\equiv\dfrac{\dot{\varphi}^{2}}{2M_{\rm Pl}^{2}H^{2}}\
,\qquad\epsilon_{\sigma}\equiv\dfrac{\dot{\sigma}^{2}}{2M_{\rm Pl}^{2}H^{2}}\
,$ (3.8)

the EoM for the fluctuation of the inflaton can be written as

$\left[\dfrac{\partial^{2}}{\partial\tau^{2}}+k^{2}-\dfrac{2}{\tau^{2}}\right]\left(a\delta\hat{\varphi}_{\bm{k}}^{(S)}\right)\simeq\dfrac{6}{\tau^{2}}\sqrt{\epsilon_{\varphi}\epsilon_{\sigma}}\left(a\delta\hat{\sigma}_{\bm{k}}^{(S)}\right)\
.$ (3.9)

The source term comes from the gravitational coupling between the inflaton and
the axion [41]. Here we ignore the mass of the inflaton and its retarded Green
function is identical to that of the axion, Eq. (3.6). Then, we obtain the
inhomogeneous solution of the above equation as

$\delta\hat{\varphi}_{\bm{k}}^{\rm(S)}=\dfrac{6\sqrt{\epsilon_{\varphi}}}{a}\int_{-\infty}^{\infty}d\tau^{\prime}\dfrac{\sqrt{\epsilon_{\sigma}(\tau^{\prime})}}{\tau{{}^{\prime}}^{2}}G_{k}(\tau,\tau^{\prime})\int_{-\infty}^{\infty}d\tau^{\prime\prime}G_{k}(\tau^{\prime},\tau^{\prime\prime})\int\dfrac{{\rm
d}^{3}p}{(2\pi)^{3}}\hat{\mathcal{S}}(\tau^{\prime\prime},\ \bm{p},\
\bm{k}-\bm{p})\ ,$ (3.10)

where we have assumed $\epsilon_{\varphi}=\text{const}$. for simplicity. The
curvature perturbation can be decomposed into two contributions
$\hat{\zeta}_{\bm{k}}=\zeta_{\bm{k}}^{(V)}+\hat{\zeta}_{\bm{k}}^{(S)}$, where
the vacuum contribution and the sourced contribution is given by

$\displaystyle\zeta_{\bm{k}}^{(V)}=-\frac{H}{\dot{\varphi}}\delta\varphi_{\bm{k}}^{(V)},\qquad\hat{\zeta}_{\bm{k}}^{(S)}=-\frac{H}{\dot{\varphi}}\delta\hat{\varphi}_{\bm{k}}^{(S)},$
(3.11)

respectively.

In passing, we quickly discuss the power spectrum of curvature perturbation
$\mathcal{P}_{\zeta}$, which is defined by

$\langle\hat{\zeta}_{\bm{k}}\hat{\zeta}_{\bm{k}^{\prime}}\rangle=(2\pi)^{3}\delta(\bm{k}+\bm{k}^{\prime})\dfrac{2\pi^{2}}{k^{3}}\mathcal{P}_{\zeta}\,.$
(3.12)

The contributions comes from the vacuum fluctuation $\delta\varphi_{k}^{(V)}$,
whose mode function is identical to Eq. (3.4) in the slow-roll regime, and the
sourced one is obtained in Eq. (3.10). Combining them and using Eq. (3.11),
one finds

$\mathcal{P}_{\zeta}=\mathcal{P}_{\zeta,v}\left[1+\mathcal{P}_{\zeta,v}f_{2}(\xi)e^{4\pi\xi}\right]\
,$ (3.13)

where $\mathcal{P}_{\zeta,v}\equiv H^{2}/(8\pi^{2}M_{\rm
Pl}^{2}\epsilon_{\varphi})$ is the dimensionless power spectrum of the vacuum
curvature perturbation. The function $f_{2}(\xi)$ was found numerically. For
large and constant $\xi$, it is given by [33]

$f_{2}(\xi)\simeq\dfrac{7.5\times 10^{-5}}{\xi^{6}}\ .$ (3.14)

When $\xi$ is dynamical, a more elaborated fitting function is known [41]. One
can also compute the bispectrum of the curvature perturbation induced by the
sourced inflaton fluctuation $\delta\hat{\varphi}_{k}^{\rm(S)}$ in the same
way [41]. However, we will skip it and calculate the trispectrum in the next
section.

##  4 Formalism of the trispectrum

In this section, we derive some expressions that will be used for the
evaluation of the trispectrum in the next section. The trispectrum of the
sourced curvature perturbation $T_{\zeta}$ is defined as

$\langle\hat{\zeta}_{\bm{k}_{1}}^{(S)}\hat{\zeta}_{\bm{k}_{2}}^{(S)}\hat{\zeta}_{\bm{k}_{3}}^{(S)}\hat{\zeta}_{\bm{k}_{4}}^{(S)}\rangle\equiv(2\pi)^{3}\delta(\bm{k}_{1}+\bm{k}_{2}+\bm{k}_{3}+\bm{k}_{4})T_{\zeta}(\bm{k}_{1},\
\bm{k}_{2},\ \bm{k}_{3},\ \bm{k}_{4})\ .$ (4.1)

In particular, we will calculate the real and imaginary parts of the above
trispectrum separately for the following reason [74]. The parity
transformation in real space is given by flipping the sign of the spatial
coordinates as $\bm{x}\to-\bm{x}$. Similarly, the parity transformation in
Fourier space corresponds to flipping the sign of the wave number vector as
$\bm{k}\to-\bm{k}$. Then, we consider the parity transformation of the
trispectrum of curvature perturbation, we obtain

$\displaystyle\langle\hat{\zeta}_{\bm{k}_{1}}^{(S)}\hat{\zeta}_{\bm{k}_{2}}^{(S)}\hat{\zeta}_{\bm{k}_{3}}^{(S)}\hat{\zeta}_{\bm{k}_{4}}^{(S)}\rangle\
\to\
\langle\hat{\zeta}_{-\bm{k}_{1}}^{(S)}\hat{\zeta}_{-\bm{k}_{2}}^{(S)}\hat{\zeta}_{-\bm{k}_{3}}^{(S)}\hat{\zeta}_{-\bm{k}_{4}}^{(S)}\rangle=\langle\hat{\zeta}_{\bm{k}_{1}}^{(S)}\hat{\zeta}_{\bm{k}_{2}}^{(S)}\hat{\zeta}_{\bm{k}_{3}}^{(S)}\hat{\zeta}_{\bm{k}_{4}}^{(S)}\rangle^{*}\
,$ (4.2)

where we use the reality condition of curvature perturbation
$\zeta_{-\bm{k}}^{(S)}=\zeta_{\bm{k}}^{(S)*}$ in the last equality. We can see
that the real part of the trispectrum ${\rm Re}[T_{\zeta}]$ is invariant under
parity transformation, so this is parity even. In contrast, the imaginary part
${\rm Im}[T_{\zeta}]$ changes the sign in the parity transformation.
Therefore, the imaginary part of trispectrum is sensitive to parity violation.

We begin with the four-point correlation function of the sourced inflaton
fluctuation $\delta\hat{\varphi}_{\bm{k}}^{(S)}$. Using Eq. (3.10), we find

$\displaystyle\langle\delta\hat{\varphi}_{\bm{k}_{1}}^{(S)}\delta\hat{\varphi}_{\bm{k}_{2}}^{(S)}\delta\hat{\varphi}_{\bm{k}_{3}}^{(S)}\delta\hat{\varphi}_{\bm{k}_{4}}^{(S)}\rangle$
$\displaystyle=\dfrac{6^{4}\epsilon_{\varphi}^{2}}{a^{4}}\prod_{i=1,2,3,4}\int_{-\infty}^{\infty}d\tau^{\prime}_{i}\dfrac{\sqrt{\epsilon_{\sigma}(\tau^{\prime}_{i})}}{\tau{{}^{\prime}}_{i}^{2}}G_{k_{i}}(\tau,\tau^{\prime}_{i})\int_{-\infty}^{\infty}d\tau^{\prime\prime}_{i}G_{k_{i}}(\tau^{\prime}_{i},\tau^{\prime\prime}_{i})\int\dfrac{d\bm{p}_{i}}{(2\pi)^{3}}$
$\displaystyle\times\langle\hat{\mathcal{S}}(\tau^{\prime\prime}_{1},\
\bm{p}_{1},\ \bm{k}_{1}-\bm{p}_{1})\hat{\mathcal{S}}(\tau^{\prime\prime}_{2},\
\bm{p}_{2},\ \bm{k}_{2}-\bm{p}_{2})\hat{\mathcal{S}}(\tau^{\prime\prime}_{3},\
\bm{p}_{3},\ \bm{k}_{3}-\bm{p}_{3})\hat{\mathcal{S}}(\tau^{\prime\prime}_{4},\
\bm{p}_{4},\ \bm{k}_{4}-\bm{p}_{4})\rangle\ .$ (4.3)

Dropping the disconnected contributions and keeping the connected ones
corresponding to the Feynman diagram shown in Fig. 1, we obtain

$\displaystyle\langle\hat{\mathcal{S}}(\tau^{\prime\prime}_{1},\ \bm{p}_{1},\
\bm{k}_{1}-\bm{p}_{1})\hat{\mathcal{S}}(\tau^{\prime\prime}_{2},\ \bm{p}_{2},\
\bm{k}_{2}-\bm{p}_{2})\hat{\mathcal{S}}(\tau^{\prime\prime}_{3},\ \bm{p}_{3},\
\bm{k}_{3}-\bm{p}_{3})\hat{\mathcal{S}}(\tau^{\prime\prime}_{4},\ \bm{p}_{4},\
\bm{k}_{4}-\bm{p}_{4})\rangle$
$\displaystyle=2^{4}(2\pi)^{12}\left(\delta(\bm{p}_{1}+\bm{p}_{2})\delta(\bm{p}_{3}+\bm{p}_{4})\delta(\bm{k}_{1}-\bm{p}_{1}+\bm{k}_{4}-\bm{p}_{4})\delta(\bm{k}_{2}-\bm{p}_{2}+\bm{k}_{3}-\bm{p}_{3})\right.$
$\displaystyle\left.+\delta(\bm{p}_{1}+\bm{p}_{2})\delta(\bm{p}_{3}+\bm{p}_{4})\delta(\bm{k}_{1}-\bm{p}_{1}+\bm{k}_{3}-\bm{p}_{3})\delta(\bm{k}_{2}-\bm{p}_{2}+\bm{k}_{4}-\bm{p}_{4})\right.$
$\displaystyle\left.+\delta(\bm{p}_{1}+\bm{p}_{3})\delta(\bm{p}_{2}+\bm{p}_{4})\delta(\bm{k}_{1}-\bm{p}_{1}+\bm{k}_{4}-\bm{p}_{4})\delta(\bm{k}_{2}-\bm{p}_{2}+\bm{k}_{3}-\bm{p}_{3})\right)$
$\displaystyle\times\mathcal{S}(\tau^{\prime\prime}_{1},\ \bm{p}_{1},\
\bm{k}_{1}-\bm{p}_{1})\mathcal{S}(\tau^{\prime\prime}_{2},\ \bm{p}_{2},\
\bm{k}_{2}-\bm{p}_{2})\mathcal{S}(\tau^{\prime\prime}_{3},\ \bm{p}_{3},\
\bm{k}_{3}-\bm{p}_{3})\mathcal{S}(\tau^{\prime\prime}_{4},\ \bm{p}_{4},\
\bm{k}_{4}-\bm{p}_{4})\ ,$ (4.4)

where the coefficient $2^{4}$ comes from the symmetric computation of the
quadratic operator $\hat{\mathcal{S}}$. To make the computation simple, in the
whole integral range we adopt the approximation function (2.12) and (2.13)
which receives its support around the horizon crossing. This approximation is
valid because the UV contribution from the deep sub-horizon regime should be
removed and the energy density of electromagnetic field dilutes at the deep
super-horizon regime as discussed in Appendix. A. Plugging Eqs. (2.12) and
(2.13), we obtain

$\displaystyle\mathcal{S}(\tau_{1},\ \bm{p}_{1},\
\bm{k}_{1}-\bm{p}_{1})\mathcal{S}(\tau_{2},\ \bm{p}_{2},\
\bm{k}_{2}-\bm{p}_{2})\mathcal{S}(\tau_{3},\ \bm{p}_{3},\
\bm{k}_{3}-\bm{p}_{3})\mathcal{S}(\tau_{4},\ \bm{p}_{4},\
\bm{k}_{4}-\bm{p}_{4})$
$\displaystyle\simeq\dfrac{\lambda^{4}}{2^{6}f^{4}}\prod_{i=1,2,3,4}\dfrac{e^{+}_{m}(\hat{\bm{p}}_{i})e^{+}_{m}(\savestack{\tmpbox}{\stretchto{\scaleto{\scalerel*[\widthof{\bm{k}_{i}-\bm{p}_{i}}]{\kern-0.6pt\bigwedge\kern-0.6pt}{\rule[-642.0pt]{4.30554pt}{642.0pt}}}{}}{0.5ex}}\stackon[1pt]{\bm{k}_{i}-\bm{p}_{i}}{\tmpbox})}{a(\tau_{i})}\sqrt{\dfrac{\xi(\tau_{i})}{-\tau_{i}}}\mathcal{A}(\xi,-p_{i}\tau_{i})\mathcal{A}(\xi,-|\bm{k}_{i}-\bm{p}_{i}|\tau_{i})\left(\sqrt{|\bm{k}_{i}-\bm{p}_{i}|}+\sqrt{p_{i}}\right)$

We introduce dimensionless time variables $x\equiv-k_{1}\tau,\ x_{i}\equiv-
k_{1}\tau_{i}$ and dimensionless momenta
$\bm{k}_{i}^{*}\equiv\bm{k}_{i}/k_{1},\ \bm{p}_{i}^{*}\equiv\bm{p}_{i}/k_{1}$
normalized by $k_{1}$. Note that these upper asterisks should be distinguished
from the lower asterisks (e.g. $\xi_{*}\equiv\xi(t_{*})$) denoting the time
when the axion passes through the inflection point. Using them, we obtain

$\displaystyle\langle\delta\hat{\varphi}_{\bm{k}_{1}}^{(S)}\delta\hat{\varphi}_{\bm{k}_{2}}^{(S)}\delta\hat{\varphi}_{\bm{k}_{3}}^{(S)}\delta\hat{\varphi}_{\bm{k}_{4}}^{(S)}\rangle$
$\displaystyle\simeq\dfrac{6^{4}\epsilon_{\varphi}^{2}\lambda^{4}H^{8}x^{4}}{4k_{1}^{12}f^{4}}\prod_{i=1,2,3,4}\Bigg{(}\int_{\infty}^{-\infty}dx^{\prime}_{i}\dfrac{\sqrt{\epsilon_{\sigma}(\tau_{i}^{\prime})}}{x_{i}{{}^{\prime}}^{2}}k_{1}G_{k_{i}}(x,x_{i}^{\prime})\int_{\infty}^{-\infty}dx^{\prime\prime}_{i}\sqrt{\dfrac{\xi(\tau_{i}^{\prime\prime})}{x^{\prime\prime}_{i}}}k_{1}x_{i}^{\prime\prime}G_{k_{i}}(x_{i}^{\prime},x_{i}^{\prime\prime})$
$\displaystyle\times\int\dfrac{d\bm{p}^{*}_{i}}{(2\pi)^{3}}e^{+}_{m}(\hat{\bm{p}_{i}})e^{+}_{m}(\savestack{\tmpbox}{\stretchto{\scaleto{\scalerel*[\widthof{\bm{k}_{i}-\bm{p}_{i}}]{\kern-0.6pt\bigwedge\kern-0.6pt}{\rule[-642.0pt]{4.30554pt}{642.0pt}}}{}}{0.5ex}}\stackon[1pt]{\bm{k}_{i}-\bm{p}_{i}}{\tmpbox})\left(\sqrt{p_{i}^{*}}+\sqrt{|\bm{k}_{i}-\bm{p}_{i}|^{*}}\right)\mathcal{A}(\xi,p^{*}_{i}x^{\prime\prime}_{i})\mathcal{A}(\xi,|\bm{k}_{i}-\bm{p}_{i}|^{*}x^{\prime\prime}_{i})$
$\displaystyle\times(2\pi)^{12}\left[\delta(\bm{p}^{*}_{1}+\bm{p}^{*}_{2})\delta(\bm{p}^{*}_{3}+\bm{p}^{*}_{4})\delta(\bm{k}^{*}_{1}-\bm{p}_{1}^{*}+\bm{k}_{4}^{*}-\bm{p}_{4}^{*})\delta(\bm{k}_{2}^{*}-\bm{p}_{2}^{*}+\bm{k}_{3}^{*}-\bm{p}_{3}^{*})\right.$
$\displaystyle\left.+(\text{permutations of}\ \bm{k}_{3}^{*}-\bm{p}_{3}^{*},\
\bm{k}_{4}^{*}-\bm{p}_{4}^{*})+(\text{permutations of}\ \bm{p}_{2}^{*},\
\bm{p}_{3}^{*})\right]\ .$ (4.5)

Since $\delta\hat{\varphi}_{\bm{k}}^{(S)}$ and $\hat{\zeta}_{\bm{k}}^{(S)}$
are related by Eq. (3.11), the trispectrum $T_{\zeta}$ can be calculated from
Eq. (4.5). By computing the integration with delta functions and using
$H/\dot{\varphi}=-1/(M_{\rm Pl}\sqrt{2\epsilon_{\varphi}})$ and
$\sqrt{\epsilon_{\sigma}(\tau_{i}^{\prime})}=\sqrt{2}f\xi(\tau_{i}^{\prime})/(\lambda
M_{\rm Pl})$, we find the asymptotic form of super-horizon limit $x\rightarrow
0$ as

$\displaystyle T_{\zeta}(\bm{k}_{1},\ \bm{k}_{2},\ \bm{k}_{3},\ \bm{k}_{4})$
$\displaystyle\simeq\dfrac{3^{4}}{2^{4}}\pi^{8}r_{v}^{4}\dfrac{\mathcal{P}_{\zeta,v}^{4}}{k_{1}^{9}}\int\dfrac{d\bm{p}^{*}_{1}}{(2\pi)^{3}}\prod_{i=1,2,3,4}e^{+}_{m}(\hat{\bm{p}_{i}})e^{+}_{m}(\savestack{\tmpbox}{\stretchto{\scaleto{\scalerel*[\widthof{\bm{k}_{i}-\bm{p}_{i}}]{\kern-0.6pt\bigwedge\kern-0.6pt}{\rule[-642.0pt]{4.30554pt}{642.0pt}}}{}}{0.5ex}}\stackon[1pt]{\bm{k}_{i}-\bm{p}_{i}}{\tmpbox})\left(\sqrt{p_{i}{{}^{*}}}+\sqrt{|\bm{k}_{i}-\bm{p}_{i}|^{*}}\right)\left(p^{*}_{i}|\bm{k}_{i}-\bm{p}_{i}|^{*}\right)$
$\displaystyle\times\mathcal{T}\left(\xi_{*},x_{*},\delta,\sqrt{p_{i}{{}^{*}}}+\sqrt{|\bm{k}_{i}-\bm{p}_{i}|^{*}}\right)N(\xi_{*},p^{*}_{i}x_{*},\delta)N(\xi_{*},|\bm{k}_{i}-\bm{p}_{i}|^{*}x_{*},\delta)$
$\displaystyle+(\text{permutations of}\ \bm{k}_{3}^{*}-\bm{p}_{3}^{*},\
\bm{k}_{4}^{*}-\bm{p}_{4}^{*})+(\text{permutations of}\ \bm{p}_{2}^{*},\
\bm{p}_{3}^{*})\ .$ (4.6)

where $r_{v}=16\epsilon_{\varphi}$ is a vacuum tensor-to-scalar ratio and in
the first line the momenta are fixed as $\bm{p}_{2}=-\bm{p}_{1},\
\bm{p}_{3}=\bm{k}_{2}+\bm{k}_{3}+\bm{p}_{1},\
\bm{p}_{4}=\bm{k}_{1}+\bm{k}_{4}-\bm{p}_{1}$. Here, $\mathcal{T}$ denotes the
function with time integration

$\displaystyle\mathcal{T}\left(\xi_{*},x_{*},\delta,\sqrt{p_{i}{{}^{*}}}+\sqrt{|\bm{k}_{i}-\bm{p}_{i}|^{*}}\right)$
$\displaystyle\equiv\int_{\infty}^{0}\dfrac{dx^{\prime}_{i}}{x_{i}^{\prime}}\dfrac{2\xi_{*}}{\left(\tfrac{x_{*}}{x_{i}^{\prime}}\right)^{\delta}+\left(\tfrac{x_{i}^{\prime}}{x_{*}}\right)^{\delta}}\sqrt{\dfrac{\pi}{2}}\dfrac{J_{3/2}(k_{i}^{*}x_{i}^{\prime})}{(k_{i}^{*}x_{i}{{}^{\prime}})^{3/2}}$
$\displaystyle\times\int_{\infty}^{x_{i}^{\prime}}dx^{\prime\prime}_{i}(x_{i}^{\prime}x_{i}^{\prime\prime})^{3/2}\dfrac{\pi}{2}\left[J_{3/2}(k_{i}^{*}x_{i}^{\prime})Y_{3/2}(k_{i}^{*}x_{i}^{\prime\prime})-Y_{3/2}(k_{i}^{*}x_{i}^{\prime})J_{3/2}(k_{i}^{*}x_{i}^{\prime\prime})\right]$
$\displaystyle\times\exp\left[-\dfrac{4\sqrt{\xi_{*}x^{\prime\prime}_{i}}}{1+\delta}\left(\dfrac{x^{\prime\prime}_{i}}{x_{*}}\right)^{\delta/2}\left(\sqrt{p_{i}{{}^{*}}}+\sqrt{|\bm{k}_{i}-\bm{p}_{i}|^{*}}\right){{}_{2}}F_{1}\left(\tfrac{1}{2},\
\tfrac{1+\delta}{4\delta};\ \tfrac{5\delta+1}{4\delta};\
-\left(\tfrac{x_{i}^{\prime\prime}}{x_{*}}\right)^{2\delta}\right)\right]\ .$
(4.7)

In the above integration, we can find that the time integration can cover most
of the interval when $p_{i}^{*}\sim|\bm{k}_{i}-\bm{p}_{i}|^{*}=\mathcal{O}(1)$
at which the two $\mathcal{A}$’s shown in Eq. (4.5) are simultaneously at
their peaks. Therefore, the trispectrum is expected to be mostly amplified at
around the equilateral limit: $k_{1}=k_{2}=k_{3}=k_{4}$. Hereafter, we
evaluate the trispctrum signal around there.

Figure 2:  Our momentum configuration in the equilateral limit. The four
momenta have the same length, $k_{1}=k_{2}=k_{3}=k_{4}$. $\bm{k}_{1}$ and
$\bm{k}_{2}$ are in the $x\text{-}y$ plane (green plane). $\bm{k}_{3}$ and
$\bm{k}_{4}$ lie in another plane at an angle of $\phi$ to the $x\text{-}y$
plane around the $x$-axis (blue plane). $\bm{k}_{1}+\bm{k}_{2}$ and
$-\bm{k}_{3}-\bm{k}_{4}$ are parallel to the $x$-axis. Each of
$\bm{k}_{1},\bm{k}_{2},-\bm{k}_{3}$ and $-\bm{k}_{4}$ makes an angle of
$\theta$ with the $x$ axis.

##  5 Results

In this section, we provide the numerical results of the trispectrum of the
sourced curvature perturbation and discuss its parity violation. We consider
two different configurations of the four momenta; the exact equilateral and
quasi-equilateral shapes.

###  5.1 Exact equilateral shape

To investigate the trispectrum, we consider four momenta $\bm{k}_{1}$,
$\bm{k}_{2}$, $\bm{k}_{3}$, and $\bm{k}_{4}$ satisfying
$\bm{k}_{1}+\bm{k}_{2}+\bm{k}_{3}+\bm{k}_{4}=\bm{0}$. Note that we cannot
generally choose the coordinate system where all of the momenta lie in the
$x$-$y$ plane unlike the case of power spectrum and bispectrum. We first
consider the exact equilateral limit, $k_{1}=k_{2}=k_{3}=k_{4}\equiv k$, and
use the parameterization of [88] as shown in Fig. 2:

$\displaystyle\bm{k}_{1}$ $\displaystyle=k(\cos\theta,\,\sin\theta,0),$ (5.1)
$\displaystyle\bm{k}_{2}$ $\displaystyle=k(\cos\theta,\,-\sin\theta,0),$ (5.2)
$\displaystyle\bm{k}_{3}$
$\displaystyle=k(-\cos\theta,\,\sin\theta\cos\phi,\,\sin\theta\sin\phi),$
(5.3) $\displaystyle\bm{k}_{4}$
$\displaystyle=k(-\cos\theta,\,-\sin\theta\cos\phi,\,-\sin\theta\sin\phi).$
(5.4)

In the equilateral limit, we can always take this coordinate system by setting
$\hat{\bm{x}}\parallel\bm{k}_{1}+\bm{k}_{2}$ and
$\hat{\bm{z}}\perp\bm{k}_{1}$. Without loss of generality, we can set
$0<\theta<\pi/2$ and $0\leq\phi<\pi$ by exchanging $\bm{k}_{i}$. Note that
this coordinate system avoids $\bm{k}_{i}\parallel\pm\hat{\bm{z}}$, where the
phases of polarization vector and tensor are not well-defined.

In the equilateral configuration, the trispectrum $T_{\zeta}(\bm{k}_{1},\
\bm{k}_{2},\ \bm{k}_{3},\ \bm{k}_{4})$ in Eq. (4.6) is labeled by the two
angles $(\theta,\phi)$ and thus we call it $T^{(\theta,\phi)}_{\zeta,\rm
eq}(k)$. Then we normalize it as

$T^{(\theta,\phi)}_{\zeta,\rm
eq}(k)=P_{\zeta,v}(k)^{3}f^{(\theta,\phi)}_{\zeta,4}(k/k_{*},\xi_{*},\delta)\
,$ (5.5)

where the power spectrum of the vacuum curvature perturbation,
$P_{\zeta,v}=2\pi^{2}\mathcal{P}_{\zeta,v}/k^{3}$, is used. The dimensionless
trispectrum is given by

$\displaystyle f^{(\theta,\phi)}_{\zeta,4}(k/k_{*},\xi_{*},\delta)$
$\displaystyle=\dfrac{3^{4}\pi^{2}r_{v}^{4}}{2^{7}}\mathcal{P}_{\zeta,v}\int\dfrac{d\bm{p}^{*}_{1}}{(2\pi)^{3}}\prod_{i=1,2,3,4}e^{+}_{m}(\hat{\bm{p}_{i}})e^{+}_{m}(\savestack{\tmpbox}{\stretchto{\scaleto{\scalerel*[\widthof{\bm{k}_{i}-\bm{p}_{i}}]{\kern-0.6pt\bigwedge\kern-0.6pt}{\rule[-642.0pt]{4.30554pt}{642.0pt}}}{}}{0.5ex}}\stackon[1pt]{\bm{k}_{i}-\bm{p}_{i}}{\tmpbox})\left(\sqrt{p_{i}{{}^{*}}}+\sqrt{|\bm{k}_{i}-\bm{p}_{i}|^{*}}\right)\left(p^{*}_{i}|\bm{k}_{i}-\bm{p}_{i}|^{*}\right)$
$\displaystyle\times\mathcal{T}\left(\xi_{*},x_{*},\delta,\sqrt{p_{i}{{}^{*}}}+\sqrt{|\bm{k}_{i}-\bm{p}_{i}|^{*}}\right)N(\xi_{*},p^{*}_{i}x_{*},\delta)N(\xi_{*},|\bm{k}_{i}-\bm{p}_{i}|^{*}x_{*},\delta)$
$\displaystyle+(\text{permutations of}\ \bm{k}_{3}^{*}-\bm{p}_{3}^{*},\
\bm{k}_{4}^{*}-\bm{p}_{4}^{*})+(\text{permutations of}\ \bm{p}_{2}^{*},\
\bm{p}_{3}^{*})\ .$ (5.6)

We numerically compute this factor. Before taking a closer look at it,
however, we first take an overview of its magnitude.

To characterize the trispectrum averaged with respect to the angular
variables, we introduce an averaged $f^{(\theta,\phi)}_{\zeta,4}$ as

$\bar{f}_{\zeta,4}(k/k_{*},\xi_{*},\delta)\equiv
P_{\zeta,v}(k)^{-3}\sqrt{\dfrac{1}{N_{\rm
bin}}\sum_{i,j}\left[T^{(\theta_{i},\phi_{j})}_{\zeta,\rm
eq}(k)\right]^{2}}\,,$ (5.7)

where we have discretized an angular space $(\theta,\phi)$ into
$(\theta_{i},\phi_{j})$ with $(i=1,...,i_{\rm max},\ j=1,...,j_{\rm max})$.
$N_{\rm bin}=i_{\rm max}\times j_{\rm max}$ is a number of bins that we
numerically evaluate. We separate $\bar{f}_{\zeta,4}$ into real and imaginary
parts:

$\bar{f}_{\zeta,4}(k/k_{*},\xi_{*},\delta)\equiv
F_{\zeta,4}(k/k_{*},\xi_{*},\delta)+iG_{\zeta,4}(k/k_{*},\xi_{*},\delta)\ .$
(5.8)

We present the spectral dependence of the averaged trispectum in Figure 3.
There we have fixed the model parameters as follows: $\xi_{*}=4,\ \delta=0.3$
in the left panel, $\xi_{*}=4.5,\ \delta=0.6$ in the right panel, and
$r_{v}=10^{-3.5},\mathcal{P}_{\zeta,v}=10^{-9}$ in the both panels. These
parameters are consistent with the constraint on the power spectrum of
axion-U(1) with $k_{*}=5\times 10^{-3}\ \text{Mpc}^{-1}$ for instance [28]. As
we have discussed in Sec. 2, $\delta$ controls the duration over which
$\xi(t)$ is significantly large. A greater $\delta$ shortens the period of
tachyonic instability caused by axion. Accordingly, as shown in Fig. 3, when
$\delta=0.6$, the width of the averaged trispectrum is narrower than when
$\delta=0.3$. We also find from this figure that the parity-odd component of
the trispectrum peaks at around $k/k_{*}\sim 10$ and the ratio between the odd
and even component is roughly $1\%$ at the peak.

Figure 3:  The real part $F_{\zeta,4}$ (blue line) and the imaginary part
$G_{\zeta,4}$ (orange line) of the normalized trispectrum in the equilateral
limit that is averaged over the angular directions defined in Eq. (5.8). The
model parameters are chosen as $\delta=0.3,\xi_{*}=4$ in the left and
$\delta=0.6,\xi_{*}=4.5$ in the right panel. In the both panels, we use
$r_{v}=10^{-3.5}$ and $i_{\rm max}=j_{\rm max}=10$. When the angular average
is taken, the parity even part (blue) is typically hundred times larger than
the parity odd part (orange).

Next, we investigate the angular dependence of the trispectrum at around the
spectral peak, $k/k_{*}=10$. In Fig. 4, we show the amplitude of ${\rm
Re}[f_{\zeta}]$ and ${\rm Im}[f_{\zeta}]$ as a function of $(\phi,\theta)$. We
choose the parameters as $\xi_{*}=4.5,\ \delta=0.6$ and $r_{v}=10^{-3.5}$. In
the case of $|\bm{k}_{1}+\bm{k}_{2}|\to 0$, so called collapsed limit that
corresponds to $\theta\to\pi/2$, the polarization vector has a planar
structure, so the parity-odd signature vanishes. Similarly, the imaginary part
vanishes in $\theta=0$ and $\phi=0,\pi$ as well. Note that the real part of
$f_{\zeta}$ becomes negative around $\pi/4\lesssim\phi\lesssim 3\pi/4$ and
$\pi/6\lesssim\theta\lesssim 5\pi/12$, but otherwise it is positive. Moreover,
corresponding to the even (odd) parity condition, ${\rm Re}[f_{\zeta}]$ (${\rm
Im}[f_{\zeta}]$) becomes symmetric (antisymmetic) with respect to the
$\phi=\pi/2$ line. These mathematical features are apparent from Fig. 4. Here,
we again confirm that the parity-even part ${\rm Re}[f_{\zeta}]$ basically
dominates the trispectrum signal.

In Fig. 5, we plot the fraction of ${\rm Im}[f_{\zeta}]/{\rm Re}[f_{\zeta}]$
as a function of $(\phi,\theta)$, and evaluate the relative size of the
parity-odd part. Note that this fraction becomes diverging on the zero-
crossing line of the parity-even part ${\rm Re}[f_{\zeta}]$.222In Fig. 5, we
have used a cutoff value of order unity to avoid the computational
divergences. Nonetheless, such areas are not much resolved in Fig. 5 due to a
lack of number of bins. We here find that the parity-odd part still has
$\gtrsim 10\%$ signal of the parity-even one in some specific parameter
regions. The $10\%$ fraction is larger than the result from a model of axion
as an inflaton analyzed in Ref. [85]. The reason is that, in contrast with the
inflaton model, in our spectator model $\xi$ is not directly connected to the
curvature perturbation and dynamically evolves in time, making the value of
$\xi_{*}$ larger than that in Ref. [85].333 In Ref. [85], $\xi\leq 2.4$ are
adopted in the equilateral-limit analysis.

In Fig. 6, we show $f_{\zeta}^{(\theta,\phi)}$ as a function of $\phi$ with
$\theta=\pi/3$. The parity-even component exhibits zero-crossings at around
$\phi\simeq\pi/4$ and $3\pi/4$, and the ratio diverges at the same points. The
parity-odd component peaks at $\phi\simeq\pi/8$ and $7\pi/8$, and the ratio
between the odd and even component reaches $\mathcal{O}(10\%)$ there. However,
again, we do not find any point where the odd component surpasses the even
one, except for the close vicinity of the zero-crossing region.

Figure 4:  Contour plots of the real part ${\rm
Re}[f_{\zeta}^{(\theta,\phi)}]$ (left panel) and the imaginary part ${\rm
Im}[f_{\zeta}^{(\theta,\phi)}]$ (right panel) of the normalized trispectrum in
the equilateral limit $f_{\zeta}^{(\theta,\phi)}$ introduced in Eq. (5.5). We
numerically compute it based on Eq. (5.6). The color denotes the amplitude of
each value. In this plot, we use $\xi_{*}=4.5,\ \delta=0.6,\
r_{v}=10^{-3.5},k/k_{*}=10$ and $i_{\rm max}=j_{\rm max}=16$.  Figure 5:
Contour plot of the ratio between the odd and even parity part of the
equilateral trispectrum, ${\rm Im}[f_{\zeta}^{(\theta,\phi)}]/{\rm
Re}[f_{\zeta}^{(\theta,\phi)}]$. In other words, the right panel of Fig. 4 is
divided by the left panel. The parameters are the same as Fig. 4. We see that
this ratio is of the order of unity in small regions highlighted by the red
and blue color.

Figure 6:  The even part ${\rm Re}[f_{\zeta}^{(\theta,\phi)}]$ (left panel),
odd part ${\rm Im}[f_{\zeta}^{(\theta,\phi)}]$ (center panel), and ratio ${\rm
Im}[f_{\zeta}^{(\theta,\phi)}]/{\rm Re}[f_{\zeta}^{(\theta,\phi)}]$ (right
panel) shown as functions of $\phi$ by fixing $\theta=\pi/3$. In other words,
they are the cross-sections between the plane of $\theta=\pi/3$ and the
previous contour plots in Fig. 4 and 5. The parameters are the same as Fig. 4.

###  5.2 Quasi-equilateral shape

In the equilateral limit, the parity-even part is dominant over the parity-odd
part. However, the equilateral configuration is just one of the numerous
possible configurations. In this section, we examine another configuration
which is slightly deviated from the exact equilateral one.

We take the coordinates of the momenta with a more generic configuration as
shown in Fig. 7:

$\displaystyle\bm{k}_{1}$
$\displaystyle=k_{1}(\cos\theta_{1},\,\sin\theta_{1},0),$ (5.9)
$\displaystyle\bm{k}_{2}$
$\displaystyle=k_{2}(\cos\theta_{2},\,-\sin\theta_{2},0),$ (5.10)
$\displaystyle\bm{k}_{3}$
$\displaystyle=k_{3}(-\cos\theta_{3},\,\sin\theta_{3}\cos\phi,\,\sin\theta_{3}\sin\phi),$
(5.11) $\displaystyle\bm{k}_{4}$
$\displaystyle=k_{4}(-\cos\theta_{4},\,-\sin\theta_{4}\cos\phi,\,-\sin\theta_{4}\sin\phi).$
(5.12)

Since two pairs of momenta lie in the same planes as before, we do not need to
re-define the angle $\phi$, but we need to generalize $\theta$ to $\theta_{i}\
(i=1,2,3,4)$ for each of momenta to satisfy the momentum conservation,
$\bm{k}_{1}+\bm{k}_{2}+\bm{k}_{3}+\bm{k}_{4}=\bm{0}$.

The momentum conservation $\bm{k}_{1}+\bm{k}_{2}+\bm{k}_{3}+\bm{k}_{4}=\bm{0}$
leads to the following conditions:

$\displaystyle k_{1}\sin\theta_{1}=k_{2}\sin\theta_{2}\ ,$ (5.13)
$\displaystyle k_{3}\sin\theta_{3}=k_{4}\sin\theta_{4}\ ,$ (5.14)
$\displaystyle
k_{1}\cos\theta_{1}+k_{2}\cos\theta_{2}=k_{3}\cos\theta_{3}+k_{4}\cos\theta_{4}\
.$ (5.15)

It is still tough to cover all considerable shapes. However, we can guess that
a large signal from the generation of the sourced scalar mode is located
around the shapes of the four momenta with roughly the same length, because
the particle production of gauge field happens around the horizon-crossing.
Therefore, we consider a small deviation from equilateral shape; namely, the
case where only one momentum has a different magnitude from the other three
momenta

$k_{1}=k_{2}=k_{3}\equiv k,\quad k_{4}\neq k\ .$ (5.16)

Substituting the condition (5.16) and defining $m\equiv k_{4}/k$, the above
conditions are reduced to

$\displaystyle\theta_{1}=\theta_{2}\equiv\theta\ ,$ (5.17)
$\displaystyle\sin\theta_{3}=m\sin\theta_{4}\ ,$ (5.18)
$\displaystyle\cos\theta=\dfrac{1}{2}\left(\cos\theta_{3}+m\cos\theta_{4}\right)\
.$ (5.19)

Then, we can express $\theta_{3}$ and $\theta_{4}$ in terms of $\theta$ and
$m$ as

$\cos\theta_{3}=\cos\theta+\dfrac{1-m^{2}}{4\cos\theta}\
,\quad\cos\theta_{4}=\dfrac{1}{m}\left(\cos\theta-\dfrac{1-m^{2}}{4\cos\theta}\right)\
.$ (5.20)

Thus, our angular variables are $\\{\theta,\phi,m\\}$. Notice that the values
of $\theta$ and $m$ should be taken to satisfy $0<\cos\theta_{3,4}<1$.

Figure 7:  Our momentum configuration in the quasi-equilateral shape. In the
same way as the equilateral configuration in Fig. 2, $\bm{k}_{1}$ and
$\bm{k}_{2}$ are in the $x\text{-}y$ plane (green plane), and $\bm{k}_{3}$ and
$\bm{k}_{4}$ lie in another plane at an angle of $\phi$ to the $x\text{-}y$
plane around the $x$-axis (blue plane). However, the fourth momentum is longer
than the other momenta, $k_{4}>k_{1}=k_{2}=k_{3}$. To keep the momentum
conservation, $\bm{k}_{1}+\bm{k}_{2}+\bm{k}_{3}+\bm{k}_{4}=0$, the angle to
the $x$-axis is generalized from the common value $\theta$ to individually
different ones $\theta_{i}\,(i=1,2,3,4)$.

We numerically compute the dimensionless trispectrum
$f_{\zeta}^{(\theta,\phi)}$ in the same way as the equilateral limit case. In
Fig. 8, we show the amplitude of ${\rm Re}[f_{\zeta}]$ and ${\rm
Im}[f_{\zeta}]$ as a function of $(\phi,\theta)$. We choose the momenta around
the inflection point $k/k_{*}=10$ and parameters are $\xi_{*}=4.5,\
\delta=0.6,\ r_{v}=10^{-3.5},$ and $m=1.3$. Regarding the even part, its
magnitude becomes smaller in comparison with the equilateral case because the
contribution from the $\bm{k}_{4}$ mode is slightly off peak and the amplitude
becomes smaller than the others. As a result, we can observe that, in this
quasi-equilateral configuration, the real and imaginary parts are of the same
order of magnitude contrary to the case of the equilateral shape. It can be
clearly seen in Fig. 9 that shows their ratio. In the darkest red and blue
regions of Fig. 9, the absolute value of the ratio, $|{\rm
Im}[f_{\zeta}^{(\theta,\phi)}]/{\rm Re}[f_{\zeta}^{(\theta,\phi)}]|$, exceeds
unity. Comparing Fig. 4 with Fig. 8, one also notices that the real part ${\rm
Re}[f_{\zeta}]$ (the imaginary part ${\rm Im}[f_{\zeta}]$) is no longer
symmetric (antisymmetric) with respect to $\phi=\pi/2$ in the quasi-
equilateral case.

In Fig. 10, we present the real and imaginary parts of
$f_{\zeta}^{(\theta,\phi)}$ and their ratio as a function of $\phi$ by fixing
$\theta=\pi/3$. The parity-even signal crosses zero at around
$\phi\simeq\pi/5,4\pi/5$, and the ratio diverges at the same points. Contrary
to the equilateral shape case, however, the absolute ratio $|{\rm
Im}[f_{\zeta}^{(\theta,\phi)}]/{\rm Re}[f_{\zeta}^{(\theta,\phi)}]|$ exceeds
unity even in regions away from the zero-crossing points. For example,
$\phi=3\pi/4$ is close to the negative maximum point of the even component,
while the ratio is unity.444 A similar amplification of the parity-odd
component in non-equilateral configurations is confirmed also at the
bispectrum level [89]. It is also interesting to note that the point where the
even component crosses zero and the odd component peaks are roughly
coincident. These results imply that the angular dependence of the trispectrum
has non-trivial structures and is important to distinguish the origin of the
parity-violating signals.

Figure 8:  Contour plot of ${\rm Re}[f_{\zeta}^{(\theta,\phi)}]$ (left panel)
and ${\rm Im}[f_{\zeta}^{(\theta,\phi)}]$ (right panel) in the quasi-
equilateral case. The parameters are the same as Fig. 4, except for $m\equiv
k_{4}/k=1.3$.  Figure 9:  Contour plot of ${\rm
Im}[f_{\zeta}^{(\theta,\phi)}]/{\rm Re}[f_{\zeta}^{(\theta,\phi)}]$ in the
quasi-equilateral case. We set the threshold at $\sim\pm 1$, because the ratio
diverges when the real part ${\rm Re}[f_{\zeta}]$ crosses zero. In this plot,
we use $\xi_{*}=4.5,\ \delta=0.6,\ r_{v}=10^{-3.5},\ m=1.3$ and $i_{\rm
max}=j_{\rm max}=16$. We see that the ratio becomes of the order of unity in
larger regions than the equilateral case.

Figure 10:  Plot of dimensionless factor of trispectrum (5.6) as a function of
$\phi$ in the quasi-equilateral case. These are even part ${\rm
Re}[f_{\zeta}^{(\theta,\phi)}]$ (left panel), odd part ${\rm
Im}[f_{\zeta}^{(\theta,\phi)}]$ (central panel), and ratio ${\rm
Im}[f_{\zeta}^{(\theta,\phi)}]/{\rm Re}[f_{\zeta}^{(\theta,\phi)}]$ (right
panel) respectively. We fixed it as $\theta=\pi/3$. In this plot, we use
$\xi_{*}=4.5,\ \delta=0.6,\ r_{v}=10^{-3.5}$.

##  6 Discussion and conclusion

In this paper, we have studied the generation of the trispectrum sourced by
the rolling axion coupled to U$(1)$ gauge field during inflation. In our
setup, one helicity mode of the gauge field is amplified by the tachyonic
instability, and the axion perturbation is sourced by the gauge field at one-
loop level through the Chern-Simons coupling. Then, the inflaton perturbation
is amplified by the axion perturbation through gravitational coupling. This
model generates sizable equilateral-type non-Gaussianity in the primordial
curvature perturbation.

We have first studied the exact equilateral shape of the curvature
trispectrum. Regarding the angular-averaged trispectrum, we have found that
both parity-even and parity-odd components have peaks at slightly smaller
scale than one crossing the horizon when the axion pass through the inflection
point in the potential, and the amplitude of parity-odd mode is two orders of
magnitude smaller than that of parity-even mode. We have also studied the
angular dependence of the trispectrum near the peak scale. We have found that
the ratio of parity-even to parity-odd mode amplitude has sharp peaks at
certain angles. This is because parity-even mode crosses zero, and the ratio
becomes large near that point. Apart from such zero-crossing points, the ratio
of parity-odd to parity-even mode is typically ${\cal O}(10^{-1})$. Therefore,
the parity-odd mode is non-negligible but not dominant compared to the parity-
even mode in the equilateral configuration.

We then have extended the analysis to quasi equilateral configurations, e.g.,
where only one of the momenta is slightly longer than the other three. In
comparison with the exact equilateral case, interestingly, the ratio of
parity-odd to parity-even mode exceeds unity in wider regions of angular
momentum variables including the area far from the zero-crossing points of the
parity-even mode. This implies the powerfullness of the quasi-equilateral
signals for assessing this model. To clarify this, and moreover, to test with
observed datasets of the CMB or large-scale structure, more comprehensive
analysis covering the whole momentum space is necessary, and we leave it to
our future work.

Finally, let us compare our results with those obtained in a previous similar
but different study [85] where the gauge field is directly coupled to an
axionic inflaton, and the speed of background axion field was assumed to be
constant. Ref. [85] fixed $\theta=\pi/3$ in equilateral form and studied the
angular dependence of $\phi$. The trispectra of ours and in Ref. [85] have
similar features such as the parity-even mode is an even function with
$\phi=\pi/2$, while the parity-odd mode is an odd function. Also, for
$\phi=0,\pi$, the parity-odd mode becomes zero because the momenta vector has
a planar structure. On the other hand, differently from our case, the parity-
even mode does not cross zero in Ref. [85]; therefore, the ratio never becomes
very large. Moreover, in our results, the ratio of parity-odd to parity-even
mode at its maximum is ${\cal O}(10^{-1})$ apart from the zero-crossing
regions of the parity-even mode, while in Ref. [85], it is only ${\cal
O}(10^{-2})$. Therefore, even excluding the fact that our results exhibit
zero-crossings and we studied the quasi-equilateral shape, we found that the
parity-odd signal is larger than what was evaluated in Ref. [85]. The reason
is that our spectator model imprints unique scale dependence in correlation
functions, and it allows us to make the value of $\xi_{*}$ larger than the
inflaton model in Ref. [85] (an upper bound $\xi\leq 2.4$ is used in Ref.
[85], while this paper obeyed observational bounds in Ref. [28]). In light of
this, our scenario may be more useful on a physical interpretation of the
large parity-odd trispectrum signal observed in the BOSS data.

## Acknowledgments

This work is supported by the Japan Society for the Promotion of Science
(JSPS) KAKENHI Grant Nos. JP20H05854 (T. F.), JP23K03424 (T. F.), JP23KJ2007
(T. M.), 19K14702 (I. O.), JP20H05859 (I. O. and M. S.) and JP23K03390 (M.
S.).

##  Appendix A WKB solution

The differential equation (2.11) is unsolvable in a closed form due to the
dynamical $\xi(\tau)$. Here we find the approximate solution of $A^{+}_{k}$ by
relying on the WKB method. The equation of motion is rewritten as

$\left[\partial_{\tau}^{2}+\omega(\tau)^{2}\right]A^{+}_{k}(\tau)=0\
,\qquad\omega(\tau)^{2}\equiv k^{2}+\dfrac{2k\xi(\tau)}{\tau}\ .$ (A.1)

The WKB approximation is valid when the mode function adiabatically evolves in
time: $|\omega^{\prime}(\tau)|\ll\omega(\tau)^{2}$. Since we are interested in
a regime where $\omega(\tau)^{2}$ becomes negative and the tachyonic
instability happens, we consider a time window $\tau>\bar{\tau}$ at which
$\omega(\bar{\tau})^{2}=0$. Then, defining
$\rho(\tau)^{2}\equiv-\omega(\tau)^{2}$, namely

$\rho(\tau)\equiv\sqrt{-\dfrac{2k\xi(\tau)}{\tau}-k^{2}}\ ,$ (A.2)

we obtain a WKB solution

$\displaystyle
A^{+}_{k}(\tau>\bar{\tau})\simeq\dfrac{C_{1}}{\sqrt{\rho(\tau)}}\exp\left[\int^{\bar{\tau}}_{\tau}\rho(\tilde{\tau})d\tilde{\tau}\right]+\dfrac{C_{2}}{\sqrt{\rho(\tau)}}\exp\left[-\int^{\bar{\tau}}_{\tau}\rho(\tilde{\tau})d\tilde{\tau}\right]$
(A.3)

in the regime $|\omega^{\prime}(\tau)|\ll\omega(\tau)^{2}$. The integration
constant $C_{1,2}$ is given by $C_{1}=1/\sqrt{2k},\ C_{2}=-i/\sqrt{2k}$ from a
boundary condition in the adiabatic vacuum. To compute the integration
function, we separate the time integral $\int_{\tau}^{\bar{\tau}}$ into two
parts: $\int_{0}^{\bar{\tau}}+\int_{\tau}^{0}$, where the first part is a
constant term and the second part leads to a time-dependent function on
$\tau$. Using the slow-roll solution for $\xi(\tau)$ (2.11), the time integral
is performed as

$\int^{0}_{\tau}\rho(\tilde{\tau})d\tilde{\tau}\simeq\dfrac{4\xi_{*}^{1/2}}{1+\delta}\left(\dfrac{\tau}{\tau_{*}}\right)^{\delta/2}(-k\tau)^{1/2}{{}_{2}}F_{1}\left(\tfrac{1}{2},\
\tfrac{1+\delta}{4\delta};\ \tfrac{5\delta+1}{4\delta};\
-\left(\tfrac{\tau}{\tau_{*}}\right)^{2\delta}\right)\ .$ (A.4)

We notice that ${{}_{2}}F_{1}$ becomes unity in the limit
$|\tau/\tau_{*}|\rightarrow 0$ and (A.4) then corresponds to what has been
derived in the previous work [41]. Finally, the growing mode of WKB solution
is expressed as

$\displaystyle A^{+}_{k}$
$\displaystyle\simeq\dfrac{N(\xi_{*},-k\tau_{*},\delta)}{\sqrt{2k}}\left(\dfrac{-k\tau}{2\xi(\tau)}\right)^{1/4}$
$\displaystyle\quad\times\exp\left[-\dfrac{4\xi_{*}^{1/2}}{1+\delta}\left(\dfrac{\tau}{\tau_{*}}\right)^{\delta/2}(-k\tau)^{1/2}{{}_{2}}F_{1}\left(\tfrac{1}{2},\
\tfrac{1+\delta}{4\delta};\ \tfrac{5\delta+1}{4\delta};\
-\left(\tfrac{\tau}{\tau_{*}}\right)^{2\delta}\right)\right]\ .$ (A.5)

In Fig. 11, we depict the time evolution of mode function and compare the
numerical solution with the approximate WKB solutions. The WKB solution (A.5)
is well fitted with the exact solution for the whole momentum modes. Comparing
that, the reduced WKB solution is not good for lower values of $-k\tau_{*}$
where the approximation $|\tau/\tau_{*}|\ll 1$ is not so accurate, while it’s
a well approximation for momentum modes where the amplification mostly
happens.

Figure 11:  The time evolution of mode function with $-k\tau_{*}=10^{-2}$
(left panel), $-k\tau_{*}=1$ (middle panel) and $-k\tau_{*}=10^{2}$ (right
panel). The horizontal axis is a dimensionless time flowing from the right to
the left. The black lines denote the exact numerical solutions. The orange
lines are the WKB approximate solution (A.5). The blue lines are the solution
(A.5) where ${{}_{2}}F_{1}$ is replaced with $1$. We set the parameter
$\xi_{*}=4$ and $\delta=0.3$.

## References

  * Lue  _et al._ [1999] A. Lue, L.-M. Wang, and M. Kamionkowski, Phys. Rev. Lett. 83, 1506 (1999), arXiv:astro-ph/9812088 .
  * Minami and Komatsu [2020] Y. Minami and E. Komatsu, Phys. Rev. Lett. 125, 221301 (2020), arXiv:2011.11254 [astro-ph.CO] .
  * Diego-Palazuelos  _et al._ [2022] P. Diego-Palazuelos _et al._ , Phys. Rev. Lett. 128, 091302 (2022), arXiv:2201.07682 [astro-ph.CO] .
  * Eskilt [2022] J. R. Eskilt, Astron. Astrophys. 662, A10 (2022), arXiv:2201.13347 [astro-ph.CO] .
  * Eskilt and Komatsu [2022] J. R. Eskilt and E. Komatsu, (2022), arXiv:2205.13962 [astro-ph.CO] .
  * Eskilt  _et al._ [2023a] J. R. Eskilt _et al._ , (2023a), arXiv:2305.02268 [astro-ph.CO] .
  * Carroll [1998] S. M. Carroll, Phys. Rev. Lett. 81, 3067 (1998), astro-ph/9806099 .
  * Panda  _et al._ [2011] S. Panda, Y. Sumitomo, and S. P. Trivedi, Phys. Rev. D 83, 083506 (2011), 1011.5877 .
  * Fujita  _et al._ [2021a] T. Fujita, Y. Minami, K. Murai, and H. Nakatsuka, Phys. Rev. D 103, 063508 (2021a), arXiv:2008.02473 [astro-ph.CO] .
  * Fujita  _et al._ [2021b] T. Fujita, K. Murai, H. Nakatsuka, and S. Tsujikawa, Phys. Rev. D 103, 043509 (2021b), arXiv:2011.11894 [astro-ph.CO] .
  * Choi  _et al._ [2021] G. Choi, W. Lin, L. Visinelli, and T. T. Yanagida, Phys. Rev. D 104, L101302 (2021), arXiv:2106.12602 [hep-ph] .
  * Obata [2022] I. Obata, JCAP 09, 062 (2022), 2108.02150 .
  * Gasparotto and Obata [2022] S. Gasparotto and I. Obata, JCAP 08, 025 (2022), arXiv:2203.09409 [astro-ph.CO] .
  * Galaverni  _et al._ [2023] M. Galaverni, F. Finelli, and D. Paoletti, Phys. Rev. D 107, 083529 (2023), arXiv:2301.07971 .
  * Gasparotto and Sfakianakis [2023] S. Gasparotto and E. I. Sfakianakis, (2023), arXiv:2306.16355 [astro-ph.CO] .
  * Murai  _et al._ [2023] K. Murai, F. Naokawa, T. Namikawa, and E. Komatsu, Phys. Rev. D 107, L041302 (2023), 2209.07804 .
  * Eskilt  _et al._ [2023b] J. R. Eskilt, L. Herold, E. Komatsu, K. Murai, T. Namikawa, and F. Naokawa, (2023b), arXiv:2303.15369 [astro-ph.CO] .
  * Finelli and Galaverni [2009] F. Finelli and M. Galaverni, Phys. Rev. D 79, 063002 (2009), 0802.4210 .
  * Liu and Ng [2017] G.-C. Liu and K.-W. Ng, Phys. Dark Univ. 16, 22 (2017), 1612.02104 .
  * Fedderke  _et al._ [2019] M. A. Fedderke, P. W. Graham, and S. Rajendran, Phys. Rev. D 100, 015040 (2019), 1903.02666 .
  * Takahashi and Yin [2021] F. Takahashi and W. Yin, JCAP 04, 007 (2021), arXiv:2012.11576 [hep-ph] .
  * Kitajima  _et al._ [2022] N. Kitajima, F. Kozai, F. Takahashi, and W. Yin, JCAP 10, 043 (2022), arXiv:2205.05083 [astro-ph.CO] .
  * Jain  _et al._ [2022] M. Jain, R. Hagimoto, A. J. Long, and M. A. Amin, JCAP 10, 090 (2022), arXiv:2208.08391 [astro-ph.CO] .
  * Gonzalez  _et al._ [2022] D. Gonzalez, N. Kitajima, F. Takahashi, and W. Yin, (2022), arXiv:2211.06849 [hep-ph] .
  * Saito  _et al._ [2007] S. Saito, K. Ichiki, and A. Taruya, JCAP 09, 002 (2007), arXiv:0705.3701 [astro-ph] .
  * Gerbino  _et al._ [2016] M. Gerbino, A. Gruppuso, P. Natoli, M. Shiraishi, and A. Melchiorri, JCAP 07, 044 (2016), arXiv:1605.09357 [astro-ph.CO] .
  * Fujita  _et al._ [2022a] T. Fujita, Y. Minami, M. Shiraishi, and S. Yokoyama, Phys. Rev. D 106, 103529 (2022a), arXiv:2208.08101 [astro-ph.CO] .
  * Campeti  _et al._ [2022] P. Campeti, O. Özsoy, I. Obata, and M. Shiraishi, JCAP 07, 039 (2022), arXiv:2203.03401 [astro-ph.CO] .
  * Shiraishi  _et al._ [2015] M. Shiraishi, M. Liguori, and J. R. Fergusson, JCAP 01, 007 (2015), arXiv:1409.0265 [astro-ph.CO] .
  * Ade  _et al._ [2016] P. A. R. Ade _et al._ (Planck), Astron. Astrophys. 594, A17 (2016), arXiv:1502.01592 [astro-ph.CO] .
  * Akrami  _et al._ [2020] Y. Akrami _et al._ (Planck), Astron. Astrophys. 641, A9 (2020), arXiv:1905.05697 [astro-ph.CO] .
  * Sorbo [2011] L. Sorbo, JCAP 1106, 003 (2011), arXiv:1101.1525 [astro-ph.CO] .
  * Barnaby  _et al._ [2011] N. Barnaby, R. Namba, and M. Peloso, JCAP 04, 009 (2011), arXiv:1102.4333 [astro-ph.CO] .
  * Cook and Sorbo [2012] J. L. Cook and L. Sorbo, Phys. Rev. D85, 023534 (2012), [Erratum: Phys. Rev.D86,069901(2012)], arXiv:1109.0022 [astro-ph.CO] .
  * Barnaby  _et al._ [2012a] N. Barnaby, E. Pajer, and M. Peloso, Phys. Rev. D 85, 023525 (2012a), arXiv:1110.3327 [astro-ph.CO] .
  * Barnaby  _et al._ [2012b] N. Barnaby, J. Moxon, R. Namba, M. Peloso, G. Shiu, and P. Zhou, Phys. Rev. D86, 103508 (2012b), arXiv:1206.6117 [astro-ph.CO] .
  * Cook and Sorbo [2013] J. L. Cook and L. Sorbo, JCAP 1311, 047 (2013), arXiv:1307.7077 [astro-ph.CO] .
  * Mukohyama  _et al._ [2014] S. Mukohyama, R. Namba, M. Peloso, and G. Shiu, JCAP 08, 036 (2014), arXiv:1405.0346 [astro-ph.CO] .
  * Bartolo  _et al._ [2015a] N. Bartolo, S. Matarrese, M. Peloso, and M. Shiraishi, JCAP 01, 027 (2015a), arXiv:1411.2521 [astro-ph.CO] .
  * Bartolo  _et al._ [2015b] N. Bartolo, S. Matarrese, M. Peloso, and M. Shiraishi, JCAP 07, 039 (2015b), arXiv:1505.02193 [astro-ph.CO] .
  * Namba  _et al._ [2016] R. Namba, M. Peloso, M. Shiraishi, L. Sorbo, and C. Unal, JCAP 01, 041 (2016), arXiv:1509.07521 [astro-ph.CO] .
  * Domcke  _et al._ [2016] V. Domcke, M. Pieroni, and P. Binétruy, JCAP 06, 031 (2016), arXiv:1603.01287 [astro-ph.CO] .
  * Garcia-Bellido  _et al._ [2016] J. Garcia-Bellido, M. Peloso, and C. Unal, JCAP 12, 031 (2016), arXiv:1610.03763 [astro-ph.CO] .
  * Obata [2017] I. Obata, JCAP 06, 050 (2017), arXiv:1612.08817 [astro-ph.CO] .
  * Domcke  _et al._ [2017] V. Domcke, F. Muia, M. Pieroni, and L. T. Witkowski, JCAP 07, 048 (2017), arXiv:1704.03464 [astro-ph.CO] .
  * Garcia-Bellido  _et al._ [2017] J. Garcia-Bellido, M. Peloso, and C. Unal, JCAP 09, 013 (2017), arXiv:1707.02441 [astro-ph.CO] .
  * Özsoy [2021] O. Özsoy, JCAP 04, 040 (2021), arXiv:2005.10280 [astro-ph.CO] .
  * Caravano  _et al._ [2021] A. Caravano, E. Komatsu, K. D. Lozanov, and J. Weller, JCAP 12, 010 (2021), arXiv:2102.06378 [astro-ph.CO] .
  * Caravano  _et al._ [2022] A. Caravano, E. Komatsu, K. D. Lozanov, and J. Weller, Phys. Rev. D 105, 123530 (2022), arXiv:2110.10695 [astro-ph.CO] .
  * Caravano  _et al._ [2023] A. Caravano, E. Komatsu, K. D. Lozanov, and J. Weller, Phys. Rev. D 108, 043504 (2023), arXiv:2204.12874 [astro-ph.CO] .
  * Niu and Rahat [2023] X. Niu and M. H. Rahat, (2023), arXiv:2307.01192 [hep-ph] .
  * Unal  _et al._ [2023] C. Unal, A. Papageorgiou, and I. Obata, (2023), arXiv:2307.02322 [astro-ph.CO] .
  * Maleknejad and Sheikh-Jabbari [2011] A. Maleknejad and M. Sheikh-Jabbari, Phys. Rev. D 84, 043515 (2011), arXiv:1102.1932 [hep-ph] .
  * Dimastrogiovanni and Peloso [2013] E. Dimastrogiovanni and M. Peloso, Phys. Rev. D 87, 103501 (2013), arXiv:1212.5184 [astro-ph.CO] .
  * Shiraishi  _et al._ [2013] M. Shiraishi, A. Ricciardone, and S. Saga, JCAP 11, 051 (2013), arXiv:1308.6769 [astro-ph.CO] .
  * Adshead  _et al._ [2013a] P. Adshead, E. Martinec, and M. Wyman, JHEP 09, 087 (2013a), arXiv:1305.2930 [hep-th] .
  * Adshead  _et al._ [2013b] P. Adshead, E. Martinec, and M. Wyman, Phys. Rev. D 88, 021302 (2013b), arXiv:1301.2598 [hep-th] .
  * Dimastrogiovanni  _et al._ [2017] E. Dimastrogiovanni, M. Fasiello, and T. Fujita, JCAP 01, 019 (2017), arXiv:1608.04216 [astro-ph.CO] .
  * Shiraishi  _et al._ [2016] M. Shiraishi, C. Hikage, R. Namba, T. Namikawa, and M. Hazumi, Phys. Rev. D 94, 043506 (2016), arXiv:1606.06082 [astro-ph.CO] .
  * Adshead and Sfakianakis [2017] P. Adshead and E. I. Sfakianakis, JHEP 08, 130 (2017), arXiv:1705.03024 [hep-th] .
  * Fujita  _et al._ [2018] T. Fujita, R. Namba, and Y. Tada, Phys. Lett. B 778, 17 (2018), arXiv:1705.01533 [astro-ph.CO] .
  * Adshead  _et al._ [2016] P. Adshead, E. Martinec, E. I. Sfakianakis, and M. Wyman, JHEP 12, 137 (2016), arXiv:1609.04025 [hep-th] .
  * Agrawal  _et al._ [2018a] A. Agrawal, T. Fujita, and E. Komatsu, Phys. Rev. D97, 103526 (2018a), arXiv:1707.03023 [astro-ph.CO] .
  * Thorne  _et al._ [2018] B. Thorne, T. Fujita, M. Hazumi, N. Katayama, E. Komatsu, and M. Shiraishi, Phys. Rev. D 97, 043506 (2018), arXiv:1707.03240 [astro-ph.CO] .
  * Agrawal  _et al._ [2018b] A. Agrawal, T. Fujita, and E. Komatsu, JCAP 1806, 027 (2018b), arXiv:1802.09284 [astro-ph.CO] .
  * Papageorgiou  _et al._ [2018] A. Papageorgiou, M. Peloso, and C. Unal, JCAP 09, 030 (2018), arXiv:1806.08313 [astro-ph.CO] .
  * Fujita  _et al._ [2019a] T. Fujita, R. Namba, and I. Obata, JCAP 04, 044 (2019a), arXiv:1811.12371 [astro-ph.CO] .
  * Fujita  _et al._ [2019b] T. Fujita, E. I. Sfakianakis, and M. Shiraishi, JCAP 05, 057 (2019b), arXiv:1812.03667 [astro-ph.CO] .
  * Papageorgiou  _et al._ [2019] A. Papageorgiou, M. Peloso, and C. Unal, JCAP 07, 004 (2019), arXiv:1904.01488 [astro-ph.CO] .
  * Campeti  _et al._ [2021] P. Campeti, E. Komatsu, D. Poletti, and C. Baccigalupi, JCAP 01, 012 (2021), arXiv:2007.04241 [astro-ph.CO] .
  * Ishiwata  _et al._ [2022] K. Ishiwata, E. Komatsu, and I. Obata, JCAP 03, 010 (2022), arXiv:2111.14429 [hep-ph] .
  * Maleknejad  _et al._ [2013] A. Maleknejad, M. M. Sheikh-Jabbari, and J. Soda, Phys. Rept. 528, 161 (2013), arXiv:1212.2921 [hep-th] .
  * Komatsu [2022] E. Komatsu, Nature Rev. Phys. 4, 452 (2022), arXiv:2202.13919 [astro-ph.CO] .
  * Shiraishi [2016] M. Shiraishi, Phys. Rev. D 94, 083503 (2016), arXiv:1608.00368 [astro-ph.CO] .
  * Philcox  _et al._ [2021] O. H. E. Philcox, J. Hou, and Z. Slepian, (2021), arXiv:2108.01670 [astro-ph.CO] .
  * Cahn  _et al._ [2023] R. N. Cahn, Z. Slepian, and J. Hou, Phys. Rev. Lett. 130, 201002 (2023), arXiv:2110.12004 [astro-ph.CO] .
  * Hou  _et al._ [2023] J. Hou, Z. Slepian, and R. N. Cahn, Mon. Not. Roy. Astron. Soc. 522, 5701 (2023), arXiv:2206.03625 [astro-ph.CO] .
  * Philcox [2022] O. H. E. Philcox, Phys. Rev. D 106, 063501 (2022), arXiv:2206.04227 [astro-ph.CO] .
  * Cabass  _et al._ [2023a] G. Cabass, M. M. Ivanov, and O. H. E. Philcox, Phys. Rev. D 107, 023523 (2023a), arXiv:2210.16320 [astro-ph.CO] .
  * Creque-Sarbinowski  _et al._ [2023] C. Creque-Sarbinowski, S. Alexander, M. Kamionkowski, and O. Philcox, (2023), arXiv:2303.04815 [astro-ph.CO] .
  * Philcox [2023] O. H. E. Philcox, (2023), arXiv:2303.12106 [astro-ph.CO] .
  * Philcox and Shiraishi [2023] O. H. E. Philcox and M. Shiraishi, (2023), arXiv:2308.03831 [astro-ph.CO] .
  * Liu  _et al._ [2020] T. Liu, X. Tong, Y. Wang, and Z.-Z. Xianyu, JHEP 04, 189 (2020), arXiv:1909.01819 [hep-ph] .
  * Cabass  _et al._ [2023b] G. Cabass, S. Jazayeri, E. Pajer, and D. Stefanyszyn, JHEP 02, 021 (2023b), arXiv:2210.02907 [hep-th] .
  * Niu  _et al._ [2023] X. Niu, M. H. Rahat, K. Srinivasan, and W. Xue, JCAP 05, 018 (2023), arXiv:2211.14324 [hep-ph] .
  * Jazayeri  _et al._ [2023] S. Jazayeri, S. Renaux-Petel, X. Tong, D. Werth, and Y. Zhu, (2023), arXiv:2308.11315 [hep-th] .
  * Stefanyszyn  _et al._ [2023] D. Stefanyszyn, X. Tong, and Y. Zhu, (2023), arXiv:2309.07769 [hep-th] .
  * Fujita  _et al._ [2022b] T. Fujita, K. Murai, I. Obata, and M. Shiraishi, JCAP 01, 007 (2022b), arXiv:2109.06457 [astro-ph.CO] .
  * Shiraishi [2012] M. Shiraishi, JCAP 06, 015 (2012), arXiv:1202.2847 [astro-ph.CO] .


# Internal review on the state of loss of vacuum research on helium cryogenic
systems

Dalesandro, Andrew A

(19 March 2023)

###### Abstract

It is commonplace in helium cryogenic systems to utilize vacuum insulation to
mitigate convective heat transfer to the low temperature fluids. While the
insulating vacuum meaningfully improves cryogenic system thermal performance,
failure of the insulating vacuum results in rapid heat transfer to the
cryogenic helium often resulting in loss of fluids from venting, potential
over-pressure of the system, and risks damage to vessels and equipment. The
nature of the loss of vacuum problem is highly transient and difficult to
predict uniformly across applications, although several efforts have incurred
in recent decades to evaluate the heat transfer mechanism and develop useful
models.

##  1 Introduction

The loss of insulating vacuum on helium cryogenic systems is often considered
a worst-case failure scenario due a variety of confluent factors and risks.
Risks include loss of life, personnel injury, irreparable damage to process
equipment or vacuum vessel cryostat, inoperability of equipment, venting loss
of helium storage, and reduced performance of equipment or insulation.

###  1.1 Vacuum differential to ambient

An external leak or failure would result in loss of insulating vacuum by air
due to a differential pressure between the insulating vacuum and ambient air
of 1 atm; resulting in a high flow of air into the vacuum space. The sudden
rush of air entering the cryostat may damage the vacuum vessel or multi-layer
insulation, reducing thermal performance for future operations.

###  1.2 Vacuum differential to process

It is also possible that the insulating vacuum can be spoiled from the helium
process itself which may result in orders of magnitude higher pressure
differential (helium process dependent) which could significantly increase the
flow rate into the vacuum vessel. In addition to spoiling the vacuum, the
leaking helium entering the cryostat may also excessively cool the cryostat
vacuum vessel, which is typically not rated for low temperatures.

###  1.3 Temperature differential to process

Since cryogenic helium systems are largely fully welded systems inside the
cryostat, loss of vacuum by air is typically the more credible failure. The
large temperature differential between warm air ( 300 K) and helium (e.g.
liquid helium 4.2 K) causes the air influx to rapidly condense and desublimate
on cold helium surfaces.

###  1.4 Low surface heat capacity

Helium cryogenic systems are primarily fabricated from austenitic stainless
steels due to good mechanical properties at low temperatures, including
relatively high ductility. Heat capacity of austenitic stainless steels is low
at low temperatures, such that the material surfaces will rapidly warm to the
temperature of the condensing air on the cold surfaces, increasing the rate of
heat transfer to helium.

###  1.5 Enthalpy of process helium

The enthalpy change of the air from ambient to desublimation is ¿ 20x the
latent heat of helium at its normal boiling point. This results in a rapid
vaporization and expansion of two phase helium systems, resulting in
pressurization of the system.

###  1.6 Expansion ratio

The expansion of helium from its normal boiling point to standard temperature
and pressure is 700x. The expansion of helium may lead to over-pressure of the
system, exercising of pressure safety relief valves, and venting of helium
into ambient spaces. Due to the high volumetric expansion the venting helium
may displace large volumes of air, inert the ambient environment, and create
an oxygen deficient environment and potential asphyxiation.

##  2 Approaches

Several attempts have been made to analyze and model the loss of vacuum
effects on cryogenic helium systems, most commonly liquid helium systems at
normal boiling point 4.2 K and superfluid helium 2 K.

###  2.1 Bulk vessels

The simplest representation of the loss of vacuum problem is the bulk vessel
approach in which liquid helium is stored within a vacuum insulating vessel
and venting of the vacuum with air allows for the study of heat transfer and
propagative effects in bulk helium.

Perhaps the most well known and cited cryogenic helium loss of vacuum study is
by Lehman and Zahn, who tested the venting of air into the vacuum space of a
4.2 K liquid helium vacuum insulated vessel, with and without multi-layer
insulation (MLI) on the inner vessel [1]. Mass flow of air into the vacuum
space and venting of helium from the inner vessel was measured and used to
calculate the overall heat transfer from air to the helium. The results
concluded a peak heat flux without MLI up to 3.8 W/cm2 while the MLI wrapped
vessel limited peak heat flux to 0.6 W/cm2, thus indicating that use of MLI on
liquid helium cooled surfaces can significantly mitigate loss of vacuum
impacts.

A similar experiment was conducted by Harrison et al by venting air to the
insulating vacuum space on a bulk helium vessel both with and without
composite insulation, but instead with superfluid helium at 1.8 K [2].
Observed peak heat flux without and with insulation was 3.1 W/cm2 and 0.44
W/cm2, respectively, indicating again that insulating the internal vessel has
significant mitigative effects on heat transfer to liquid helium due to loss
of vacuum, and agreeing rather well with the Lehman and Zahn results.

Another set of researchers at KIT developed an analytical model for sudden
loss of insulating vacuum by air [3] and designed a liquid helium pressurized
cryostat test facility named PICARD to help validate the model [4] testing
loss of vacuum heat transfer of a normal boiling liquid helium volume of 100
liters, operating pressure of 16 barg, and a helium mass flow vent rates
$\leq$ 4 kg/s [5],[6]. Results from PICARD concluded that relief valves sized
for 4 W/cm2 tend to be oversized, resulting in chatter and possible
overpressure due to reduced relief capacity. Further results also concluded
that while peak heat flux can exceed 5 W/cm2, the nominal heat flux decays
quickly to $<$ 2 W/cm2 due to reduced heat transfer rates from gas
desublimation on the cold surfaces [7],[8].

###  2.2 Representative geometry

Due to the complex dynamics at work during a loss of insulating vacuum on a
liquid helium cryogenic system and fast temperature, pressure, and heat
transfer transients, attempts have been made to use more representative
geometry to test loss of vacuum effects for specific applications rather than
relying on results from bulk vessel tests.

Superconducting radio frequency (SRF) cavities for particular acceleration are
one of the most prominent applications for large scale superfluid helium usage
and over-pressure safety control tends to be strongly influenced by the very
high geometric aspect ratios of these devices (length $>>$ diameter).
Experiments on scaled SRF cavities in which the beamline vacuum is vented to
air have concluded that this high aspect ratio reduces the rate of heat
transfer to the superfluid by cryo-pumping the inflowing air on the cold pipe
inner surface, effectively slowing the speed of air in the beam-line to
$\sim$10 m/s, well below the air’s speed of sound, even at low temperatures
[9].

###  2.3 Simplified geometry

Instead of designing and developing experiment test geometry comparable to
production cryogenic system geometry to mitigate uncertainty in the transient
effects, another approach is to reduce the number of variables and study the
loss of insulating vacuum on simplified geometry.

One of the most systematic approaches to the study of loss of insulating
vacuum in normal boiling point liquid helium were performed by Bosque [10] and
Dhuley [11] in which the rate of air deposition during a loss of vacuum event
is measured across a liquid wetted copper disc. These tests concluded that the
air condensation rate quickly peaks as the vacuum is vented, but also decays
quickly as the total deposition of air steadily increases, while the rate of
heat transfer to the LHe bath is limited by film boiling. Dhuley and Van
Sciver [12],[13] continued these experiments with an evacuated vertical tube
immersed in LHe, and instrumented longitudinally with pressure and temperature
transmitters to accurately measure the air deposition propagation rate, and
calculated an effective heat flux $\sim$2.3 W/cm2 based on the pool boiling
heat transfer curve for 4.2 K saturated LHe [14]. The authors followed up
their experimental results with development of a numerical model based on
parametric data across multiple vacuum venting air flow rates to predict and
explain that the slow air propagation and deposition speed is a function of
the reduced heat transfer rate by desublimation of air on the helium cooled
tube surface and exponential decrease in air temperature along the length of
the channel [15],[16]. The model predicts that given a high aspect ratio (e.g.
in an evacuated tube or a in a SRF beamline) there is a finite length in which
the air ingress venting the evacuated space will freeze out entirely with only
minimal thermal impact to the rest of the cold tube. Further studies analyzing
the dynamic effects of the air deposition heat transfer concluded that while
the peak heat flux from air can exceed 10 W/cm2, the peak heat flux is
actually limited by film boiling in LHe to $\leq$ 2.5 W/cm2 [17]-[20].

Continuing the work of Dhuley et al, Garceau et al performed comparable
experiments with an extended length of evacuated tubes and performed tests
with superfluid helium [21]-[24]. Bao et al developed a representative
theoretical model of Garceau’s experimental results [25].

###  2.4 Field tests

While several options exist to model and simulate loss of insulating vacuum in
liquid helium systems, the most representative way to understand the transient
effects and heat transfer rate for a specific system design is to perform a
field test on a full-scale prototype of an as-built system.

For SRF cryo-module cavity applications, two such sets of loss of vacuum field
experiments have been performed at full-scale. Tests at CEBAF at Thomas
Jefferson National Laboratory in 1994 studied loss of beamline vacuum on a
quarter SRF cryomodule at both 2 K and 4 K LHe and measured peak heat flux of
3.5 W/cm2 and 2.8 W/cm2, respectively, with a maximum sustained heat flux of 2
W/cm2 [26]. Similarly, tests at XFEL at Deutsches Elektronen-Synchrotron
(DESY) in 2008 measured the heat flux to helium from loss of insulating vacuum
and loss of beam vacuum at 2 K and at 4.4 K with a full-scale, complete XFEL
SRF cryomodule [27]. The measured heat flux due to loss of insulating vacuum
at 2 K and 4.4 K was $\leq$ 0.4 W/cm2 and $\leq$ 0.65 W/cm2, respectively
while the loss of beam vacuum heat flux was $\leq$ 2.3 W/cm2. The experiments
also concluded substantial delay in air pressure propagation during venting
the beam vacuum of $\sim$3 m/s based on a $>$ 4 seconds delayed in pressure
pickup between pressure transmitters at the near and far end of $\sim$12 m SRF
cryomodule beam line.

##  3 Codes and standards

Loss of vacuum influences on system design typically manifests in overpressure
protection. Design codes and standards provide methodology for sizing relief
devices and vent lines to prevent over-pressure of piping and pressure
vessels. Notable codes include Compressed Gas Association S1-3 [28] and
International Standards Organization ISO-EN-21013 [29]. ISO-EN-21013 provides
guidance to estimate heat load to LHe for MLI and uninsulated vessels. CGA
S1-3 specifies calculating heat load based on the surface area of the helium
containing vessel, the thickness of thermal insulation, and thermal
conductivity of the insulation, which for vacuum insulated vessels is not
always appropriate. For dynamic thermofluid properties of helium, reference is
made to National Bureau of Standards Technology of Liquid Helium [30].

##  4 Conclusion

Loss of insulating vacuum, and loss of beam vacuum for SRF cavity
applications, is a reasonably plausible failure mode for liquid helium systems
which often is considered a worst-case failure. There have been significant
efforts to study and understand the dynamic heat flux and transient effects of
air deposition and propagative affects, especially for high aspect ratios
geometries. Generally for LHe systems, loss of insulating vacuum can be
predicted between $\sim$2-4 W/cm2 for uninsulated vessels and $\sim$0.6 W/cm2
for MLI insulated vessels, while loss of beam vacuum failures can be predicted
between $\sim$1.4-2.3 W/cm2 since beam lines are typically uninsulated.

## References

  * [1] W. Lehmann and G. Zahn, “Safety Aspects for the LHe Cryostats and LHe Transport Containers,” Proc. of the ICEC 7, London, IPC Science and Technology Press, pp. 569-579, 1978.
  * [2] S. M. Harrison, “Loss of vacuum experiments on a superfluid helium vessel,” IEEE Transactions on Applied Superconductivity 12(1), 1343 – 1346, 2002.
  * [3] C. Heidt, S. Grohmann, and M. Süßer, “Modeling the pressure increase in liquid helium cryostats after failure of the insulating vacuum,” AIP Conference Proceedings 1573 (1), 1574-1580, 2014.
  * [4] C. Heidt, H. Schön, M. Stamm, and S. Grohmann, “Commissioning of the cryogenic safety test facility PICARD,” IOP Conference Series: Materials Science and Engineering 101 (1), 012161, 2015.
  * [5] C. Heidt, A. Henriques, M. Stamm, and S. Grohmann, “First experimental data of the cryogenic safety test facility PICARD,” IOP Conference Series: Materials Science and Engineering 171, 012044, 2017.
  * [6] C. Heidt, “Experimental Investigation and Modelling of Incidents in Liquid Helium Cryostats,” PhD Dissertation, Karlsruhe Institute of Technology, 2018.
  * [7] C. Weber, A. Henriques, C. Zoller, and S. Grohmann, “Safety studies on vacuum insulated liquid helium cryostats,” IOP Conference Series: Materials Science and Engineering 278, 012169, 2017.
  * [8] C. Weber, A. Henriques, and S. Grohmann, “Study on the heat transfer of helium cryostats following loss of insulating vacuum,” IOP Conference Series: Materials Science and Engineering 502, 012170, 2019.
  * [9] Andrew A. Dalesandro, Ram C. Dhuley, Jay C. Theilacker, and Steven W. Van Sciver, “Results from sudden catastrophic loss of vacuum on scaled superconducting radio frequency cryomodule”, AIP Conference Proceedings 1573, 1822-1829, 2014.
  * [10] E. S. Bosque, R. C. Dhuley, and S. W. Van Sciver, “Transient heat transfer in Helium II due to sudden vacuum break”, AIP Conference Proceedings 1573, 260-267, 2014.
  * [11] R. C. Dhuley, E. S. Bosque, and S. W. Van Sciver, “Cryodeposition of nitrogen gas on a surface cooled by helium II,” AIP Conference Proceedings 1573, 626-632, 2014.
  * [12] Ram C. Dhuley and Steven W. Van Sciver, “Sudden vacuum loss in long liquid helium cooled tubes”, IEEE Transactions on Applied Superconductivity 25(3), 9000305, 2015.
  * [13] R. C. Dhuley and S. W. Van Sciver, “Epoxy encapsulation of a CernoxTM SD thermometer for measuring the temperature of surfaces in liquid helium”, Cryogenics 77, 49-52, 2016.
  * [14] W. G. Steward, “Transient helium heat transfer: Phase I– Static coolant,” International Journal of Heat and Mass Transfer 21, 863–874, 1978.
  * [15] R. C. Dhuley and S. W. Van Sciver, “Propagation of nitrogen gas in a liquid helium cooled vacuum tube following sudden vacuum loss- Part I: Experimental investigations and analytical modeling”, International Journal of Heat and Mass Transfer 96, 573-581, 2016.
  * [16] R. C. Dhuley and S. W. Van Sciver, “Propagation of nitrogen gas in a liquid helium cooled vacuum tube following sudden vacuum loss- Part II: Analysis of propagation speed”, International Journal of Heat and Mass Transfer 98, 728-737, 2016.
  * [17] R. C. Dhuley and S. W. Van Sciver, “Heat transfer in a liquid helium cooled vacuum tube following sudden vacuum loss”, IOP Conference Series: Material Science and Engineering 101, 012006, 2015.
  * [18] R. C. Dhuley and S. W. Van Sciver, “Nitrogen gas propagation in a liquid helium cooled vacuum tube following a sudden vacuum loss,” IOP Conference Series: Material Science and Engineering 171, 012084, 2015.
  * [19] R. C. Dhuley, S. Posen, O. Prokofiev, M. I. Geelhoed, and J. C. T. Thangaraj, “First demonstration of a cryocooler conduction-cooled superconducting radiofrequency cavity operating at practical cw accelerating gradients”, Superconductor Science and Technology 33, 06LT01, 2020.
  * [20] R. C. Dhuley, R. Kostin, O. Prokofiev, M. I. Geelhoed, T. H. Nicol, S. Posen, J. C. T. Thangaraj, T. K. Kroc, and R. D. Kephart, “Thermal link design for conduction cooling of a superconducting radio frequency cavity using cryocoolers”, IEEE Transactions on Applied Superconductivity 29(5), 0500205, 2019.
  * [21] N. Garceau, S. Bao, and W. Guo, “Effect of mass flow rate on gas propagation after vacuum break in a liquid helium cooled tube,” IOP Conference Series: Material Science and Engineering 755, 012112, 2015.
  * [22] N. Garceau, S. Bao, and W. Guo, “Heat and mass transfer during a sudden loss of vacuum in a liquid helium cooled tube – Part I: Interpretation of experimental observations,” International Journal of Heat and Mass Transfer 129, 1144-1150, 2019.
  * [23] N. Garceau, S. Bao, W. Guo, and S. W. Van Sciver, “The design and testing of a liquid helium cooled tube system for simulating sudden vacuum loss in particle accelerators,” Cryogenics 100, 92-96, 2019.
  * [24] N. Garceau, S. Bao, W. Guo, “Heat and mass transfer during a sudden loss of vacuum in a liquid helium cooled tube - Part III: Heat deposition in He II,” International Journal of Heat and Mass Transfer 181, 121885, 2021.
  * [25] S. Bao, N. Garceau, and W. Guo, “Heat and mass transfer during a sudden loss of vacuum in a liquid helium cooled tube – Part II: Theoretical modeling,” International Journal of Heat and Mass Transfer 146, 118883, 2020.
  * [26] M. Wiseman, K. Crawford, M. Drury, K. Jordan, J. Preble, Q. Saulter, and W. Schneider, “Loss of cavity vacuum experiment at CEBAF,” Advances in cryogenic engineering, Volume 39(A), 997-1003, 1994.
  * [27] T. Boekman et al., “Experimental tests of fault conditions during the cryogenic operation of a XFEL prototype cryomodule,” Proc. ICEC 22- ICMC 2008, 22:723–728, 2009.
  * [28] Compressed Gas Association, “S.1: Pressure Relief Device Standards”
  * [29] Cryogenic vessels — Pressure-relief accessories for cryogenic service — Part 3: Sizing and capacity determination.
  * [30] B. W. Birmingham and R. H. Kropschot, “Technology of liquid helium”, 1968.