# Gravitational corrections to the Einstein-Scalar-QCD model

Huan Souza<EMAIL_ADDRESS>Faculdade de Física, Universidade Federal
do Pará, 66075-110, Belém, Pará, Brazil.  L. Ibiapina Bevilaqua
<EMAIL_ADDRESS>Escola de Ciências e Tecnologia, Universidade
Federal do Rio Grande do Norte
Caixa Postal 1524, 59072-970, Natal, Rio Grande do Norte, Brazil.  A. C. Lehum
<EMAIL_ADDRESS>Faculdade de Física, Universidade Federal do Pará, 66075-110,
Belém, Pará, Brazil.

###### Abstract

This study employs the effective field theory approach to quantum gravity to
investigate a non-Abelian gauge theory involving scalar particles coupled to
gravity. The study demonstrates explicitly that the Slavnov-Taylor identities
are maintained at one-loop order, which indicates that the universality of the
color charge is preserved. Additionally, the graviton corrections to the two-
loop gluon self-energy and its renormalization are computed.

##  I Introduction

Although we are still in need of a consistent and generally accepted
description of quantum gravity at high energies, if we restrict ourselves to
low energies compared to the Planck scale, we can nevertheless draw some
trustful conclusions about the gravitational phenomena at quantum level using
the viewpoint and methods of effective field theories Donoghue:1994dn ;
Burgess:2003jk ; Shapiro . Thus, the well known nonrenormalizability of
Einstein’s theory coupled to other fields 'tHooft:1974bx ; PhysRevLett.32.245
; Deser:1974cy  is not an impediment to study the influence of gravity in the
renormalization of other fields and parameters in a meaningful way. The
central idea is that we add to the action the high-order terms needed to
renormalize the parameters of the lower-order terms and the new parameters
introduced will be irrelevant to the low-energy behavior of the theory.

As it is well known, the renormalized quantities of a theory depend on an
arbitrary scale and the renormalization group is the theoretical tool to study
this dependence and allows us to describe how the coupling constants change
with this scale, establishing the so-called running of the coupling constants
Srednicki:2007qs . If this dependence is such that the coupling constant gets
weaker as we go to higher energies the theory is said to be asymptotically
free Gross:1973id ; Politzer:1973fx ; Gross:1974jv . The possibility that
gravitational corrections could render all gauge coupling constants
asymptotically free was suggested by Robinson and Wilczek, who used the
effective field theory approach of quantum gravity to reach this conclusion
Robinson:2005fj . However, this result was soon contested by Pietrykowski
Pietrykowski:2006xy , who showed that the result was gauge dependent.
Subsequently, many works investigate the use of the renormalization group in
quantum gravity as an effective field theory (See for instance Refs.
Felipe:2012vq ; Felipe:2013vq ; Ebert:2007gf ; Nielsen:2012fm ; Toms:2008dq ;
Toms:2010vy ; Ellis:2010rw ; Anber:2010uj ; Bevilaqua:2015hma ;
Bevilaqua:2021uzk ; Bevilaqua:2021uev ). In a previous work Bevilaqua:2015hma
, we used dimensional regularization to compute gravitational effects on the
beta function of the scalar quantum electrodynamics at one-loop order and
found that all gravitational contributions cancel out. The situation is
different at two-loop order, in which we do find nonzero gravitational
corrections to the beta function for both scalar and fermionic QED, as shown
in a latter work Bevilaqua:2021uzk . However, those corrections give a
positive contribution to the beta function and thus the electrical charge is
not asymptotically free neither has a nontrivial fixed point.

The use of renormalization group in the context of non-renormalizable field
theories raise some subtle questions. The universality of the coupling
constants in effective field theories was discussed by Anber et al. in
Anber:2010uj , where it was suggested that an operator mixing could make the
coupling constants dependent on the process under consideration and therefore
non-universal. That would imply that, unlike renormalizable field theories,
the concept of running coupling may not be useful in the effective field
theory approach to quantum gravity. This is indeed the case for the quartic
self-interaction of scalars in scalar-QED, as discussed in Bevilaqua:2015hma
but, as shown in Bevilaqua:2015hma  for scalar-QED and in Bevilaqua:2021uev
for fermionic-QED it seems not to be the case for the gauge coupling because
of the Ward identity. The central role of the gauge symmetry in the
universality of the gauge coupling for QED led us to explore this issue in the
non-Abelian case. Using dimensional regularization, we showed that the
Slavnov-Taylor identities are satisfied in a non-Abelian gauge theory coupled
to fermions and gravity Souza:2022ovu . In the same work, we have also
calculated the gravitatinal correction for the beta function at one-loop thus
verifying directly the absence of contributions from the gravitational sector.

In previous studies, the coupling of non-Abelian gauge theories to gravity has
been investigated Souza:2022ovu ; Buchbinder:1983nug ; Tang:2008ah ;
Tang:2011gz . In this research, we extend our previous analysis by
investigating the asymptotic behavior of a non-Abelian gauge theory coupled to
complex scalars and gravity. This exploration is motivated by the significant
role scalar theories play in the advancement of high-energy theory. Over the
years, scalar models have been proposed to tackle issues such as
renormalization group theory for non-renormalizable theories Barvinsky:1993zg
, the study of dilatons Shapiro:1995yc , and potential candidates for dark
matter Cohen:2011ec ; Arkani-Hamed:2008hhe .. In fact, Ref. Calmet:2021iid
argue that quantum gravity might have crucial implications in a theory of dark
matter. Additionally, a recent study and:2022ttn  investigated the interaction
between SU(2) Yang-Mills waves and gravitational waves. The results revealed
that while the problem can be perturbatively studied in the symmetric phase,
non-perturbative approaches are necessary in the broken phase. Hence, the
examination of a non-Abelian gauge theory coupled to complex scalars and
gravity is of particular interest due to the fundamental role scalar theories
have played in addressing diverse problems in high-energy theory.

The paper is structured as follows. Section II introduces the Lagrangian and
propagators of the model. In Section III, the one-loop renormalization of the
model is presented, highlighting the preservation of gauge invariance of the
gravitational interaction and respect for the Slavnov-Taylor identities.
Section IV utilizes the Tarasov algorithm to compute the two-loop counterterm
for the gluon wave-function. Finally, concluding remarks are provided in
Section V. The minimal subtraction (MS) scheme is used throughout this work to
handle the UV divergences, with $(+---)$ being the spacetime signature and
natural units of $c=\hbar=1$ are adopted.

##  II The Einstein-Scalar-QCD model

To get an effective field theory description for our model, we add higher
order terms to the Lagrangian of a non-Abelian gauge theory with complex
scalars coupled to gravity:

$\displaystyle\mathcal{L}=$
$\displaystyle\sqrt{-g}\sum_{f}\Big{\\{}\frac{2}{\kappa^{2}}R-\frac{1}{4}g^{\mu\alpha}g^{\nu\beta}G_{\mu\nu}^{a}G_{\alpha\beta}^{a}+g^{\mu\nu}(D_{\mu}\phi^{i})^{\dagger}D_{\nu}\phi^{i}-m_{i}(\phi^{i})^{\dagger}\phi^{i}+\lambda((\phi^{i})^{\dagger}\phi^{i})^{2}+\mathcal{L}_{HO}\Big{\\}},$
(1)

where the index $i=1,2,\cdots,N_{s}$ runs over the scalars flavors,
$G^{a}_{\mu\nu}=\nabla_{\mu}A_{\nu}^{a}-\nabla_{\nu}A_{\mu}^{a}+gf^{abc}A^{b}_{\mu}A^{c}_{\nu}$
is the non-Abelian field-strength with $f^{abc}$ being the structure constants
of the $SU(N)$ group, and $D_{\mu}=\partial_{\mu}-igt^{a}A^{a}_{\mu}$ is the
covariant derivative. The higher order terms $\mathcal{L}_{HO}$ are written as

$\mathcal{L}_{HO}=\frac{\tilde{\lambda}_{1}}{M_{P}^{2}}\left[\mathrm{Re}((\phi^{i})^{\dagger}\partial_{\mu}\phi^{i})\right]^{2}+\frac{\tilde{\lambda}_{2}}{M_{P}^{2}}\left[\mathrm{Im}((\phi^{i})^{\dagger}\partial_{\mu}\phi^{i})\right]^{2}-\frac{\tilde{e}_{3}}{4}G_{a}^{\mu\nu}\frac{\Box}{M_{P}^{2}}G^{a}_{\mu\nu}.$
(2)

To obtain the usual quadratic term for the gravitational field, we need to
expand $g_{\mu\nu}$ around the flat metric as

$g_{\mu\nu}=\eta_{\mu\nu}+\kappa h_{\mu\nu},$ (3)

such that

$g^{\mu\nu}=\eta^{\mu\nu}-\kappa
h^{\mu\nu}+\cdots\qquad\text{and}\qquad\sqrt{-g}=1+\frac{\kappa}{2}h+\cdots,$
(4)

where $h=\eta^{\mu\nu}h_{\mu\nu}$. The affine connection is written as

$\Gamma^{\lambda}_{\leavevmode\nobreak\
\mu\nu}=\frac{1}{2}\kappa(\eta^{\lambda\sigma}-\kappa
h^{\lambda\sigma})(\partial_{\mu}h_{\sigma\nu}+\partial_{\nu}h_{\sigma\mu}-\partial_{\sigma}h_{\mu\nu}).$
(5)

Organizing the Lagrangian as,

$\displaystyle\mathcal{L}$ $\displaystyle=$
$\displaystyle\mathcal{L}_{h}+\mathcal{L}_{f}+\mathcal{L}_{A};$ (6a)
$\displaystyle\mathcal{L}_{h}$ $\displaystyle=$
$\displaystyle\frac{2}{\kappa^{2}}\sqrt{-g}R;$ (6b)
$\displaystyle\mathcal{L}_{s}$ $\displaystyle=$
$\displaystyle\sqrt{-g}[g^{\mu\nu}(D_{\mu}\phi^{i})^{\dagger}D_{\nu}\phi^{i}-m_{i}(\phi^{i})^{\dagger}\phi^{i}+\lambda((\phi^{i})^{\dagger}\phi^{i})^{2}];$
(6c) $\displaystyle\mathcal{L}_{A}$ $\displaystyle=$
$\displaystyle-\frac{\sqrt{-g}}{4}g^{\mu\alpha}g^{\nu\beta}G_{\mu\nu}^{a}G^{a}_{\alpha\beta}.$
(6d)

Using Eqs. (3)-(5), we write the pure gravity sector (6b) in terms of
$h_{\mu\nu}$. Moreover, it is convinient to organize $\mathcal{L}_{h}$ in
powers of $h$ as follows:

$\displaystyle\mathcal{L}_{h}$ $\displaystyle=$
$\displaystyle\mathcal{L}_{h}^{0}+\kappa\mathcal{L}_{h}^{1}+\cdots$ (7a)
$\displaystyle\mathcal{L}_{h}^{0}$ $\displaystyle=$
$\displaystyle-\frac{1}{4}\partial_{\mu}h\partial^{\mu}h+\frac{1}{2}\partial_{\mu}h^{\sigma\nu}\partial^{\mu}h_{\sigma\nu};$
(7b) $\displaystyle\mathcal{L}_{h}^{1}$ $\displaystyle=$
$\displaystyle\frac{1}{2}h^{\alpha}_{\leavevmode\nobreak\
\beta}\partial^{\mu}h^{\beta}_{\leavevmode\nobreak\
\alpha}\partial_{\mu}h-\frac{1}{2}h^{\alpha}_{\leavevmode\nobreak\
\beta}\partial_{\alpha}h^{\mu}_{\leavevmode\nobreak\
\nu}\partial^{\beta}h^{\nu}_{\leavevmode\nobreak\
\mu}-h^{\alpha}_{\leavevmode\nobreak\
\beta}\partial_{\mu}h^{\nu}_{\leavevmode\nobreak\
\alpha}\partial^{\mu}h^{\beta}_{\leavevmode\nobreak\ \nu}$ (7c)
$\displaystyle+\frac{1}{4}h\partial^{\beta}h^{\mu}_{\leavevmode\nobreak\
\nu}\partial_{\beta}h^{\nu}_{\leavevmode\nobreak\
\mu}+h^{\beta}_{\leavevmode\nobreak\
\mu}\partial_{\nu}h^{\alpha}_{\leavevmode\nobreak\
\beta}\partial^{\mu}h^{\nu}_{\leavevmode\nobreak\
\alpha}-\frac{1}{8}h\partial^{\nu}h\partial_{\nu}h,$

where the indices are raised and lowered with the flat metric (here and
henceforth, we are following the results in Ref. Choi:1994ax ).

For the matter sector (6c), the expansion around the flat metric give us

$\displaystyle\mathcal{L}_{s}$ $\displaystyle=$
$\displaystyle(D^{\mu}\phi^{i})^{\dagger}D_{\mu}\phi^{i}-m^{2}_{i}((\phi^{i})^{\dagger}\phi^{i})-\frac{\lambda}{4}((\phi^{i})^{\dagger}\phi^{i})^{2}-\kappa
h^{\mu\nu}(D_{\mu}\phi^{i})^{\dagger}D_{\nu}\phi^{i}$ (8a)
$\displaystyle+\frac{\kappa}{2}h\left[(D^{\mu}\phi^{i})^{\dagger}D_{\mu}\phi^{i}-m^{2}_{i}(\phi^{i})^{\dagger}\phi^{i}-\frac{\lambda}{4}((\phi^{i})^{\dagger}\phi^{i})^{2}\right],$

which we organize as follows

$\displaystyle\mathcal{L}_{s}$ $\displaystyle=$
$\displaystyle\mathcal{L}_{s}^{0}+\kappa\mathcal{L}_{s}^{1}+\cdots$ (9a)
$\displaystyle\mathcal{L}_{s}^{0}$ $\displaystyle=$
$\displaystyle(D^{\mu}\phi^{i})^{\dagger}D_{\mu}\phi^{i}-m_{i}^{2}((\phi^{i})^{\dagger}\phi^{i})-\frac{\lambda}{4}((\phi^{i})^{\dagger}\phi^{i})^{2}$
(9b) $\displaystyle\mathcal{L}_{s}^{1}$ $\displaystyle=$
$\displaystyle-h^{\mu\nu}(D_{\mu}\phi^{i})^{\dagger}D_{\nu}\phi^{i}+\frac{1}{2}h\left[(D^{\mu}\phi^{i})^{\dagger}D_{\mu}\phi^{i}-m_{i}^{2}(\phi^{i})^{\dagger}\phi^{i}-\frac{\lambda}{4}((\phi^{i})^{\dagger}\phi^{i})^{2}\right];$
(9c)

and finally, for the gauge sector,

$\displaystyle\mathcal{L}_{A}$ $\displaystyle=$
$\displaystyle\mathcal{L}_{A}^{0}+\kappa\mathcal{L}_{A}^{1}+\cdots$ (10a)
$\displaystyle\mathcal{L}_{A}^{0}$ $\displaystyle=$
$\displaystyle-\frac{1}{4}G_{\mu\nu}^{a}G^{\mu\nu}_{a}$ (10b)
$\displaystyle\mathcal{L}_{A}^{1}$ $\displaystyle=$
$\displaystyle\frac{1}{2}h^{\tau}_{\leavevmode\nobreak\
\nu}G^{\mu\nu}_{a}G_{\mu\tau}^{a}+\frac{1}{2}h\mathcal{L}_{A}^{0}.$ (10c)

As usual for gauge theories, in order to quantize this model, we have to deal
with the excess of degrees of freedom in $A_{\mu}^{a}$ and $h_{\mu\nu}$ due to
their symmetries. In our calculations, we have followed the Faddeev-Popov
procedure that introduces gauge-fixing terms in the action that will modify
the propagators of both $A_{\mu}^{a}$ and $h_{\mu\nu}$. Moreover, we must also
introduce ghost fields for both vector and tensor fields. However, the ghost
field associated with the graviton will not appear in this text because, since
we are working with the one-graviton exchange approximation, the new term
containing the ghosts added to the action will not contribute to the
renormalization of the gauge coupling constant. Therefore, whenever we refer
to ghost field in what follows, we mean the one associated with $A_{\mu}^{a}$.
The propagators for scalars, ghosts, gluons and gravitons are given,
respectively, by

$\displaystyle\Delta_{s}(p)$ $\displaystyle=$
$\displaystyle\frac{i}{p^{2}-m_{a}^{2}};$ (11a) $\displaystyle\Delta_{ab}(p)$
$\displaystyle=$ $\displaystyle\frac{i}{p^{2}}\delta_{ab};$ (11b)
$\displaystyle\Delta^{\mu\nu}_{ab}(p)$ $\displaystyle=$
$\displaystyle\frac{i}{p^{2}}\left(\eta^{\mu\nu}-(1-\xi_{A})\frac{p^{\mu}p^{\nu}}{p^{2}}\right)\delta_{ab};$
(11c) $\displaystyle\Delta^{\alpha\beta\mu\nu}(p)$ $\displaystyle=$
$\displaystyle\frac{i}{p^{2}}\left(P^{\alpha\beta\mu\nu}-(1-\xi_{h})\frac{Q^{\alpha\beta\mu\nu}}{p^{2}}\right).$
(11d)

The gauge-fixing parameters $\xi_{A}$ and $\xi_{h}$ will be carried out
through the whole calculation, since we do not want to choose any specific
gauge. The projectors $P^{\alpha\beta\mu\nu}$ and $Q^{\alpha\beta\mu\nu}$ in
the graviton propagator are given by

$\displaystyle P^{\alpha\beta\mu\nu}$ $\displaystyle=$
$\displaystyle\frac{1}{2}\left(\eta^{\alpha\mu}\eta^{\beta\nu}+\eta^{\alpha\nu}\eta^{\beta\mu}-\eta^{\alpha\beta}\eta^{\mu\nu}\right);$
$\displaystyle Q^{\alpha\beta\mu\nu}$ $\displaystyle=$
$\displaystyle(\eta^{\alpha\mu}p^{\beta}p^{\nu}+\eta^{\alpha\nu}p^{\beta}p^{\mu}+\eta^{\beta\mu}p^{\alpha}p^{\nu}+\eta^{\beta\nu}p^{\alpha}p^{\mu}).$
(12)

##  III The one-loop renormalization

The Slavnov-Taylor identities are a set of relations that must be satisfied by
the n-point functions to ensure the gauge independence of the observables of
the theory. In this section we want to explicitly show that the Slavnov-Taylor
identities are respected at one-loop order for our model. To simplify our
computations, we will consider here that all the masses are the same, so we
drop the index $i$. As we will see, this will not affect our final result.

We start by computing the n-point functions. Namely, the self-energy of
scalar, vector and ghost fields ($\Sigma_{s},\Pi^{\mu\nu}_{ab}$ and
$\Sigma_{ab}$, respectively), also the scalar-gluon, ghost-gluon and gluon-
gluon three-point functions ($\Gamma^{\mu}_{a}$, $\Gamma^{\mu}_{abc}$ and
$\Pi^{\mu\nu\alpha}_{abc}$, respectively), the gluon four-point function
($\Gamma^{\mu\nu\rho\sigma}_{abcd}$), and finally the scalar-gluon four-point
function ($\Pi^{\mu\nu}_{abcd}$). All the computations were done using the
Mathematica packages: FeynRules to generate the models feynrules , FeynArts to
draw the diagrams Hahn:2000kx , and FeynCalc to simplify and compute the
amplitudes Shtabovenko:2020gxv .

At one-loop, the self-energy of the scalar field, Fig. 1, results in

$\displaystyle-i\Sigma_{s}(p)$
$\displaystyle=ip^{2}\left(\frac{C_{A}\left(\xi_{A}-3\right)g^{2}-(\xi_{h}-2)\kappa^{2}m^{2}}{16\pi^{2}\epsilon}+Z_{2s}^{(1)}\right)$
(13) $\displaystyle+im^{2}\left(\frac{-C_{A}\xi_{A}g^{2}+4\lambda
N_{s}-(\xi_{h}-2)\kappa^{2}m^{2}}{16\pi^{2}\epsilon}-Z_{m_{s}}^{(1)}\right)+\mathrm{finite},$

where $C_{A}=N$ for the $SU(N)$ group. By imposing finiteness to
$\Sigma_{s}(p)$, we find the following one-loop counterterms:

$\displaystyle Z_{2s}^{(1)}$ $\displaystyle=$
$\displaystyle\frac{\kappa^{2}m^{2}\left(\xi_{h}-2\right)-C_{A}\left(\xi_{A}-3\right)g^{2}}{16\pi^{2}\epsilon},$
(14a) $\displaystyle Z_{m}^{(1)}$ $\displaystyle=$
$\displaystyle\frac{-C_{A}\xi_{A}g^{2}+4\lambda
N_{s}-(\xi_{h}-2)\kappa^{2}m^{2}}{16\pi^{2}\epsilon}.$ (14b)

Figure 1: Feynman diagrams for the scalar self-energy. Continuous, wiggly,
dotted, and dashed lines represent the scalar, gluon, ghost, and graviton
propagators, respectively.

For the gluon self-energy, it is convenient to write the one-loop correction
(corresponding to the diagrams in Fig. 2) as

$\displaystyle\Pi^{\mu\nu}_{ab}(p)=\left(p^{2}\eta^{\mu\nu}-p^{\mu}p^{\nu}\right)\Pi(p)\delta_{ab},$
(15)

where the function $\Pi(p)$ is found to be

$\displaystyle\Pi(p)=-iZ_{3}^{(1)}-ip^{2}\tilde{Z}_{3}^{(1)}+\frac{i\kappa^{2}p^{2}\left(2-3\xi_{h}\right)}{96\pi^{2}\epsilon}-\frac{iC_{A}g^{2}\left(2N_{s}+3\xi_{A}-13\right)}{96\pi^{2}\epsilon}+\mathrm{finite},$
(16)

and, imposing the finiteness on $\Pi(p)$, we find

$\displaystyle Z_{3}^{(1)}$ $\displaystyle=$
$\displaystyle-\frac{C_{A}g^{2}\left(2N_{s}+3\xi_{A}-13\right)}{96\pi^{2}\epsilon},$
(17a) $\displaystyle\tilde{Z}_{3}^{(1)}$ $\displaystyle=$
$\displaystyle-\frac{\kappa^{2}(3\xi_{h}-2)}{96\pi^{2}\epsilon}\leavevmode\nobreak\
.$ (17b)

We can see from Eq. (16) that $Z_{3}$ is the relevant counterterm to the beta
function of the color charge, since it is the renormalizing factor for the
quadratic term $G^{\mu\nu}_{a}G_{\mu\nu}^{a}$, while $\tilde{Z}_{3}$
renormalizes a higher derivative term like $G^{\mu\nu}_{a}\Box
G_{\mu\nu}^{a}$. Notice also that the UV divergent part of Eq. (16) is not
dependent on the masses of the scalars.

Figure 2: Feynman diagrams for the gluon self-energy. Figure 3: Feynman
diagrams for the ghost self-energy.

Contributions to the ghost self-energy up to one-loop order are depicted in
Fig. 3. The resulting expression is

$\displaystyle-i\Sigma_{ab}$ $\displaystyle=$
$\displaystyle\left(\frac{ip^{2}C_{A}\left(\xi_{A}-3\right)g^{2}}{64\pi^{2}\epsilon}+ip^{2}Z_{2c}^{(1)}\right)\delta_{ab}+\mathrm{finite},$
(18)

and, imposing finiteness, we find

$\displaystyle Z_{2_{c}}^{(1)}$ $\displaystyle=$
$\displaystyle-\frac{C_{A}g^{2}\left(\xi_{A}-3\right)}{64\pi^{2}\epsilon}.$
(19)

Notice that in Fig. 3 the gravitational interactions are not shown. Although
in the action there is a coupling of $h^{\mu\nu}$ to the kinetic term of the
ghosts associated with the gluons, the gravitational contributions to the
ghost self-energy will be renormalized by a higher-order term and is therefore
irrelevant for our purposes here. One way to see why this is happens is to
observe that both the ghosts and the graviton are massless, so the only
contribution proportional to $\kappa^{2}$ must be of the order $p^{4}$.

Figure 4: Feynman diagrams for the vertex interaction between gluons and
ghosts up to one-loop order.

For the 3-point functions, let’s first consider the ghost-ghost-gluon vertex
(Fig. 4), where again all the gravitational corrections are renormalized by
higher-order terms and are therefore omitted here. Also, in the following
expressions, we will use $p_{1}$ and $p_{2}$ to represent incoming external
momenta, and $p_{3}$ and $p_{4}$ for outgoing momenta. The expression obtained
for these diagrams is

$\Gamma_{abc}^{\mu}=-gp_{3}^{\mu}f_{abc}\left(\frac{C_{A}g^{2}\xi_{A}}{32\pi^{2}\epsilon}+Z_{1c}^{(1)}\right)+\text{finite},$
(20)

and the subtraction of the UV pole will give us

$Z_{1_{c}}^{(1)}=-\frac{C_{A}g^{2}\xi_{A}}{32\pi^{2}\epsilon}.$ (21)

Figure 5: Feynman diagrams to the vertex interaction between quarks top and
gluons up to one-loop order.

For the other 3-point function, the scalar-scalar-gluon vertex, the
gravitational interaction will be present in some diagrams, as we can see in
Fig. 5, where the relevant contributions to this function up to one-loop order
are shown. The resulting expression is

$\displaystyle-i\Gamma^{\mu}_{abc}$ $\displaystyle=$ $\displaystyle
gf_{abc}(p^{\mu}_{2}-p^{\mu}_{3})\left(\frac{C_{A}\left(9-5\xi_{A}\right)g^{2}+4\kappa^{2}m^{2}\left(\xi_{h}-2\right)}{64\pi^{2}\epsilon}-Z_{1}^{(1)}\right)$
(22) $\displaystyle+O(p^{3})+\text{finite},$

from which, through MS, we find

$Z_{1}^{(1)}=\frac{C_{A}\left(9-5\xi_{A}\right)g^{2}+4\kappa^{2}m^{2}\left(\xi_{h}-2\right)}{64\pi^{2}\epsilon}.$
(23)

The 3-point function describing the vertex with three gluons in shown in Fig.
6. We have used the projection

Figure 6: Feynman diagrams to the gluons vertex interaction at one-loop order.

$\Pi^{\mu\nu\alpha}_{abc}=\eta^{\mu\nu}\Pi^{\alpha}_{abc}\qquad\Rightarrow\qquad\Pi^{\alpha}_{abc}=\frac{1}{4}\eta_{\mu\nu}\Pi^{\mu\nu\alpha}_{abc}$
(24)

and used the fact that $p_{3}=p_{1}+p_{2}$, to get

$\displaystyle-i\Pi^{\alpha}_{abc}$ $\displaystyle=$
$\displaystyle\frac{g^{3}f_{abc}C_{A}\left(-9\xi_{A}-4N_{s}+17\right)(p_{1}-p_{2})^{\alpha}}{256\pi^{2}\epsilon}-\frac{3}{4}Z_{3g}^{(1)}g\left(p_{1}-p_{2}\right){}^{\alpha}f_{abc}$
(25) $\displaystyle+O(p^{2})+\mathrm{finite},$

Through MS, we impose finiteness and find

$Z_{3g}^{(1)}=-\frac{g^{2}C_{A}\left(9\xi_{A}-17-4N_{s}\right)}{192\pi^{2}\epsilon}.$
(26)

Now, we consider the scattering of four gluons (Fig. 7 showed at the end of
the paper for convenience). Since the interaction of four gluons has no
derivatives, the $Z_{4g}$ counterterm will renormalize terms proportional to
$p^{0}$ and therefore we can set external momentum equals to zero if we
restrict ourselves to the computation of this counterterm. Also, for
simplicity, we have used the scalar projection

$\Gamma_{abcd}=\frac{1}{16}\eta_{\mu\nu}\eta_{\rho\sigma}\Gamma^{\mu\nu\rho\sigma}_{abcd},$
(27)

to obtain the expression for the gluon 4-point function

$\displaystyle-i\Gamma_{abcd}$ $\displaystyle=$
$\displaystyle-\left(\frac{iC_{A}g^{4}\left(N_{s}+3\xi_{A}-2\right)}{32\pi^{2}\epsilon}+\frac{3}{2}iZ_{4g}^{(1)}g^{2}\right)\Bigr{(}\text{tr}(t_{a}t_{b}t_{c}t_{d})-2\text{tr}(t_{a}t_{c}t_{b}t_{d})-2\text{tr}(t_{b}t_{c}t_{a}t_{d})$
(28)
$\displaystyle+\text{tr}(t_{b}t_{a}t_{c}t_{d})+\text{tr}(t_{c}t_{a}t_{b}t_{d})+\text{tr}(t_{c}t_{b}t_{a}t_{d})\Bigr{)},$

Then, again imposing finiteness through MS, we have

$Z^{(1)}_{1_{4g}}=-\frac{C_{A}g^{2}\left(N_{s}+3\xi_{A}-2\right)}{48\pi^{2}\epsilon}.$
(29)

The other 4-point function involves two scalars and two gluons (Fig. 8, again
showed at the end of the paper for convenience). For this vertex, we use the
following projection

$\Pi^{\mu\nu}_{abcd}=\eta^{\mu\nu}\Pi_{abcd}\qquad\Rightarrow\qquad\Pi_{abcd}=\frac{1}{4}\eta_{\mu\nu}\Pi^{\mu\nu}_{abcd}$
(30)

and then we have

$\displaystyle\Pi_{abcd}$
$\displaystyle=\left(\frac{ig^{2}-3C_{A}\left(\xi_{A}-1\right)g^{2}-2(\xi_{h}-2)\kappa^{2}m^{2}}{16\pi^{2}\epsilon}-2iZ_{2g}^{(1)}g^{2}\right)\Bigr{(}2\text{tr}(t_{a}t_{b}t_{c}t_{d})-\text{tr}(t_{a}t_{c}t_{b}t_{d})$
(31)
$\displaystyle-\text{tr}(t_{b}t_{a}t_{c}t_{d})-\text{tr}(t_{b}t_{c}t_{a}t_{d})-\text{tr}(t_{c}t_{a}t_{b}t_{d})+2\text{tr}(t_{c}t_{b}t_{a}t_{d})\Bigr{)}.$

and the counterterm is found to be

$Z_{2g}^{(1)}=-\frac{3C_{A}\left(\xi_{A}-1\right)g^{2}-2(\xi_{h}-2)\kappa^{2}m^{2}}{32\pi^{2}\epsilon}.$
(32)

From Eqs. (14a), (17a), (19), (21), (23), (26), (29) we conclude that

$Z_{1}^{(1)}-Z_{2s}^{(1)}=Z_{3g}^{(1)}-Z_{3}^{(1)}=\frac{1}{2}\left(Z_{4g}^{(1)}-Z_{3}^{(1)}\right)=\frac{1}{2}\left(Z_{2g}^{(1)}-Z_{2s}^{(1)}\right)=Z_{1c}^{(1)}-Z_{2c}^{(1)}=-\frac{C_{A}g^{2}(3+\xi_{A})}{64\pi^{2}\epsilon}$
(33)

so the Slavnov-Taylor identities Slavnov:1972fg ; Taylor:1971ff  are indeed
respected and thus gravitational interaction does not spoil the gauge
symmetry. This result allows us to define a global color charge.

Moreover, we can show that the beta function is independent of $\kappa$ and
$m$, as the expression the one-loop beta function of the color charge can be
found through the relations between the renormalized coupling constants and
the counterterms given by

$\displaystyle g$ $\displaystyle=$
$\displaystyle\mu^{-2\epsilon}\frac{Z_{2s}Z_{3}^{1/2}}{Z_{1}}g_{0};$ (34a)
$\displaystyle g$ $\displaystyle=$
$\displaystyle\mu^{-2\epsilon}\frac{Z_{3}^{3/2}}{Z_{3g}}g_{0};$ (34b)
$\displaystyle g$ $\displaystyle=$
$\displaystyle\mu^{-2\epsilon}\frac{Z_{3}}{Z_{4g}^{1/2}}g_{0};$ (34c)
$\displaystyle g$ $\displaystyle=$
$\displaystyle\mu^{-2\epsilon}\frac{Z_{2c}Z_{3}^{1/2}}{Z_{1c}}g_{0};$ (34d)
$\displaystyle g$ $\displaystyle=$
$\displaystyle\mu^{-2\epsilon}\frac{Z_{2}^{1/2}Z_{3}^{1/2}}{Z_{2g}^{1/2}}g_{0}.$
(34e)

Therefore, the beta function for the color charge is

$\displaystyle\beta(g)$ $\displaystyle=$
$\displaystyle\lim_{\epsilon\rightarrow
0}\mu\frac{dg}{d\mu}=\lim_{\epsilon\rightarrow
0}\mu\frac{d}{d\mu}\left[g_{0}\left(1-Z_{1}^{(1)}+Z_{2s}^{(1)}+\frac{Z_{3}^{(1)}}{2}\right)\mu^{-2\epsilon}\right]$
(35) $\displaystyle=$
$\displaystyle-\frac{g^{3}}{(4\pi)^{2}}\left(\frac{11}{3}C_{A}-\frac{2}{6}N_{s}\right).$

The observed outcome is gauge-independent, a characteristic that was
previously established via a functional approach in Ref.Folkerts:2011jz . This
property has also been verified in the context of the Effective Field Theory
of gravity when coupled with fermionic QCD in Souza:2022ovu .

As we can see, it does not depend on the mass, so our choice to make all
masses the same does not affect our result for the beta function at one-loop
order. On the other hand, as discussed in Bevilaqua:2021uev , at two-loop we
would expect a $\sum_{i}\kappa^{2}m_{i}^{2}$ term.

It is needed to stress here the importance of a regularization scheme that
preserves the symmetries of the model. In fact, the authors in
Ref.Folkerts:2011jz  showed that in the weak-gravity limit there is no
gravitational contribution at one-loop order if the regularization scheme
preserves the symmetries of the model, such as dimensional regularization. On
the other hand, if the regularization scheme does not preserve all the
symmetries, there will be a negative contribution to the beta function (as
seen in Robinson:2005fj ).

##  IV Two-loop Gluon self-energy

This section presents the computation of the two-loop gluon self-energy and
its renormalization. TARCER Mertig:1998vk , in combination with previously
cited Mathematica packages, is utilized for this computation. TARCER
implements the Tarasov algorithm for the reduction of two-loop scalar
propagator type integrals with external momentum and arbitrary masses
Tarasov:1997kx . The Feynman and harmonic gauges ($\xi_{A}=\xi_{h}=1$) are
used for simplicity, and the analysis is limited to the case in which there is
only one scalar particle ($N_{s}=1$).

The Feynman diagrams we need to compute are showed in Fig. 9. Due to gauge
invariance, our result can be expressed as

$\Pi^{(2)}_{\mu\nu}=\left(p^{2}g_{\mu\nu}-p_{\mu}p_{\nu}\right)\Pi^{(2)},$
(36)

where the function $\Pi^{(2)}$ is a scalar function that can be expressed in
terms of a set of basic integrals. To present the results in a simplified
manner, we will adopt a notation similar to the one used in the original
TARCER paper Mertig:1998vk  for the basic integrals that will be utilized,

$\displaystyle\textbf{A}_{\nu}(m)=\frac{1}{\pi^{D/2}}\int\frac{d^{D}k}{[k^{2}-m^{2}]^{\nu}}$
(37a)
$\displaystyle\textbf{B}_{\nu_{1},\nu_{2}}(m_{1},m_{2})=\frac{1}{\pi^{D/2}}\int\frac{d^{D}k}{[k^{2}-m_{1}^{2}]^{\nu_{1}}[(k-p)^{2}-m_{2}^{2}]^{\nu_{2}}}$
(37b)
$\displaystyle\textbf{J}_{\nu_{1},\nu_{2},\nu_{3}}(m_{1},m_{2},m_{3})=\frac{1}{\pi^{D}}\int\frac{d^{D}k_{1}d^{D}k_{2}}{[k_{1}^{2}-m_{1}^{2}]^{\nu_{1}}[k_{5}^{2}-m_{2}^{2}]^{\nu_{2}}[k_{4}^{2}-m_{3}^{2}]^{\nu_{3}}}$
(37c)
$\displaystyle\textbf{F}_{\nu_{1},...,\nu_{5}}(m_{1},...,m_{5})=\frac{1}{\pi^{D}}\int\frac{d^{D}k_{1}d^{D}k_{2}}{[k_{1}^{2}-m_{1}^{2}]^{\nu_{1}}[k_{2}^{2}-m_{2}^{2}]^{\nu_{2}}[k_{3}^{2}-m_{3}^{2}]^{\nu_{3}}[k_{4}^{2}-m_{4}^{2}]^{\nu_{4}}[k_{5}^{2}-m_{5}^{2}]^{\nu_{5}}},$
(37d)

in which $p$ is the external momentum and we introduced $k_{3}=k_{1}-p$,
$k_{4}=k_{2}-p$, and $k_{5}=k_{1}-k_{2}$.

Therefore, we can write

$\displaystyle\Pi^{(2)}$ $\displaystyle=$ $\displaystyle
c_{1}\leavevmode\nobreak\ \textbf{A}_{1}(m)\leavevmode\nobreak\
\textbf{B}_{1,1}(0,0)+c_{2}\leavevmode\nobreak\
\textbf{A}_{1}(m)\leavevmode\nobreak\
\textbf{B}_{1,1}(m,m)+c_{3}\leavevmode\nobreak\
\textbf{B}_{1,1}(0,0)\leavevmode\nobreak\
\textbf{B}_{1,1}(m,m)+c_{4}\left(\textbf{A}_{1}(m)\right)^{2}$ (38)
$\displaystyle
c_{5}\left(\textbf{B}_{1,1}(0,0)\right)^{2}+c_{6}\left(\textbf{B}_{1,1}(m,m)\right)^{2}+c_{7}\leavevmode\nobreak\
\textbf{J}_{1,1,1}(0,0,0)+c_{8}\leavevmode\nobreak\
\textbf{J}_{1,1,1}(m,m,0)+c_{9}\leavevmode\nobreak\ \textbf{J}_{2,1,1}(m,m,0)$
$\displaystyle c_{10}\leavevmode\nobreak\
\textbf{F}_{1,1,1,1,1}(0,m,0,m,m)+c_{11}\leavevmode\nobreak\
\textbf{F}_{1,1,1,1,1}(m,0,m,0,m).$

All of the aforementioned integrals are established and can be found in
Refs.Martin:2005qm ; Martin:2003qz , and the coefficients $c_{i}$ are
presented in appendix A. As we are only concerned with the renormalization of
the gluon wave-function, we expand Eq.(38) around $p=0$ and retain only terms
proportional to $p^{0}$. Higher powers in the external momentum will be
renormalized by higher-order terms. Thus, we obtain:

$\displaystyle\Pi^{(2)}$ $\displaystyle=$ $\displaystyle-\frac{i\lambda
C_{A}\leavevmode\nobreak\
g^{2}}{384\pi^{4}\epsilon}-\frac{i\kappa^{2}m^{2}C_{A}\leavevmode\nobreak\
g^{2}}{256\pi^{4}\epsilon}+\frac{iC_{A}^{2}\leavevmode\nobreak\
g^{4}\log\left(m^{2}\right)}{384\pi^{4}\epsilon}-\frac{iC_{A}^{2}\leavevmode\nobreak\
g^{4}\log\left(-p^{2}\right)}{64\pi^{4}\epsilon}-\frac{i\lambda
C_{A}\leavevmode\nobreak\ g^{2}}{384\pi^{4}\epsilon}+\frac{5i\gamma
C_{A}^{2}\leavevmode\nobreak\ g^{4}}{384\pi^{4}\epsilon}$ (39)
$\displaystyle+\frac{17iC_{A}^{2}\leavevmode\nobreak\
g^{4}}{576\pi^{4}\epsilon}+\frac{5i\log(4\pi)C_{A}^{2}\leavevmode\nobreak\
g^{4}}{384\pi^{4}\epsilon}+\frac{5iC_{A}^{2}\leavevmode\nobreak\
g^{4}}{768\pi^{4}\epsilon^{2}}+O(p)+\mathrm{finite}.$

Now, we should compute the 1-loop diagrams with counterterms insertion in Fig.
10. By doing so, we obtain

$\displaystyle\Pi^{(2)}_{\mu\nu
CT}=(p^{2}g_{\mu\nu}-p_{\mu}p_{\nu})\Pi^{(2)}_{CT},$ (40)

where

$\displaystyle\Pi^{(2)}_{CT}$ $\displaystyle=$
$\displaystyle-\frac{iC_{A}^{2}\leavevmode\nobreak\
g^{4}\log\left(m^{2}\right)}{384\pi^{4}\epsilon}+\frac{iC_{A}^{2}\leavevmode\nobreak\
g^{4}\log\left(-p^{2}\right)}{64\pi^{4}\epsilon}-\frac{5iC_{A}^{2}\leavevmode\nobreak\
g^{4}}{384\pi^{4}\epsilon^{2}}+\frac{i\lambda C_{A}\leavevmode\nobreak\
g^{2}}{192\pi^{4}\epsilon}-\frac{5i\gamma C_{A}^{2}\leavevmode\nobreak\
g^{4}}{384\pi^{4}\epsilon}$ (41)
$\displaystyle-\frac{59iC_{A}^{2}\leavevmode\nobreak\
g^{4}}{2304\pi^{4}\epsilon}-\frac{5i\log(4\pi)C_{A}^{2}\leavevmode\nobreak\
g^{4}}{384\pi^{4}\epsilon}+O(p)+\mathrm{finite}.$

Therefore, we obtain that the two-loop gluon wave-function counterterm is
given by

$Z_{3}^{(2)}=\frac{C_{A}^{2}g^{4}}{256\pi^{4}\epsilon}-\frac{5C_{A}^{2}g^{4}}{768\pi^{4}\epsilon^{2}}-\frac{\kappa^{2}m^{2}C_{A}g^{2}}{256\pi^{4}\epsilon}.$
(42)

##  V Concluding remarks

In summary, we have evaluated the n-point functions for the Einstein-Scalar-
QCD model and demonstrated that there are no gravitational corrections to the
beta function of the color charge at one-loop order. Additionally, we have
explicitly verified that the Slavnov-Taylor identities are preserved at this
order of perturbation theory, indicating that the universality of the color
charge is maintained. Lastly, we have computed the counterterm for the gluon
wave-function at two-loop order.

It is important to contextualize our results and compare them with previous
research. To this end, we will follow the discussion in Donoghue:2019clr  and
highlight some distinctions between our findings and theirs. One such
difference lies in the adoption of a distinct regularization scheme. In
reference Tang:2008ah , it is argued that there are three primary concerns
that should be considered when working with quantum gravity: gauge invariance,
gauge conditions introduced in the quantization process, and the ability of
the method to regulate any type of divergence. It was further argued that
although dimensional regularization (DR) satisfies the first two requirements,
it cannot handle more than logarithmic divergences. Therefore, Tang and Wu
employed the Loop Regularization method (LP) in their studies Tang:2008ah ;
Tang:2011gz  to regulate the divergences. This method is capable of dealing
with the quadratic divergences that appear in the Feynman diagrams. The
authors used LP to compute the beta functions of the Einstein-Yang-Mills
theory and compared the results with those obtained using DR. They found that
while using DR leads to no gravitational contribution at one-loop, the use of
LP leads to a contribution that is proportional to $\mu^{2}$.

It is a fundamental requirement that physical results should not depend on the
choice of the regularization scheme. Anber pointed out in Anber:2010uj  that
the quadratic divergences are not relevant when using the S-matrix, which is a
physical quantity. Moreover, Toms demonstrated in Toms:2011zza  that it is
possible to define the electrical charge in quantum gravity using the
background field method in a physically meaningful way that is not influenced
by the quadratic divergences. Therefore, such contributions should be regarded
as unphysical and should not be included in the evaluation of the running
coupling.

An intriguing avenue for further investigation pertains to the existence of a
non-Abelian scalar particle serving as a potential dark matter candidate, as
well as the implications of quantum gravity for dark matter. In the study
conducted in Ref.Calmet:2021iid , the potential ramifications of quantum
gravity on dark matter models were explored. It was demonstrated that quantum
gravity would give rise to a fifth force-like interaction, setting a lower
limit on the masses of bosonic dark matter candidates. The authors also argued
that, due to the influence of quantum gravity, these potential candidates
would decay. However, given the ongoing observation of dark matter in the
present universe, the authors were able to calculate an upper bound on the
mass of a scalar singlet dark matter particle. In our future work, we intend
to investigate the mass range for a non-Abelian scalar dark matter candidate,
as presented in our study. In such a scenario, the fifth force-like
interaction would also be non-Abelian in nature. This particular scenario was
discussed in Arkani-Hamed:2008hhe .

In our future endeavors, we plan to investigate the dynamics of the
renormalized coupling constant in non-Abelian gauge theories, considering the
presence of fermions and scalars coupled to gravity at the two-loop level.
This investigation will involve an expansion of our research to incorporate
modified theories of gravity, such as quadratic gravity Odintsov:1991nd ;
Salvio:2014soa ; Donoghue:2018izj ; Donoghue:2021cza . Drawing on the
qualitative analysis presented in Souza:2022ovu , we expect that modified
theories of gravity, characterized by unconventional properties such as
repulsive gravity under specific regimes, could potentially impact the
behavior of the beta function. These modified gravity theories introduce
additional gravitational interactions and might influence the running of the
coupling constant in non-Abelian gauge theories, leading to intriguing and
novel phenomena.

###### Acknowledgements.

The work of HS is partially supported by Coordenação de Aperfeiçoamento de
Pessoal de Nível Superior (CAPES).

##  Appendix A Two-loop coefficients

In this section we present the two-loop coefficients for the two-loop gluon
self-energy from Eq. (38).

$\displaystyle c_{1}$ $\displaystyle=$
$\displaystyle-\frac{i\left(D^{4}-10D^{3}+35D^{2}-50D+24\right)C_{A}g_{s}^{2}}{960(D-4)(D-3)(D-1)^{2}m^{4}}(-4C_{A}g_{s}^{2}(20\left(2D^{2}-3D-11\right)m^{2}$
(43a)
$\displaystyle+\left(2D^{2}-11D+12\right)p^{2})-5\left(D^{2}-8D+12\right)\kappa^{2}m^{2}\left((D-8)p^{2}-48m^{2}\right));$
$\displaystyle c_{2}$ $\displaystyle=$
$\displaystyle-\frac{iC_{A}g_{s}^{2}}{16(D-4)(D-3)(D-1)^{2}m^{2}p^{2}}(-64(D-1)^{2}\left(D^{2}-7D+12\right)\lambda
m^{2}$
$\displaystyle+8(D-1)C_{A}g_{s}^{2}\left(4\left(D^{3}-8D^{2}+19D-16\right)m^{2}+(D-2)Dp^{2}\right)+2D^{6}\kappa^{2}m^{4}-18D^{5}\kappa^{2}m^{4}$
$\displaystyle-D^{5}\kappa^{2}m^{2}p^{2}+22D^{4}\kappa^{2}m^{4}-64D^{4}\lambda
m^{2}+23D^{4}\kappa^{2}m^{2}p^{2}+262D^{3}\kappa^{2}m^{4}+576D^{3}\lambda
m^{2}$
$\displaystyle-196D^{3}\kappa^{2}m^{2}p^{2}-1124D^{2}\kappa^{2}m^{4}-1728D^{2}\lambda
m^{2}+696D^{2}\kappa^{2}m^{2}p^{2}+8D^{2}\kappa^{2}p^{4}+1712D\kappa^{2}m^{4}$
$\displaystyle+1984D\lambda
m^{2}-1048D\kappa^{2}m^{2}p^{2}-24D\kappa^{2}p^{4}-928\kappa^{2}m^{4}-768\lambda
m^{2}+544\kappa^{2}m^{2}p^{2}+16\kappa^{2}p^{4});$ $\displaystyle c_{3}$
$\displaystyle=$
$\displaystyle-\frac{i\left(D^{3}-8D^{2}+19D-12\right)C_{A}g_{s}^{2}\left(2C_{A}g_{s}^{2}+\kappa^{2}\left(2(D-2)m^{2}-(D-4)p^{2}\right)\right)}{2(D-4)(D-3)(D-1)^{2}};$
(43b) $\displaystyle c_{4}$ $\displaystyle=$
$\displaystyle\frac{i\left(3D^{4}-40D^{3}+180D^{2}-320D+192\right)C_{A}g_{s}^{2}}{960(D-6)(D-5)(D-4)^{2}(D-3)(D-2)(D-1)^{2}(3D-4)m^{4}p^{4}}(-1920(D-1)^{2}(D^{4}-14D^{3}$
(43c) $\displaystyle+71D^{2}-154D+120)\lambda
m^{2}p^{2}+4\left(D^{2}-3D+2\right)C_{A}g_{s}^{2}(\left(2D^{3}-19D^{2}+54D-45\right)(D-4)^{2}p^{4}$
$\displaystyle+32\left(4D^{5}-48D^{4}+113D^{3}+616D^{2}-3099D+3470\right)m^{4}+4(8D^{5}-40D^{4}-281D^{3}+2224D^{2}$
$\displaystyle-4899D+3924)m^{2}p^{2})+5(D-5)m^{2}p^{2}(\left(D^{2}-3D+2\right)((D^{5}-23D^{4}+200D^{3}-820D^{2}+1584D$
$\displaystyle-1056)\kappa^{2}p^{2}-384\left(D^{3}-8D^{2}+19D-12\right)\lambda)+4(5D^{7}-113D^{6}+1052D^{5}-5122D^{4}+13896D^{3}$
$\displaystyle-20896D^{2}+16032D-4800)\kappa^{2}m^{2}));$ $\displaystyle
c_{5}$ $\displaystyle=$
$\displaystyle\frac{iC_{A}g_{s}^{2}}{128(D-4)(D-1)^{2}}(64\left(D^{3}-5D^{2}+2D+2\right)C_{A}g_{s}^{2}+(-24D^{5}+497D^{4}-3680D^{3}+12984D^{2}$
(43d) $\displaystyle-21560D+11840);\kappa^{2}p^{2})$ $\displaystyle c_{6}$
$\displaystyle=$
$\displaystyle\frac{iC_{A}g_{s}^{2}}{64(D-4)(D-1)^{2}p^{2}}(\kappa^{2}(16\left(D^{3}-10D^{2}+36D-36\right)m^{4}-8\left(D^{3}-10D^{2}+48D-48\right)m^{2}p^{2}$
(43e)
$\displaystyle+\left(D^{3}-10D^{2}+64D-64\right)p^{4})-128(D-1)C_{A}g_{s}^{2}\left(2m^{2}-p^{2}\right));$
$\displaystyle c_{7}$ $\displaystyle=$
$\displaystyle-\frac{iC_{A}g_{s}^{2}}{48(D-6)(D-4)^{2}(D-1)(3D-4)p^{2}}(24(9D^{6}-189D^{5}+1364D^{4}-4756D^{3}+9280D^{2}$
(43f)
$\displaystyle-10336D+4992)C_{A}g_{s}^{2}+(6D^{8}-35D^{7}-2454D^{6}+39327D^{5}-240012D^{4}+695044D^{3}$
$\displaystyle-915664D^{2}+366464D+98304)\kappa^{2}p^{2});$ $\displaystyle
c_{8}$ $\displaystyle=$
$\displaystyle-\frac{iC_{A}g_{s}^{2}}{480(D-4)(D-2)(D-1)m^{2}p^{4}}(4(D-2)C_{A}g_{s}^{2}(32\left(12D^{4}-92D^{3}-41D^{2}+1577D-2776\right)m^{4}$
(43g)
$\displaystyle+4\left(24D^{4}-172D^{3}+273D^{2}+193D-516\right)m^{2}p^{2}+\left(6D^{4}-67D^{3}+271D^{2}-468D+288\right)p^{4})$
$\displaystyle+5\kappa^{2}m^{2}p^{2}(4\left(6D^{6}-213D^{5}+2417D^{4}-12716D^{3}+34112D^{2}-45272D+23616\right)m^{2}$
$\displaystyle+\left(3D^{6}-63D^{5}+518D^{4}-2092D^{3}+4296D^{2}-3968D+1024\right)p^{2}));$
$\displaystyle c_{9}$ $\displaystyle=$
$\displaystyle\frac{iC_{A}g_{s}^{2}}{480(D-4)(D-3)(D-2)(D-1)m^{2}p^{4}}(4(D-2)C_{A}g_{s}^{2}(240\left(7D^{2}-57D+100\right)m^{4}p^{2}$
(43h)
$\displaystyle-(D-4)^{2}\left(2D^{2}-9D+9\right)p^{6}+128\left(4D^{4}-32D^{3}-7D^{2}+548D-1041\right)m^{6}$
$\displaystyle-4\left(6D^{4}-39D^{3}-22D^{2}+517D-876\right)m^{2}p^{4})+5\kappa^{2}m^{2}p^{2}(16(2D^{6}-69D^{5}+789D^{4}-4236D^{3}$
$\displaystyle+11684D^{2}-16012D+8664)m^{4}-4(D^{6}-44D^{5}+543D^{4}-3040D^{3}+8736D^{2}-12616D$
$\displaystyle+7296)m^{2}p^{2}-\left(D^{6}-25D^{5}+246D^{4}-1220D^{3}+3224D^{2}-4416D+2496\right)p^{4}));$
$\displaystyle c_{10}$ $\displaystyle=$
$\displaystyle\frac{i\kappa^{2}m^{2}C_{A}g_{s}^{2}\left(\left(D^{2}-6D+4\right)p^{2}-4(D-2)m^{2}\right)}{2(D-1)};$
(43i) $\displaystyle c_{11}$ $\displaystyle=$
$\displaystyle-\frac{iC_{A}g_{s}^{2}\left(C_{A}g_{s}^{2}\left(8m^{2}-p^{2}\right)+(D-2)\kappa^{2}m^{2}\left((D-4)p^{2}-8m^{2}\right)\right)}{2(D-1)}.$
(43j)

## References

  * (1) J. F. Donoghue, General relativity as an effective field theory: The leading quantum corrections, Phys. Rev. D 50, 3874-3888 (1994) doi:10.1103/PhysRevD.50.3874 [arXiv:gr-qc/9405057 [gr-qc]].
  * (2) C. P. Burgess, Quantum gravity in everyday life: General relativity as an effective field theory, Living Rev. Rel. 7, 5-56 (2004) doi:10.12942/lrr-2004-5 [arXiv:gr-qc/0311082 [gr-qc]].
  * (3) I. L. Buchbinder, S. Odintsov, and L. Shapiro, Effective Action in Quantum Gravity (CRC Press, Boca Raton, 1992).
  * (4) G. ’t Hooft and M. J. G. Veltman, One loop divergencies in the theory of gravitation, Annales Poincare Phys. Theor. A 20, 69 (1974).
  * (5) S. Deser and P. van Nieuwenhuizen, Nonrenormalizability of the Quantized Einstein-Maxwell System, Phys. Rev. Lett. 32, 245-247 (1974) doi:10.1103/PhysRevLett.32.245
  * (6) S. Deser and P. van Nieuwenhuizen, Nonrenormalizability of the Quantized Dirac-Einstein System, Phys. Rev. D 10, 411 (1974) doi:10.1103/PhysRevD.10.411
  * (7) M. Srednicki, Quantum Field Theory (Cambridge University Press, New York, 2007).
  * (8) D. J. Gross and F. Wilczek, Ultraviolet Behavior of Nonabelian Gauge Theories, Phys. Rev. Lett. 30, 1343-1346 (1973) doi:10.1103/PhysRevLett.30.1343
  * (9) H. D. Politzer, Reliable Perturbative Results for Strong Interactions?, Phys. Rev. Lett. 30, 1346-1349 (1973) doi:10.1103/PhysRevLett.30.1346
  * (10) D. J. Gross and A. Neveu, Dynamical Symmetry Breaking in Asymptotically Free Field Theories, Phys. Rev. D 10, 3235 (1974) doi:10.1103/PhysRevD.10.3235
  * (11) S. P. Robinson and F. Wilczek, Gravitational correction to running of gauge couplings, Phys. Rev. Lett. 96, 231601 (2006) doi:10.1103/PhysRevLett.96.231601 [arXiv:hep-th/0509050 [hep-th]].
  * (12) A. R. Pietrykowski, Gauge dependence of gravitational correction to running of gauge couplings, Phys. Rev. Lett. 98, 061801 (2007) doi:10.1103/PhysRevLett.98.061801 [arXiv:hep-th/0606208 [hep-th]].
  * (13) J. C. C. Felipe, L. C. T. Brito, M. Sampaio and M. C. Nemes, Quantum gravitational contributions to the beta function of quantum electrodynamics, Phys. Lett. B 700, 86 (2011). [arXiv:1103.5824 [hep-th]].
  * (14) J. C. C. Felipe, L. A. Cabral, L. C. T. Brito, M. Sampaio and M. C. Nemes, Ambiguities in the gravitational correction of quantum electrodynamics running coupling, Mod. Phys. Lett. A 28, 1350078 (2013). [arXiv:1205.6779 [hep-th]].
  * (15) D. Ebert, J. Plefka and A. Rodigast, Absence of gravitational contributions to the running Yang-Mills coupling, Phys. Lett. B 660, 579-582 (2008) doi:10.1016/j.physletb.2008.01.037 [arXiv:0710.1002 [hep-th]].
  * (16) N. K. Nielsen, The Einstein-Maxwell system, Ward identities, and the Vilkovisky construction, Annals Phys. 327, 861-892 (2012) doi:10.1016/j.aop.2011.12.010 [arXiv:1109.2699 [hep-th]].
  * (17) D. J. Toms, Cosmological constant and quantum gravitational corrections to the running fine structure constant, Phys. Rev. Lett. 101, 131301 (2008) doi:10.1103/PhysRevLett.101.131301 [arXiv:0809.3897 [hep-th]].
  * (18) D. J. Toms, Quantum gravitational contributions to quantum electrodynamics, Nature 468, 56-59 (2010) doi:10.1038/nature09506 [arXiv:1010.0793 [hep-th]].
  * (19) J. Ellis and N. E. Mavromatos, On the Interpretation of Gravitational Corrections to Gauge Couplings, Phys. Lett. B 711, 139-142 (2012) doi:10.1016/j.physletb.2012.04.005 [arXiv:1012.4353 [hep-th]].
  * (20) M. M. Anber, J. F. Donoghue and M. El-Houssieny, Running couplings and operator mixing in the gravitational corrections to coupling constants, Phys. Rev. D 83, 124003 (2011) doi:10.1103/PhysRevD.83.124003 [arXiv:1011.3229 [hep-th]].
  * (21) L. Ibiapina Bevilaqua, A. C. Lehum and A. J. da Silva, Effective field theory of quantum gravity coupled to scalar electrodynamics, Class. Quant. Grav. 33, no.9, 095008 (2016) doi:10.1088/0264-9381/33/9/095008 [arXiv:1506.00027 [hep-th]].
  * (22) L. I. Bevilaqua, M. Dias, A. C. Lehum, C. R. Senise, A. J. da Silva and H. Souza, Gravitational corrections to two-loop beta function in quantum electrodynamics, Phys. Rev. D 104, no.12, 125001 (2021) doi:10.1103/PhysRevD.104.125001 [arXiv:2105.12577 [hep-th]].
  * (23) L. I. Bevilaqua, A. C. Lehum and H. Souza, Universality of gauge coupling constant in the Einstein-QED system, Phys. Rev. D 104, no.12, 125019 (2021) doi:10.1103/PhysRevD.104.125019 [arXiv:2105.12732 [hep-th]].
  * (24) H. Souza, L. Ibiapina Bevilaqua and A. C. Lehum, Gravitational corrections to a non-Abelian gauge theory, Phys. Rev. D 106, no.4, 045010 (2022) doi:10.1103/PhysRevD.106.045010 [arXiv:2206.02941 [hep-th]].
  * (25) I. L. Buchbinder and S. D. Odintsov, One loop renormalization of the Yang-Mills field theory in a curved space-time, Sov. Phys. J. 26 (1983), 359-361 doi:10.1007/BF01882976
  * (26) Y. Tang and Y. L. Wu, Gravitational Contributions to the Running of Gauge Couplings, Commun. Theor. Phys. 54, 1040-1044 (2010) doi:10.1088/0253-6102/54/6/15 [arXiv:0807.0331 [hep-ph]].
  * (27) Y. Tang and Y. L. Wu, Gravitational Contributions to Gauge Green’s Functions and Asymptotic Free Power-Law Running of Gauge Coupling, JHEP 11, 073 (2011) doi:10.1007/JHEP11(2011)073 [arXiv:1109.4001 [hep-ph]].
  * (28) A. O. Barvinsky, A. Y. Kamenshchik and I. P. Karmazin, The Renormalization group for nonrenormalizable theories: Einstein gravity with a scalar field, Phys. Rev. D 48, 3677-3694 (1993) doi:10.1103/PhysRevD.48.3677 [arXiv:gr-qc/9302007 [gr-qc]].
  * (29) I. L. Shapiro and H. Takata, One loop renormalization of the four-dimensional theory for quantum dilaton gravity, Phys. Rev. D 52, 2162-2175 (1995) doi:10.1103/PhysRevD.52.2162 [arXiv:hep-th/9502111 [hep-th]].
  * (30) T. Cohen, J. Kearney, A. Pierce and D. Tucker-Smith, Singlet-Doublet Dark Matter, Phys. Rev. D 85, 075003 (2012) doi:10.1103/PhysRevD.85.075003 [arXiv:1109.2604 [hep-ph]].
  * (31) N. Arkani-Hamed, D. P. Finkbeiner, T. R. Slatyer and N. Weiner, A Theory of Dark Matter, Phys. Rev. D 79, 015014 (2009) doi:10.1103/PhysRevD.79.015014 [arXiv:0810.0713 [hep-ph]].
  * (32) X. Calmet and F. Kuipers, Implications of quantum gravity for dark matter, Int. J. Mod. Phys. D 30, no.14, 2142004 (2021) doi:10.1142/S0218271821420049 [arXiv:2107.13529 [hep-ph]].
  * (33) N. R. G. and and A. Dasgupta, Interaction of Gravitational Waves with Yang-Mills fields, [arXiv:2212.02416 [gr-qc]].
  * (34) S. Y. Choi, J. S. Shim and H. S. Song, Factorization and polarization in linearized gravity, Phys. Rev. D 51, 2751-2769 (1995) doi:10.1103/PhysRevD.51.2751 [arXiv:hep-th/9411092 [hep-th]].
  * (35) A. Alloul, N. D. Christensen, C. Degrande, C. Duhr and B. Fuks, FeynRules 2.0 - A complete toolbox for tree-level phenomenology, Comput. Phys. Commun. 185, 2250-2300 (2014) doi:10.1016/j.cpc.2014.04.012 [arXiv:1310.1921 [hep-ph]].
  * (36) T. Hahn, Generating Feynman diagrams and amplitudes with FeynArts 3, Comput. Phys. Commun. 140, 418-431 (2001) doi:10.1016/S0010-4655(01)00290-9 [arXiv:hep-ph/0012260 [hep-ph]].
  * (37) V. Shtabovenko, R. Mertig and F. Orellana, FeynCalc 9.3: New features and improvements, Comput. Phys. Commun. 256, 107478 (2020) doi:10.1016/j.cpc.2020.107478 [arXiv:2001.04407 [hep-ph]].
  * (38) A. A. Slavnov, Ward Identities in Gauge Theories, Theor. Math. Phys. 10, 99-107 (1972) doi:10.1007/BF01090719
  * (39) J. C. Taylor, Ward Identities and Charge Renormalization of the Yang-Mills Field, Nucl. Phys. B 33, 436-444 (1971) doi:10.1016/0550-3213(71)90297-5
  * (40) S. Folkerts, D. F. Litim and J. M. Pawlowski, Asymptotic freedom of Yang-Mills theory with gravity, Phys. Lett. B 709, 234-241 (2012) doi:10.1016/j.physletb.2012.02.002 [arXiv:1101.5552 [hep-th]].
  * (41) R. Mertig and R. Scharf, TARCER: A Mathematica program for the reduction of two loop propagator integrals, Comput. Phys. Commun. 111, 265-273 (1998) doi:10.1016/S0010-4655(98)00042-3 [arXiv:hep-ph/9801383 [hep-ph]].
  * (42) O. V. Tarasov, Generalized recurrence relations for two loop propagator integrals with arbitrary masses, Nucl. Phys. B 502, 455-482 (1997) doi:10.1016/S0550-3213(97)00376-3 [arXiv:hep-ph/9703319 [hep-ph]].
  * (43) S. P. Martin and D. G. Robertson, TSIL: A Program for the calculation of two-loop self-energy integrals, Comput. Phys. Commun. 174, 133-151 (2006) doi:10.1016/j.cpc.2005.08.005 [arXiv:hep-ph/0501132 [hep-ph]].
  * (44) S. P. Martin, Evaluation of two loop selfenergy basis integrals using differential equations, Phys. Rev. D 68, 075002 (2003) doi:10.1103/PhysRevD.68.075002 [arXiv:hep-ph/0307101 [hep-ph]].
  * (45) J. F. Donoghue, A Critique of the Asymptotic Safety Program, Front. in Phys. 8, 56 (2020) doi:10.3389/fphy.2020.00056 [arXiv:1911.02967 [hep-th]].
  * (46) D. J. Toms, Quadratic divergences and quantum gravitational contributions to gauge coupling constants, Phys. Rev. D 84, 084016 (2011) doi:10.1103/PhysRevD.84.084016
  * (47) K. S. Stelle, Renormalization of Higher Derivative Quantum Gravity, Phys. Rev. D 16, 953-969 (1977) doi:10.1103/PhysRevD.16.953
  * (48) S. D. Odintsov and I. L. Shapiro, General relativity as the low-energy limit in higher derivative quantum gravity, Class. Quant. Grav. 9, 873-882 (1992) doi:10.1088/0264-9381/9/4/006
  * (49) A. Salvio and A. Strumia, Agravity, JHEP 06, 080 (2014) doi:10.1007/JHEP06(2014)080 [arXiv:1403.4226 [hep-ph]].
  * (50) J. F. Donoghue and G. Menezes, Gauge Assisted Quadratic Gravity: A Framework for UV Complete Quantum Gravity, Phys. Rev. D 97, no.12, 126005 (2018) doi:10.1103/PhysRevD.97.126005 [arXiv:1804.04980 [hep-th]].
  * (51) J. F. Donoghue and G. Menezes, On quadratic gravity, Nuovo Cim. C 45, no.2, 26 (2022) doi:10.1393/ncc/i2022-22026-7 [arXiv:2112.01974 [hep-th]].

Figure 7: Feynman diagrams to the scattering between gluons up to one-loop
order and one graviton exchange. Figure 8: Feynman diagrams to the scattering
between gluons and quarks up to one-loop order and one graviton exchange.
Figure 9: Feynman diagrams to the gluon self-energy involving only one
graviton exchange at two-loop order. Figure 10: Gluon self-energy 1-loop
diagrams with counterterms insertions.


CTPU-PTC-24-02
EPHOU-24-001

Finite modular axion and radiative moduli stabilization
Tetsutaro Higaki${}^{\ a}$<EMAIL_ADDRESS>and Junichiro
Kawamura${}^{\ b}$<EMAIL_ADDRESS>and Tatsuo Kobayashi${}^{\ c}$
<EMAIL_ADDRESS>

a Department of Physics, Keio University, Yokohama, 223-8522, Japan

b Center for Theoretical Physics of the Universe, Institute for Basic Science
(IBS), Daejeon 34051, Korea

c Department of Physics, Hokkaido University, Sapporo 060-0810, Japan

We propose a simple setup which can stabilize a modulus field of the finite
modular symmetry by the Coleman-Weinberg potential. Our scenario leads to a
large hierarchy suppressing instanton-like corrections $e^{2\pi i\tau}$ and to
a light axion identified as $\mathrm{Re}\tau$, where $\tau$ is the modulus
field. This stabilization mechanism provides the axion solution to the strong
CP problem. The potential has a minimum at a large $\mathrm{Im}\tau$ which
suppresses explicit $U(1)_{\mathrm{PQ}}$ violation terms proportional to
$e^{-2\pi{\mathrm{Im}\tau}}$, and hence the quality of the axion is ensured by
the residual symmetry associated with the $T$-transformation, $\tau\to\tau+1$,
around the fixed point $\tau\sim i\infty$.

##  1 Introduction

The finite modular symmetry has been intensively studied to explain the flavor
structure of the Standard Model (SM) quarks and leptons [1, 2, 3, 4, 5, 6, 7,
8, 9]. (See for reviews [10, 11].) Under the symmetry, a Yukawa coupling
constant behaves as the so-called modular form which is a holomorphic function
of the modulus $\tau$. Interestingly, the finite modular symmetries
$\Gamma_{N}$, $N\in\mathbb{N}$, are isomorphic to the non-Abelian discrete
symmetries [12], which can explain the flavor structure of the SM fermions
[13, 14, 15, 16, 17, 18, 19, 20, 21, 22]. Such modular dependence of the
Yukawa coupling to a modulus was found in heterotic orbifold models [23, 24,
25, 26, 27, 28] and magnetized D-brane models [29, 30, 31, 32, 33, 34, 35] in
superstring theory. Moreover, the finite modular symmetry can explain the
hierarchies of the quark and lepton masses and mixing if the modulus is
stabilized near one of the fixed points where the residual symmetry remains
unbroken [36, 37, 38, 39, 40, 41, 42, 43, 44, 45, 46, 47].

In the models with finite modular flavor symmetries, the value of modulus
$\tau$ plays an essential role to fit to the data and is often treated as a
free parameter. Whereas in a UV theory, its value should be fixed dynamically
through moduli stabilization. There are studies on the moduli stabilization
with the three-form fluxes [48, 49], the negative power of modular form [50,
51] 111 The positive power of the trivial singlet is also studied in Ref.
[51]. See also Ref. [52]. , and the general superpotential invariant under the
modular symmetry $SL(2,\mathbb{Z})$ [53, 54, 55, 56, 57], utilizing the Klein
$j$ function [58]. The latter two mechanisms will need non-perturbative
dynamics to generate those potentials. Recently, effect of the radiative
correction to the tree-level stabilization is discussed in Ref. [59].

In this work, we point out that the modulus can be stabilized only by the
Coleman-Weinberg (CW) potential generated by couplings with matter fields.
This mechanism provides a perturbative way to stabilize a modulus through the
supersymmetry breaking. For illustration, we shall consider a simple model
with vector-like quarks which transforms as non-trivial singlets under the
finite modular flavor symmetry. The resultant potential has a global minimum
at $\mathrm{Im}\,\tau\gg 1$, where the residual $\mathbb{Z}^{T}_{N}$ symmetry
associated with the $T$ transformation, $\tau\to\tau+1$, is unbroken. That can
lead to a large hierarchy suppressing instanton-like correction, i.e.,
$e^{-2\pi{\mathrm{Im}\,\tau}}\ll 1$. An interesting feature of this mechanism
is that the $\mathrm{Re}\,\tau$ direction can be identified as the QCD axion
due to the accidental Peccei-Quinn (PQ) symmetry originated from the
$\mathbb{Z}^{T}_{N}$ residual symmetry. We shall discuss the condition that
the axion solution to the strong CP problem is not spoiled due to the CW
potential, i.e. the axion quality is good enough.

This paper is organized as follows. The finite modular symmetry is briefly
reviewed in Sec. 2. In Sec. 3, we discuss the mechanism of the modular
stabilization by utilizing the CW potential, and then we argue existence of
the QCD axion in the setup in Sec. 4. Section 5 concludes.

##  2 Finite modular symmetry

We consider the series of groups

$\displaystyle\Gamma(N):=\left\\{\begin{pmatrix}a&b\\\
c&d\end{pmatrix}\in\mathrm{SL}(2,\mathbb{Z}),\quad\begin{pmatrix}a&b\\\
c&d\end{pmatrix}\equiv\begin{pmatrix}1&0\\\
0&1\end{pmatrix}\quad\mathrm{mod}~{}N\right\\},$ (2.1)

where $N\in\mathbb{N}$ is called level and
$\Gamma:=SL(2,\mathbb{Z})=\Gamma(1)$ is a group of $2\times 2$ matrices with
determinant unity. The modular transformation $\gamma\in\Gamma(N)$ for a
complex parameter $\tau$ is defined as

$\displaystyle\tau\to\gamma\tau=\frac{a\tau+b}{c\tau+d}.$ (2.2)

The modular group is generated by three generators

$\displaystyle S=\begin{pmatrix}0&1\\\ -1&0\end{pmatrix},\quad
T=\begin{pmatrix}1&1\\\ 0&1\end{pmatrix},\quad R=\begin{pmatrix}-1&0\\\
0&-1\end{pmatrix}.$ (2.3)

The finite modular symmetry is defined as the quotient group
$\Gamma_{N}:=\overline{\Gamma}/\Gamma(N)$ where
$\overline{\Gamma}:=\Gamma(1)/\mathbb{Z}^{R}_{2}$. Here $\mathbb{Z}^{R}_{2}$
is the $\mathbb{Z}_{2}$ symmetry generated by $R$. Under $\Gamma_{N}$ with
$N<6$, the generators satisfy

$\displaystyle S^{2}=(ST)^{3}=T^{N}=1.$ (2.4)

There are Abelian discrete symmetries $\mathbb{Z}^{S}_{2}$,
$\mathbb{Z}^{ST}_{3}$ and $\mathbb{Z}^{T}_{N}$ in associated with the $S$,
$ST$ and $T$ generators which are unbroken at $\tau=i$, $w:=e^{2\pi i/3}$ and
$i\infty$, respectively. The finite modular symmetries are isomorphic to the
non-Abelian discrete symmetries, e.g. $\Gamma_{3}\simeq A_{4}$. A modular form
$Y^{({k})}_{r}(\tau)$ with a modular weight $k$ and representation $r$
transforms under $\Gamma_{N}$ as

$\displaystyle Y^{({k})}_{r}(\tau)\to(c\tau+d)^{k}\rho(r)Y^{({k})}_{r}(\tau),$
(2.5)

where $\rho(r)$ is the representation matrix of $\Gamma_{N}$. We assume that
the chiral superfield $Q$ with weight $-k_{Q}$ and representation $r_{Q}$
transforms as

$\displaystyle Q\to(c\tau+d)^{-k_{Q}}\rho(r_{Q})Q.$ (2.6)

More detailed discussions can be found in e.g. Refs. [10, 11].

##  3 Radiative moduli stabilization

We consider a simple supersymmetric model with the following Kähler potneital
and superpotnetial:

$\displaystyle K=$ $\displaystyle\
-h\log(-i\tau+i\tau^{\dagger})+\sum_{i}\left(\frac{Q_{i}^{\dagger}Q_{i}}{(-i\tau+i\tau^{\dagger})^{k_{Q_{i}}}}+\frac{\overline{Q}_{i}^{\dagger}\overline{Q}_{i}}{(-i\tau+i\tau^{\dagger})^{k_{\overline{Q}_{i}}}}\right),$
(3.1) $\displaystyle W=$ $\displaystyle\
\sum_{i}M_{Q_{i}}Y^{({k_{i}})}_{r_{i}}(\tau)\overline{Q}_{i}Q_{i},$

where $h\in\mathbb{N}$. The reduced Planck mass $M_{p}=2.4\times
10^{18}~{}\mathrm{GeV}$ is set to unity. The chiral superfield $Q_{i}$
($\overline{Q}_{i}$) has a modular weight $k_{Q_{i}}$ ($k_{\overline{Q}_{i}}$)
and is assumed to be a singlet under $\Gamma_{N}$ for simplicity. Here,
$Q_{i}$ and $\overline{Q}_{i}$ are vector-like pairs under the SM gauge group.
The weight of the modular form $k_{i}=k_{Q_{i}}+k_{\overline{Q}_{i}}-h$ and
the representation $r_{i}$ are chosen so that the combination $e^{K}|W|^{2}$
in supergravity is invariant under $\Gamma_{N}$. We introduce the vector-like
mass parameter $M_{Q_{i}}$ which could be replaced by vacuum expectation
values of fields. Throughout this work, we assume that the SM fields are
trivial singlets under the modular symmetry $\Gamma_{N}$ for simplicity and
they have no effects in the following analysis.

We shall show that the modulus $\tau$ can be stabilized without introducing a
superpotential for the modulus $\tau$, but with the 1-loop CW potential

$\displaystyle V_{\mathrm{CW}}=\frac{1}{32\pi^{2}}\sum_{i}$
$\displaystyle\left[\left(m_{i}^{2}+m_{Q_{i}}^{2}(\tau)\right)^{2}\left(\log\left(\frac{m_{i}^{2}+m_{Q_{i}}^{2}(\tau)}{\mu^{2}}\right)-\frac{3}{2}\right)\right.$
(3.2) $\displaystyle\hskip
85.35826pt\left.-(m_{Q_{i}}^{2}(\tau))^{2}\left(\log\left(\frac{m_{Q_{i}}^{2}(\tau)}{\mu^{2}}\right)-\frac{3}{2}\right)\right],$

with

$\displaystyle
m_{Q_{i}}^{2}(\tau):=M_{Q_{i}}^{2}\left(-i\tau+i\tau^{\dagger}\right)^{k_{i}}\left|{Y^{({k_{i}})}_{r_{i}}(\tau)}\right|^{2}.$
(3.3)

Here, we consider the $\overline{\mathrm{MS}}$ scheme with $\mu$ being the
renormalization scale. The first (second) term in the parenthesis corresponds
to the scalar (fermion) contribution. The supersymmetric mass $m_{Q_{i}}$ is
multiplied by the factor $(-i\tau+i\tau^{\dagger})^{k_{i}/2}$ due to the
canonical normalization. We introduce the soft supersymmetry breaking mass
squared $m_{i}^{2}$ for the scalar component of $Q_{i}$ which is assumed to be
independent of $\tau$ for simplicity. This is the case if supersymmetry is
broken by a mechanism irrelevant to the modulus $\tau$. We shall first see the
simplest case analytically at large $\mathrm{Im}\,\tau$, and then see the
results numerically later. Throughout this paper, we assume that the mass
parameters are common for $Q_{i}$’s, i.e. $M_{Q_{i}}=:M_{Q}$, $m_{i}=:m_{0}$
and $k_{i}=:k$ for simplicity. Note that the tree-level scalar potential in
supergravity is irrelevant when $Q_{i}=\overline{Q}_{i}=0$.

###  3.1 Simplified analysis

We argue the model of radiative stabilization of the modulus $\tau$ where
there is only one pair of $(\overline{Q},Q)$, so we omit the index $i$. The
singlet modular form of $\Gamma_{N}$ can be expanded as

$\displaystyle
Y^{({k})}_{1_{t}}(\tau)=q^{t/N}\sum_{n=0}^{\infty}c_{n}q^{n},\quad q:=e^{2\pi
i\tau},$ (3.4)

where $r=1_{t}$ is a singlet whose charge is $0\leq t<N$ under the
$\mathbb{Z}^{T}_{N}$ symmetry 222 In the usual notation [1], $1_{0}=1$,
$1_{1}=1^{\prime}$ and $1_{2}=1^{{\prime\prime}}$ for $A_{4}$. . Assuming
$\left|{q}\right|=e^{-2\pi\mathrm{Im}\,\tau}\ll 1$ and $m_{0}^{2}\ll
M_{Q}^{2}$, the CW potential is given by

$\displaystyle
V_{\mathrm{CW}}=\frac{m_{0}^{2}\tilde{M}_{Q}^{2}}{16\pi^{2}}x^{k}e^{-px}\left(\log\frac{\tilde{M}_{Q}^{2}}{e\mu^{2}}+k\log
x-px\right)+\mathcal{O}\left({\left|{q}\right|,m_{0}^{4}}\right),$ (3.5)

where $x:=2\mathrm{Im}\,\tau$, $p:=2\pi t/N$ and $\tilde{M}_{Q}:=c_{0}M_{Q}$.
At the leading order, there is no potential for $\mathrm{Re}\,\tau$ which will
be important for the axion interpretation discussed in the next section. The
first derivative of the potential is given by

$\displaystyle\frac{dV_{\mathrm{CW}}}{dx}\sim\frac{m_{0}^{2}\tilde{M}_{Q}^{2}}{16\pi^{2}}x^{k-1}e^{-px}\left(k-px\right)\left(\log\frac{\tilde{M}_{Q}^{2}}{\mu^{2}}+k\log
x-px\right).$ (3.6)

The point $x=k/p$ is a maximum for sufficiently large $k/p$, and hence this
potential can have a minimum at

$\displaystyle
x_{0}=-\frac{k}{p}\mathcal{W}\left(-\frac{p}{k}\left(\frac{\mu}{\tilde{M}_{Q}}\right)^{2/k}\right),$
(3.7)

where the second parenthesis of Eq. (3.6) is vanishing. Here, $\mathcal{W}(z)$
is the Lambert function satisfying $\mathcal{W}(z)e^{\mathcal{W}(z)}=z$ which
has two real values for $-e^{-1}<z<0$, and the larger one can be a solution
where our approximation $\left|{q}\right|\ll 1$ is valid. Taking
$\mu=\tilde{M}_{Q}$ and $p=2\pi/3~{}$, the location of the minimum is
$x_{0}=3.9,~{}7.9,~{}12$ and $16$ for $k=6,8,10$ and $12$, respectively.
$x_{0}$ increases slightly by decreasing $\mu$. Thus we find that there is the
minimum at $\mathrm{Im}\,\tau\gg 1$ in the CW potential. Note that this
minimum exists only for $t>0$ where the modular form $Y^{({k})}_{1_{t}}$ is a
non-trivial singlet. For the CW potential from the trivial-singlet modular
form, the potential increases as $\mathrm{Im}\,\tau$ increases in contrast to
that from the non-trivial one, and thus there is no minimum at
$\mathrm{Im}\,\tau\gg 1$ 333 The potentials discussed in Refs. [53, 54, 55,
56, 57, 58] have a similar behavior, that is, $V\to\infty$ for
$\mathrm{Im}\,\tau\to\infty$.  See also the potential in Ref. [52]. Therefore
the dominant CW potential should be generated from the non-trivial singlet one
for the existence of the minimum at $\mathrm{Im}\,\tau\gg 1$, so that the
potential tends to approach to zero for a sufficiently large ${\rm Im}\tau$.

###  3.2 Numerical analysis

Figure 1:  The shape of the CW potential along the $\mathrm{Im}\,\tau$
direction with $\mathrm{Re}\,\tau=0$ (red solid) and $\mathrm{Re}\,\tau=-0.5$
(blue dashed).

We shall study the potential Eq. (3.2) numerically without relying on the
$q$-expansion in Eq. (3.4). For illustration, we consider $\Gamma_{3}\simeq
A_{4}$ as the finite modular group, and two pairs of vector-like quarks whose
superpotential is given by

$\displaystyle
W=M_{Q}\sum_{i=1,2}Y^{({12})}_{1_{i}}(\tau)\overline{Q}_{i}Q_{i}.$ (3.8)

We can easily choose the representations and modular weights such that only
those two terms are allowed but mixing terms like $\bar{Q}_{1}Q_{2}$ are
forbidden. An explicit assignment of representations and modular weights will
be shown in Eq. (4.13) when we shall discuss the vector-like quark decays. The
modular forms of weight 12 are given by

$\displaystyle
Y^{({12})}_{1_{1}}(\tau)=(Y_{1}^{2}+2Y_{2}Y_{3})^{2}(Y_{3}^{2}+2Y_{1}Y_{2}),\quad
Y^{({12})}_{1_{2}}(\tau)=(Y_{1}^{2}+2Y_{2}Y_{3})(Y_{3}^{2}+2Y_{1}Y_{2})^{2},$
(3.9)

where the functions $Y_{1,2,3}$ are defined as [1]

$\displaystyle Y_{1}(\tau)=$ $\displaystyle\
\frac{i}{2\pi}\left[\frac{\eta^{\prime}(\tau/3)}{\eta(\tau/3)}+\frac{\eta^{\prime}((\tau+1)/3)}{\eta((\tau+1)/3)}+\frac{\eta^{\prime}((\tau+2)/3)}{\eta((\tau+2)/3)}-27\frac{\eta^{\prime}(3\tau)}{\eta(3\tau)}\right],$
(3.10) $\displaystyle Y_{2}(\tau)=$ $\displaystyle\
\frac{-i}{\pi}\left[\frac{\eta^{\prime}(\tau/3)}{\eta(\tau/3)}+w^{2}\frac{\eta^{\prime}((\tau+1)/3)}{\eta((\tau+1)/3)}+w\frac{\eta^{\prime}((\tau+2)/3)}{\eta((\tau+2)/3)}\right],$
(3.11) $\displaystyle Y_{3}(\tau)=$ $\displaystyle\
\frac{-i}{\pi}\left[\frac{\eta^{\prime}(\tau/3)}{\eta(\tau/3)}+w\frac{\eta^{\prime}((\tau+1)/3)}{\eta((\tau+1)/3)}+w^{2}\frac{\eta^{\prime}((\tau+2)/3)}{\eta((\tau+2)/3)}\right],$
(3.12)

with $\eta(\tau)$ being the Dedekind eta function. The $q$-expansion of the
modular forms are given by 444 There is an ambiguity of the normalization of
the modular forms which can not be determined from the symmetry. We simply
employ the normalization used in Ref. [1]. See also Ref. [60].

$\displaystyle
Y^{({12})}_{1_{1}}=-12q^{1/3}\left(1+472q+\mathcal{O}\left({q^{2}}\right)\right),\quad
Y^{({12})}_{1_{2}}=144q^{2/3}\left(1+224q+\mathcal{O}\left({q^{2}}\right)\right).$
(3.13)

The modular weights are chosen to $k=12$ so that the strong CP problem is
solved with keeping the modulus mass heavy, as discussed in the next section
555 Since $1_{1}$ and $1_{2}$ can be in the $2^{\prime}$ representation under
$A_{4}^{\prime}$, the two fields could be embedded into a doublet. . Here, we
introduce the two vector-like pairs, so that the mixed anomaly of $A_{4}$ with
the QCD is canceled [61, 62]. The singlets are trivial under the
$S$-transformation, and transform as $Y^{({k})}_{1_{t}}\to
w^{t}Y^{({k})}_{1_{t}}$ under the $T$-transformation. Hence, in our setup,
$\det(\rho(T))=w^{1+2}=1$ and the anomaly is canceled. The anomaly associated
with the modular weight can be canceled by the Green-Schwarz mechanism with an
another gauge coupling modulus or by adding $\log Y^{({k})}_{1_{0}}(\tau)$ to
the gauge kinetic function with a certain coefficient, where
$Y^{({k})}_{1_{0}}(\tau)$ is a trivial singlet modular form and could be
explained by threshold corrections from heavy modes [63].

Figure 1 shows the shape of the potential along the $\mathrm{Im}\,\tau$
direction. Here, we take $m^{2}/M_{Q}^{2}=10^{-8}$ and $\mu/M_{Q}=10^{-2}$ so
that the minimum resides at $\mathrm{Im}\,\tau\simeq 13$. The red solid (blue
dashed) line corresponds to $\mathrm{Re}\,\tau=0$ ($-0.5$). Note that the
fundamental domain is $\mathrm{Im}\,\tau\geq 1$ ($\sqrt{3}/2$) at
$\mathrm{Re}\,\tau=0~{}(-0.5)$. In the left panel, the origin of the
horizontal axis is $\sqrt{3}/2$ and the global picture of the potential is
shown. The right panel is the shape of the potential around the global minimum
at large $\mathrm{Im}\,\tau$. For $\mathrm{Im}\,\tau\lesssim 1$, where the
$q$-expansion is not efficient, there is a local minimum at $\tau\sim w$, but
$V_{\mathrm{CW}}(\tau=w)=0$, since $Y^{({12})}_{1_{1,2}}(w)=0$, is shallower
than the global minimum at $\mathrm{Im}\,\tau\sim 12$ where the potential
value is negative. The potential clearly depends on $\mathrm{Re}\,\tau$ at
small $\mathrm{Im}\,\tau$, while it is independent of $\mathrm{Re}\,\tau$
within numerical precision. Thus, the CW potential in the $\mathrm{Re}\,\tau$
direction is approximately flat near the global minimum at
$\mathrm{Im}\,\tau\gg 1$. As a result, our scenario of modulus stabilization
can lead a light axion as well as a large hierarchy by
$|q|=e^{-2\pi\mathrm{Im}\,\tau}\sim 10^{-36}$. We will confirm this flatness
analytically in the next section to discuss the quality of the axion solution
to the strong CP problem.

##  4 Axion solution and its quality

We have seen that the CW potential Eq. (3.2) induced by the superpotential Eq.
(3.8) is very flat in the $\mathrm{Re}\,\tau$ direction. The stabilization
mechanism can realize the axion solution to the strong CP problem by assuming
that the matter fields are vector-like quarks as in the KSVZ axion model [64,
65]. The effective $\theta$-angle $\overline{\theta}$ in the QCD is given by

$\displaystyle\overline{\theta}(\tau)=\theta_{0}+\mathrm{Arg}\left(Y^{({12})}_{1_{1}}(\tau)Y^{({12})}_{1_{2}}(\tau)\right)=\theta_{0}+\phi+\mathcal{O}\left({\left|{q}\right|}\right),$
(4.1)

where $\theta_{0}$ is a constant and $\phi:=2\pi\mathrm{Re}\,\tau$. The $\phi$
dependence appears as a result of the chiral anomaly of the QCD with the
approximate $U(1)_{\rm PQ}$, which originates from the residual
$\mathbb{Z}^{T}_{3}$ symmetry as shown below. For $\mathrm{Im}\,\tau\gg 1$,
the Yukawa coupling is given by

$\displaystyle Y^{({k})}_{1_{t}}\overline{Q}_{t}Q_{t}\propto
e^{i(t/3)\phi}\;\overline{Q}_{t}Q_{t},$ (4.2)

so we can find an accidental $U(1)_{\mathrm{PQ}}$ symmetry as

$\displaystyle\phi\to\phi+\alpha,\quad\overline{Q}_{t}Q_{t}\to
e^{-it\alpha/3}\overline{Q}_{t}Q_{t},$ (4.3)

where $\alpha$ is a real transformation parameter. Hence the PQ charge of
$\overline{Q}_{t}Q_{t}$ is $-t$, and the $U(1)_{\rm PQ}$ is anomalous. This
symmetry is ensured from the residual $\mathbb{Z}^{T}_{3}$ symmetry
corresponding to $\alpha=2\pi$, and the explicit breaking effects of the
$U(1)_{\mathrm{PQ}}$ are from $q=e^{2\pi i\tau}$ which is invariant under
$\mathbb{Z}^{T}_{3}$ but violates $U(1)_{\mathrm{PQ}}$ 666 The corrections to
the axion from $\log Y^{({k})}_{1_{0}}$ in the gauge kinetic function which
may exist to cancel the weight anomaly will also appear at
$\mathcal{O}\left({q}\right)$. . Therefore, $\mathrm{Re}\,\tau$ can be
identified as the QCD axion. Since the kinetic term of the modulus is given by
777 We write explicitly the reduced Planck mass $M_{p}$ in this section.

$\displaystyle\frac{hM_{p}^{2}}{(2\mathrm{Im}\,\tau)^{2}}\partial_{\mu}\tau^{\dagger}\partial^{\mu}\tau,$
(4.4)

the axion decay constant is $f_{a}=\sqrt{h}M_{p}/(2\pi
x_{0})\sim\mathcal{O}\left({10^{16}}\right)$ GeV which requires a fine-tuning
of the initial condition or early matter domination [66], to avoid the
overabundance of the axion.

The axion potential is given by

$\displaystyle
V_{a}=-\Lambda_{\mathrm{QCD}}^{4}\cos\left(\overline{\theta}(\phi)\right)+\Delta
V.$ (4.5)

Here, $\Delta V:=V_{\mathrm{CW}}|_{x=x_{0}}$ is the axion-dependent part of
the CW potential which can spoil the axion solution. The shift of the
$\theta$-angle due to $\Delta V$ is estimated as

$\displaystyle\Delta\theta=\frac{1}{\Lambda_{\mathrm{QCD}}^{4}}\left(\frac{\partial\Delta
V}{\partial\phi}\right)\simeq\frac{m_{0}^{2}M_{Q}^{2}}{8\pi^{2}\Lambda_{\mathrm{QCD}}^{4}}x_{0}^{k}\epsilon^{4+N}\Omega\sin{\theta_{0}},$
(4.6)

with

$\displaystyle\Omega:=\left(\beta_{1}\frac{kN-4\pi x_{0}}{kN-2\pi
x_{0}}-\beta_{2}\right)\left(\log\frac{m^{2}_{0}+m_{Q_{2}}^{2}}{e\mu^{2}}+\frac{m_{Q_{2}}^{2}}{m_{0}^{2}}\log\left(1+\frac{m^{2}_{0}}{m_{Q_{2}}^{2}}\right)\right),$
(4.7)

where $\epsilon:=e^{-\pi x_{0}/N}$. In our model, $k=12$ and $N=3$. Here,
$m^{2}_{Q_{2}}$ is the mass squared of $Q_{2}$ defined in Eq. (3.3) with
taking $x=x_{0}$ and dropping the negligible corrections at
$\mathcal{O}\left({q}\right)$. $\beta_{t}$ is the ratio $c_{1}/c_{0}$ in the
expansion of the modular forms $Y^{({k})}_{1_{t}}$ in Eq. (3.4), which is
given by $\beta_{1}=472$ and $\beta_{2}=224$ for $Y^{({12})}_{1_{1}}$ and
$Y^{({12})}_{1_{2}}$, respectively. Note that the
$\mathcal{O}\left({q}\right)$ correction proportional to $m_{Q_{1}}^{2}$ is
canceled at the minimum $x=x_{0}$ and thus the leading term is
$\epsilon^{4+N}$ which is proportional to
$m_{Q_{2}}^{2}\sim\mathcal{O}\left({\epsilon^{4}x_{0}^{k}M_{Q}^{2}}\right)$.
This shift of the angle should be less than
$\mathcal{O}\left({10^{-10}}\right)$ to be consistent with the measurement of
the neutron electric dipole moment [67]. A possible modulus $\tau$-dependence
of $\Lambda_{\mathrm{QCD}}$ from the gauge kinetic function will also be
suppressed by $\mathcal{O}\left({q}\right)$, and thus is sub-dominant for the
stabilization. Similarly the $\phi$ dependence of
$\partial_{\phi}\overline{\theta}$ is also negligible. In the $A_{4}$ model
discussed in the previous section, we obtain

$\displaystyle\left|{\Delta\theta}\right|\simeq 2\times
10^{-10}\times\sin\theta_{0}\left(\frac{\left|{\Omega}\right|}{10^{4}}\right)\left(\frac{m_{0}}{10^{7}~{}\mathrm{GeV}}\right)^{2}\left(\frac{M_{Q}}{M_{p}}\right)^{2}\left(\frac{x_{0}}{28}\right)^{8}\left(\frac{\epsilon}{10^{-12}}\right)^{7}.$
(4.8)

Thus the axion quality is ensured if $\epsilon\lesssim 10^{-12}$, where
$\mathrm{Im}\,\tau\sim 13$ as shown in Fig. 1. The required tuning becomes
mild if $\sin\theta_{0}\ll 1$ which might happen when the (generalized) CP is
conserved at some level as discussed in e.g. Ref. [56] 888See Refs. [68, 69]
for recent studies on the strong CP problem relevant to the modular symmetry..
In such a case, the strong CP problem is partially solved by the modular axion
$\phi\sim\mathrm{Re}\,\tau$.

Figure 2:  Masses of the modulus $X\sim\mathrm{Im}\,\tau$ (red) and vector-
like quarks $Q_{1}$ (blue dashed) and $Q_{2}$ (green dot-dashed). In this
plot, $m_{0}=10^{7}~{}\mathrm{GeV}$, $M_{Q}=10^{4}\times\mu=M_{p}$, and $h=3$.

Before closing, we discuss the masses of the modulus
$X\propto\mathrm{Im}\,\tau$ and the vector-like quarks $Q_{1,2}$. The modulus
mass may need to be heavier than
$\mathcal{O}\left({10~{}\mathrm{TeV}}\right)$, so that the modulus decay
occurs before the big-bang nucleosynthesis and the moduli problem is avoided
[70, 71, 72, 73]. After canonically normalizing the kinetic term Eq. (4.4),
the modulus mass is given by

$\displaystyle
m_{X}^{2}=\frac{m_{0}^{2}\tilde{M}_{Q}^{2}}{8\pi^{2}hM_{p}^{2}}x_{0}^{k}\epsilon^{2}(k-px_{0})^{2}.$
(4.9)

Hence the modulus mass is related to $\Delta\theta$ as

$\displaystyle m_{X}\simeq$ $\displaystyle\
c_{0}\left|{k-px_{0}}\right|\sqrt{\frac{\Delta\theta}{h\Omega\sin{\theta_{0}}}}\frac{\Lambda_{\mathrm{QCD}}^{2}}{\epsilon^{1+N/2}M_{p}}$
(4.10) $\displaystyle\sim$ $\displaystyle\
42~{}\mathrm{TeV}\times\frac{c_{0}\left|{k-2\pi
x_{0}/N}\right|}{\sqrt{h\sin{\theta_{0}}}}\left(\frac{\Delta\theta}{10^{-10}}\right)^{1/2}\left(\frac{10^{4}}{\Omega}\right)^{1/2}\left(\frac{\Lambda_{\mathrm{QCD}}}{100~{}\mathrm{MeV}}\right)^{2}\left(\frac{10^{-12}}{\epsilon}\right)^{5/2},$

in the $A_{4}$ model. Thus the modulus mass decreases as the axion quality is
higher by $\epsilon\to 0$. For $\Delta\theta/\sin\theta_{0}\sim 10^{-10}$,
$\epsilon\lesssim 10^{-12}$ is required to avoid the moduli problem 999 The
lighter modulus would not be a problem if there is the thermal inflation [74,
75, 76], large Hubble induced mass during the inflation exists [77, 78],
and/or an inflation scale is lower than the modulus mass. .

The mass of the lighter vector-like quark $Q_{2}$, whose $T$ charge is $2$, is
given by

$\displaystyle m_{Q_{2}}\sim\tilde{M}_{Q}x_{0}^{k/2}\epsilon^{2}\sim
1~{}\mathrm{TeV}\times\left(\frac{\tilde{M}_{Q}}{M_{p}}\right)\left(\frac{x_{0}}{28}\right)^{6}\left(\frac{\epsilon}{10^{-12}}\right)^{2},$
(4.11)

thus $M_{Q}\sim M_{p}$ is necessary for
$m_{Q_{2}}\gtrsim\mathcal{O}\left({\mathrm{TeV}}\right)$ to be consistent with
the LHC constraints [79, 80, 81, 82, 83]. Conversely, if we require
$m_{Q_{2}}\gtrsim\mathrm{TeV}$ and $\tilde{M}_{Q}\lesssim M_{p}$, the
hierarchy parameter is bounded from below as $\epsilon\gtrsim 10^{-8}\times
x_{0}^{-k/4}$, and thus larger weight is preferred to allow
$\epsilon\sim\mathcal{O}\left({10^{-12}}\right)$ without having the too light
vector-like quark. Figure 2 shows the masses of the modulus
$X\propto\mathrm{Im}\,\tau$ and the vector-like quarks. As $\mathrm{Im}\,\tau$
increases, the axion quality becomes better, whereas the modulus and vector-
like quark masses decrease, and hence there is the $lower$ bound
$\Delta\theta\gtrsim 10^{-13}$ for sufficiently heavy vector-like quark
$Q_{2}$.

The vector-like quarks should be able to decay into SM particles. For
illustration, if the vector-like quarks have the same quantum number as the
down-type quarks, denoted by $d$, we can write down the superpotential

$\displaystyle
W_{\mathrm{mix}}=\sum_{i=1,2}{\lambda}_{i}Y^{({14})}_{1_{1}}H_{d}qQ_{i},$
(4.12)

where $H_{d}$ and $q$ are respectively the down-type Higgs doublet and the
doublet SM quark. Here, the representations and weights of quarks are set to

$\displaystyle r_{q}=r_{d}=r_{\overline{Q}_{1}}=1_{0},\quad
r_{Q_{1,2}}=r_{\overline{Q}_{2}}=1_{2},\quad k_{q}=h-k_{d}=8,\quad
k_{Q_{1,2}}-h=k_{\overline{Q}_{1,2}}=6,$ (4.13)

and those of $H_{d}$ is the trivial-singlet with weight $0$. We omit the
flavor indices of the SM quarks. With this assignment, the interactions in
Eqs. (3.8) and (4.12), as well as the SM Yukawa coupling $yH_{d}qd$ with $y$
being a constant, are realized, while the $\overline{Q}_{i}d$ is vanishing
because of the negative weight. The size of mixing $s_{Q_{i}}$ is estimated as

$\displaystyle s_{Q_{i}}\sim$ $\displaystyle\
\frac{{\lambda}_{i}x_{0}^{7}Y^{({14})}_{1_{1}}{\langle{H_{d}}\rangle}}{M_{Q}x_{0}^{6}Y^{({12})}_{1_{i}}}\sim\frac{{\lambda}_{i}x_{0}{\langle{H_{d}}\rangle}}{M_{Q}\epsilon^{i-1}}$
(4.14) $\displaystyle\sim$ $\displaystyle\
10^{-3}\times{\lambda}_{2}\left(\frac{x}{28}\right)\left(\frac{{\langle{H_{d}}\rangle}}{100~{}\mathrm{GeV}}\right)\left(\frac{M_{p}}{M_{Q}}\right)\left(\frac{10^{-12}}{\epsilon}\right),$

where $i=2$ in the second line, and we used
$Y^{({k})}_{1_{t}}\sim\epsilon^{t}$. Thus the vector-like quarks decay
promptly at the collider scale. It is noted that the CW potential induced by
the mixing term is negligible because ${\langle{H_{d}}\rangle}\ll M_{Q}$.

##  5 Summary and discussions

In this work, we point out that the modulus of the finite modular symmetry can
be stabilized by the Coleman-Weinberg (CW) potential. For illustration, we
studied the model with $\Gamma_{3}\simeq A_{4}$ symmetry and two vector-like
quark pairs which transform as non-trivial singlets, namely $1_{1}$ and
$1_{2}$. This is the minimal possibility to cancel the mixed QCD anomaly of
the finite modular symmetry, but to induce that of the $U(1)_{\mathrm{PQ}}$
symmetry. The CW potential has the global minimum at $\mathrm{Im}\,\tau\gg 1$
if the corresponding modular form is non-trivial singlet under $A_{4}$, where
the residual symmetry $\mathbb{Z}^{T}_{3}$ remains unbroken. Since the
potential along the $\mathrm{Re}\,\tau$ direction has extremely flat due to
this residual symmetry, we can regard $\mathrm{Re}\,\tau$ as the QCD axion to
solve the strong CP problem. The accidental $U(1)_{\mathrm{PQ}}$ arises due to
the $\mathbb{Z}^{T}_{3}$ symmetry and hence the PQ breaking effects are
controlled by $\Gamma_{3}$ which is a different situation argued in Ref. [84].
We examined the condition to ensure the quality of this finite modular axion,
and correlate with the masses of the modulus $X\sim\mathrm{Im}\,\tau$ and
vector-like quarks. Interestingly, the modulus $X$ and the lighter quark
$Q_{2}$ are expected to be $\mathcal{O}\left({\mathrm{TeV}}\right)$ scale for
$\Delta\theta<10^{-10}$ and thus could be probed by cosmology and collider
experiments. This mechanism can be generalized to other modular forms as long
as the CW potential is dominated by that with non-zero $T$-charge and the
potential converges to zero for $\mathrm{Im}\,\tau\gg 1$.

The modular $A_{4}$ symmetry in our scenario can not be used to explain the
hierarchical structure of the SM fermions, as studied in Refs. [36, 37, 38,
39, 40, 41, 42, 43, 44, 45, 46, 47], since the hierarchy parameter $\epsilon$
is $\mathcal{O}\left({10^{-12}}\right)$ and is too small. Such a tiny
$\epsilon$ is required due to the correlation of the axion quality
$\Delta\theta$ with the modulus mass in Eq. (4.10). This relation can be
relaxed if the $\tau$-independent angle $\theta_{0}$ is small or the level $N$
is large. In the latter case, there exist a small hierarchy $|e^{2\pi
i\tau/N}|$ in the SM sector and a large one $|e^{2\pi i\tau}|$ in the vector-
like quark sector.101010 In general, for example, if one could break the
periodicity $\tau\sim\tau+1$ to $\tau\sim\tau+1/n$ ($n>1$) in the vector-like
$(Q,\bar{Q})$ sector with keeping the periodicity $\tau\sim\tau+1$ in the SM
fermion sector, then modular forms may appear $Y^{k}_{r}(\tau)$ in the SM
fermion sector and $Y^{k^{\prime}}_{r^{\prime}}(n\tau)$ in the $(Q,\bar{Q})$
sector, and one could realize both a small hierarchy $|e^{2\pi i\tau}|$ in the
SM sector and a large one $|e^{2\pi in\tau}|$ in the vector-like quark pairs.
. Thus, we could construct a model which explains the SM fermion hierarchies
by the modulus stabilized by the mechanism proposed in this paper. An explicit
model building is subject of our future work.

We have proposed a new scenario to lead to a large hierarchy
$e^{-2\pi\mathrm{Im}\,\tau}\sim\mathcal{O}\left({10^{-36}}\right)$ and a light
axion by stabilizing the modulus at $\mathrm{Im}\,\tau\gg 1$. Although we have
applied it to solve the strong CP problem by assuming $Q$ and $\bar{Q}$ as
vector-like matter fields with the QCD colors, our scenario would be useful to
explain other large hierarchies, which would be required for the proton
stability, tiny neutrino masses, CP/flavor violations, quintessence and so on,
together with light axions by assuming $Q$ and $\bar{Q}$ as visible or hidden
matter fields. We would study it elsewhere.

In this work, we assume that the soft scalar masses of the vector-like quarks
are independent of the modulus $\tau$ and supersymmetry is dominantly broken
by other sector, so that the supersymmetry breaking sector does not contribute
to the stabilization of the modulus $\tau$. In addition, there would be
contributions from gravitational instantons [85] which would change the
potential structure. Those effects and other dynamical effects including
stabilization mechanisms of other moduli, and uplifting of the vacuum energy
if necessary, may become important in complete models such as superstring, but
beyond the scope of this paper.

## Acknowledgement

We thank Y.Abe and N.Aso for helpful discussions in the early stage of this
work. The work of J.K. is supported in part by the Institute for Basic Science
(IBS-R018-D1). This work is supported in part by he Grant-in-Aid for
Scientific Research from the Ministry of Education, Science, Sports and
Culture (MEXT), Japan No. JP22K03601 (T.H.) and JP23K03375 (T.K.).

## References

  * [1] F. Feruglio, Are neutrino masses modular forms?, pp. 227–266.  2019\.  arXiv:1706.08749.
  * [2] T. Kobayashi, K. Tanaka, and T. H. Tatsuishi, Neutrino mixing from finite modular groups, Phys. Rev. D 98 (2018), no. 1 016004, [arXiv:1803.10391].
  * [3] J. T. Penedo and S. T. Petcov, Lepton Masses and Mixing from Modular $S_{4}$ Symmetry, Nucl. Phys. B 939 (2019) 292–307, [arXiv:1806.11040].
  * [4] P. P. Novichkov, J. T. Penedo, S. T. Petcov, and A. V. Titov, Modular A5 symmetry for flavour model building, JHEP 04 (2019) 174, [arXiv:1812.02158].
  * [5] G.-J. Ding, S. F. King, and X.-G. Liu, Neutrino mass and mixing with $A_{5}$ modular symmetry, Phys. Rev. D 100 (2019), no. 11 115005, [arXiv:1903.12588].
  * [6] X.-G. Liu and G.-J. Ding, Neutrino Masses and Mixing from Double Covering of Finite Modular Groups, JHEP 08 (2019) 134, [arXiv:1907.01488].
  * [7] P. P. Novichkov, J. T. Penedo, and S. T. Petcov, Double cover of modular $S_{4}$ for flavour model building, Nucl. Phys. B 963 (2021) 115301, [arXiv:2006.03058].
  * [8] X.-G. Liu, C.-Y. Yao, and G.-J. Ding, Modular invariant quark and lepton models in double covering of $S_{4}$ modular group, Phys. Rev. D 103 (2021), no. 5 056013, [arXiv:2006.10722].
  * [9] X.-G. Liu, C.-Y. Yao, B.-Y. Qu, and G.-J. Ding, Half-integral weight modular forms and application to neutrino mass models, Phys. Rev. D 102 (2020), no. 11 115035, [arXiv:2007.13706].
  * [10] T. Kobayashi and M. Tanimoto, Modular flavor symmetric models, 7, 2023\.  arXiv:2307.03384.
  * [11] G.-J. Ding and S. F. King, Neutrino Mass and Mixing with Modular Symmetry, arXiv:2311.09282.
  * [12] R. de Adelhart Toorop, F. Feruglio, and C. Hagedorn, Finite Modular Groups and Lepton Mixing, Nucl. Phys. B 858 (2012) 437–467, [arXiv:1112.1340].
  * [13] G. Altarelli and F. Feruglio, Discrete Flavor Symmetries and Models of Neutrino Mixing, Rev. Mod. Phys. 82 (2010) 2701–2729, [arXiv:1002.0211].
  * [14] H. Ishimori, T. Kobayashi, H. Ohki, Y. Shimizu, H. Okada, and M. Tanimoto, Non-Abelian Discrete Symmetries in Particle Physics, Prog. Theor. Phys. Suppl. 183 (2010) 1–163, [arXiv:1003.3552].
  * [15] D. Hernandez and A. Y. Smirnov, Lepton mixing and discrete symmetries, Phys. Rev. D 86 (2012) 053014, [arXiv:1204.0445].
  * [16] S. F. King and C. Luhn, Neutrino Mass and Mixing with Discrete Symmetry, Rept. Prog. Phys. 76 (2013) 056201, [arXiv:1301.1340].
  * [17] S. F. King, A. Merle, S. Morisi, Y. Shimizu, and M. Tanimoto, Neutrino Mass and Mixing: from Theory to Experiment, New J. Phys. 16 (2014) 045018, [arXiv:1402.4271].
  * [18] M. Tanimoto, Neutrinos and flavor symmetries, AIP Conf. Proc. 1666 (2015), no. 1 120002.
  * [19] S. F. King, Unified Models of Neutrinos, Flavour and CP Violation, Prog. Part. Nucl. Phys. 94 (2017) 217–256, [arXiv:1701.04413].
  * [20] S. T. Petcov, Discrete Flavour Symmetries, Neutrino Mixing and Leptonic CP Violation, Eur. Phys. J. C 78 (2018), no. 9 709, [arXiv:1711.10806].
  * [21] F. Feruglio and A. Romanino, Lepton flavor symmetries, Rev. Mod. Phys. 93 (2021), no. 1 015007, [arXiv:1912.06028].
  * [22] T. Kobayashi, H. Ohki, H. Okada, Y. Shimizu, and M. Tanimoto, An Introduction to Non-Abelian Discrete Symmetries for Particle Physicists.  1, 2022.
  * [23] J. Lauer, J. Mas, and H. P. Nilles, Duality and the Role of Nonperturbative Effects on the World Sheet, Phys. Lett. B 226 (1989) 251–256.
  * [24] J. Lauer, J. Mas, and H. P. Nilles, Twisted sector representations of discrete background symmetries for two-dimensional orbifolds, Nucl. Phys. B 351 (1991) 353–424.
  * [25] S. Ferrara, . D. Lust, and S. Theisen, Target Space Modular Invariance and Low-Energy Couplings in Orbifold Compactifications, Phys. Lett. B 233 (1989) 147–152.
  * [26] A. Baur, H. P. Nilles, A. Trautner, and P. K. S. Vaudrevange, Unification of Flavor, CP, and Modular Symmetries, Phys. Lett. B 795 (2019) 7–14, [arXiv:1901.03251].
  * [27] H. P. Nilles, S. Ramos-Sánchez, and P. K. S. Vaudrevange, Eclectic Flavor Groups, JHEP 02 (2020) 045, [arXiv:2001.01736].
  * [28] H. P. Nilles, S. Ramos–Sánchez, and P. K. S. Vaudrevange, Eclectic flavor scheme from ten-dimensional string theory - II detailed technical analysis, Nucl. Phys. B 966 (2021) 115367, [arXiv:2010.13798].
  * [29] T. Kobayashi and S. Nagamoto, Zero-modes on orbifolds : magnetized orbifold models by modular transformation, Phys. Rev. D 96 (2017), no. 9 096011, [arXiv:1709.09784].
  * [30] T. Kobayashi, S. Nagamoto, S. Takada, S. Tamba, and T. H. Tatsuishi, Modular symmetry and non-Abelian discrete flavor symmetries in string compactification, Phys. Rev. D 97 (2018), no. 11 116002, [arXiv:1804.06644].
  * [31] T. Kobayashi and S. Tamba, Modular forms of finite modular subgroups from magnetized D-brane models, Phys. Rev. D 99 (2019), no. 4 046001, [arXiv:1811.11384].
  * [32] H. Ohki, S. Uemura, and R. Watanabe, Modular flavor symmetry on a magnetized torus, Phys. Rev. D 102 (2020), no. 8 085008, [arXiv:2003.04174].
  * [33] S. Kikuchi, T. Kobayashi, S. Takada, T. H. Tatsuishi, and H. Uchida, Revisiting modular symmetry in magnetized torus and orbifold compactifications, Phys. Rev. D 102 (2020), no. 10 105010, [arXiv:2005.12642].
  * [34] S. Kikuchi, T. Kobayashi, H. Otsuka, S. Takada, and H. Uchida, Modular symmetry by orbifolding magnetized $T^{2}\times T^{2}$: realization of double cover of $\Gamma_{N}$, JHEP 11 (2020) 101, [arXiv:2007.06188].
  * [35] K. Hoshiya, S. Kikuchi, T. Kobayashi, Y. Ogawa, and H. Uchida, Classification of three-generation models by orbifolding magnetized $T^{2}\times T^{2}$, PTEP 2021 (2021), no. 3 033B05, [arXiv:2012.00751].
  * [36] S. T. Petcov and M. Tanimoto, $A_{4}$ modular flavour model of quark mass hierarchies close to the fixed point $\tau=\omega$, Eur. Phys. J. C 83 (2023), no. 7 579, [arXiv:2212.13336].
  * [37] S. T. Petcov and M. Tanimoto, A4 modular flavour model of quark mass hierarchies close to the fixed point $\tau$ = i$\infty$, JHEP 08 (2023) 086, [arXiv:2306.05730].
  * [38] S. Kikuchi, T. Kobayashi, K. Nasu, S. Takada, and H. Uchida, Quark mass hierarchies and CP violation in $A_{4}\times A_{4}\times A_{4}$ modular symmetric flavor models, JHEP 07 (2023) 134, [arXiv:2302.03326].
  * [39] Y. Abe, T. Higaki, J. Kawamura, and T. Kobayashi, Quark masses and CKM hierarchies from $S_{4}^{\prime}$ modular flavor symmetry, Eur. Phys. J. C 83 (2023), no. 12 1140, [arXiv:2301.07439].
  * [40] S. Kikuchi, T. Kobayashi, K. Nasu, S. Takada, and H. Uchida, Quark hierarchical structures in modular symmetric flavor models at level 6, Phys. Rev. D 107 (2023), no. 5 055014, [arXiv:2301.03737].
  * [41] F. Feruglio, V. Gherardi, A. Romanino, and A. Titov, Modular invariant dynamics and fermion mass hierarchies around $\tau=i$, JHEP 05 (2021) 242, [arXiv:2101.08718].
  * [42] P. P. Novichkov, J. T. Penedo, and S. T. Petcov, Fermion mass hierarchies, large lepton mixing and residual modular symmetries, JHEP 04 (2021) 206, [arXiv:2102.07488].
  * [43] Y. Abe, T. Higaki, J. Kawamura, and T. Kobayashi, Quark and lepton hierarchies from S4’ modular flavor symmetry, Phys. Lett. B 842 (2023) 137977, [arXiv:2302.11183].
  * [44] Y. Abe, T. Higaki, J. Kawamura, and T. Kobayashi, Fermion hierarchies in SU(5) grand unification from ${\Gamma}_{6}^{\prime}$ modular flavor symmetry, JHEP 08 (2023) 097, [arXiv:2307.01419].
  * [45] I. de Medeiros Varzielas, M. Levy, J. T. Penedo, and S. T. Petcov, Quarks at the modular S4 cusp, JHEP 09 (2023) 196, [arXiv:2307.14410].
  * [46] S. Kikuchi, T. Kobayashi, K. Nasu, S. Takada, and H. Uchida, $Sp(6,Z)$ modular symmetry in flavor structures: quark flavor models and Siegel modular forms for $\widetilde{\Delta}(96)$, arXiv:2310.17978.
  * [47] S. Kikuchi, T. Kobayashi, and K. Nasu, CP phase in modular flavor models and discrete Froggatt-Nielsen models, arXiv:2312.11809.
  * [48] K. Ishiguro, T. Kobayashi, and H. Otsuka, Landscape of Modular Symmetric Flavor Models, JHEP 03 (2021) 161, [arXiv:2011.09154].
  * [49] K. Ishiguro, H. Okada, and H. Otsuka, Residual flavor symmetry breaking in the landscape of modular flavor models, JHEP 09 (2022) 072, [arXiv:2206.04313].
  * [50] T. Kobayashi, Y. Shimizu, K. Takagi, M. Tanimoto, and T. H. Tatsuishi, $A_{4}$ lepton flavor model and modulus stabilization from $S_{4}$ modular symmetry, Phys. Rev. D 100 (2019), no. 11 115045, [arXiv:1909.05139]. [Erratum: Phys.Rev.D 101, 039904 (2020)].
  * [51] T. Kobayashi, Y. Shimizu, K. Takagi, M. Tanimoto, T. H. Tatsuishi, and H. Uchida, $CP$ violation in modular invariant flavor models, Phys. Rev. D 101 (2020), no. 5 055046, [arXiv:1910.11553].
  * [52] Y. Abe, T. Higaki, F. Kaneko, T. Kobayashi, and H. Otsuka, Moduli inflation from modular flavor symmetries, JHEP 06 (2023) 187, [arXiv:2303.02947].
  * [53] P. P. Novichkov, J. T. Penedo, and S. T. Petcov, Modular flavour symmetries and modulus stabilisation, JHEP 03 (2022) 149, [arXiv:2201.02020].
  * [54] J. M. Leedom, N. Righi, and A. Westphal, Heterotic de Sitter beyond modular symmetry, JHEP 02 (2023) 209, [arXiv:2212.03876].
  * [55] V. Knapp-Perez, X.-G. Liu, H. P. Nilles, S. Ramos-Sanchez, and M. Ratz, Matter matters in moduli fixing and modular flavor symmetries, Phys. Lett. B 844 (2023) 138106, [arXiv:2304.14437].
  * [56] F. Feruglio, A. Strumia, and A. Titov, Modular invariance and the QCD angle, JHEP 07 (2023) 027, [arXiv:2305.08908].
  * [57] S. F. King and X. Wang, Modulus stabilisation in the multiple-modulus framework, arXiv:2310.10369.
  * [58] M. Cvetic, A. Font, L. E. Ibanez, D. Lust, and F. Quevedo, Target space duality, supersymmetry breaking and the stability of classical string vacua, Nucl. Phys. B 361 (1991) 194–232.
  * [59] T. Kobayashi, K. Nasu, R. Sakuma, and Y. Yamada, Radiative correction on moduli stabilization in modular flavor symmetric models, Phys. Rev. D 108 (2023), no. 11 115038, [arXiv:2310.15604].
  * [60] S. T. Petcov, On the Normalisation of the Modular Forms in Modular Invariant Theories of Flavour, arXiv:2311.04185.
  * [61] T. Araki, Anomaly of Discrete Symmetries and Gauge Coupling Unification, Prog. Theor. Phys. 117 (2007) 1119–1138, [hep-ph/0612306].
  * [62] T. Araki, T. Kobayashi, J. Kubo, S. Ramos-Sanchez, M. Ratz, and P. K. S. Vaudrevange, (Non-)Abelian discrete anomalies, Nucl. Phys. B 805 (2008) 124–147, [arXiv:0805.0207].
  * [63] L. E. Ibanez and D. Lust, Duality anomaly cancellation, minimal string unification and the effective low-energy Lagrangian of 4-D strings, Nucl. Phys. B 382 (1992) 305–361, [hep-th/9202046].
  * [64] J. E. Kim, Weak Interaction Singlet and Strong CP Invariance, Phys. Rev. Lett. 43 (1979) 103.
  * [65] M. A. Shifman, A. I. Vainshtein, and V. I. Zakharov, Can Confinement Ensure Natural CP Invariance of Strong Interactions?, Nucl. Phys. B 166 (1980) 493–506.
  * [66] M. Kawasaki, T. Moroi, and T. Yanagida, Can decaying particles raise the upper bound on the Peccei-Quinn scale?, Phys. Lett. B 383 (1996) 313–316, [hep-ph/9510461].
  * [67] C. Abel et al., Measurement of the Permanent Electric Dipole Moment of the Neutron, Phys. Rev. Lett. 124 (2020), no. 8 081803, [arXiv:2001.11966].
  * [68] T. Kobayashi and H. Otsuka, Common origin of the strong CP and CKM phases in string compactifications, Phys. Lett. B 807 (2020) 135554, [arXiv:2002.06931].
  * [69] Y. H. Ahn and S. K. Kang, Simple modular invariant model for quark, lepton, and flavored QCD axion, Phys. Rev. D 108 (2023), no. 9 095034, [arXiv:2306.14467].
  * [70] G. D. Coughlan, W. Fischler, E. W. Kolb, S. Raby, and G. G. Ross, Cosmological Problems for the Polonyi Potential, Phys. Lett. B 131 (1983) 59–64.
  * [71] A. S. Goncharov, A. D. Linde, and M. I. Vysotsky, COSMOLOGICAL PROBLEMS FOR SPONTANEOUSLY BROKEN SUPERGRAVITY, Phys. Lett. B 147 (1984) 279–283.
  * [72] J. R. Ellis, D. V. Nanopoulos, and M. Quiros, On the Axion, Dilaton, Polonyi, Gravitino and Shadow Matter Problems in Supergravity and Superstring Models, Phys. Lett. B 174 (1986) 176–182.
  * [73] B. de Carlos, J. A. Casas, F. Quevedo, and E. Roulet, Model independent properties and cosmological implications of the dilaton and moduli sectors of 4-d strings, Phys. Lett. B 318 (1993) 447–456, [hep-ph/9308325].
  * [74] K. Yamamoto, Saving the Axions in Superstring Models, Phys. Lett. B 161 (1985) 289–293.
  * [75] G. Lazarides, C. Panagiotakopoulos, and Q. Shafi, Baryogenesis and the Gravitino Problem in Superstring Models, Phys. Rev. Lett. 56 (1986) 557.
  * [76] D. H. Lyth and E. D. Stewart, Thermal inflation and the moduli problem, Phys. Rev. D 53 (1996) 1784–1798, [hep-ph/9510204].
  * [77] A. D. Linde, Relaxing the cosmological moduli problem, Phys. Rev. D 53 (1996) R4129–R4132, [hep-th/9601083].
  * [78] K. Nakayama, F. Takahashi, and T. T. Yanagida, On the Adiabatic Solution to the Polonyi/Moduli Problem, Phys. Rev. D 84 (2011) 123523, [arXiv:1109.2073].
  * [79] ATLAS Collaboration, G. Aad et al., Search for single vector-like B quark production and decay via B $\to$ bH($b\overline{b}$) in pp collisions at $\sqrt{s}$ = 13 TeV with the ATLAS detector, JHEP 23 (2020) 168, [arXiv:2308.02595].
  * [80] ATLAS Collaboration, G. Aad et al., Search for pair-production of vector-like quarks in pp collision events at s=13 TeV with at least one leptonically decaying Z boson and a third-generation quark with the ATLAS detector, Phys. Lett. B 843 (2023) 138019, [arXiv:2210.15413].
  * [81] ATLAS Collaboration, G. Aad et al., Search for single production of a vectorlike $T$ quark decaying into a Higgs boson and top quark with fully hadronic final states using the ATLAS detector, Phys. Rev. D 105 (2022), no. 9 092012, [arXiv:2201.07045].
  * [82] CMS Collaboration, A. M. Sirunyan et al., A search for bottom-type, vector-like quark pair production in a fully hadronic final state in proton-proton collisions at $\sqrt{s}=$ 13 TeV, Phys. Rev. D 102 (2020) 112004, [arXiv:2008.09835].
  * [83] CMS Collaboration, A. Tumasyan et al., Search for heavy resonances decaying to a pair of Lorentz-boosted Higgs bosons in final states with leptons and a bottom quark pair at $\sqrt{s}$= 13 TeV, JHEP 05 (2022) 005, [arXiv:2112.03161].
  * [84] B. Heidenreich, J. McNamara, and M. Reece, Non-standard axion electrodynamics and the dual Witten effect, arXiv:2309.07951.
  * [85] S. B. Giddings and A. Strominger, Axion Induced Topology Change in Quantum Gravity and String Theory, Nucl. Phys. B 306 (1988) 890–907.