Novel loop-diagrammatic approach to QCD $\theta$ parameter and application to the left-right model

When the QCD axion is absent in full theory, the strong $CP$ problem has to be explained by an additional mechanism, e.g., the left-right symmetry. Even though tree-level QCD $\bar\theta$ parameter is restricted by the mechanism, radiative corrections to $\bar\theta$ are mostly generated, which leads to a dangerous neutron electric dipole moment (EDM). The ordinary method for calculating the radiative $\bar\theta$ utilizes an equation $\bar \theta= - \text{arg}\, \text{det}\, m_q^{\rm loop}$ based on the chiral rotations of complex quark masses. In this paper, we point out that when full theory includes extra heavy quarks, the ordinary method is unsettled for the extra quark contributions and does not contain its full radiative corrections. We formulate a novel method to calculate the radiative corrections to $\bar\theta$ through a direct loop-diagrammatic approach, which should be more robust than the ordinary one. As an application, we investigate the radiative $\bar\theta$ in the minimal left-right symmetric model. We first confirm a seminal result that two-loop level radiative $\bar\theta$ completely vanishes (corresponding to one-loop corrections to the quark mass matrices). Furthermore, we estimate the size of a non-vanishing radiative $\bar\theta$ at three-loop level. It is found that the resultant induced neutron EDM is comparable to the current experimental bound, and the expected size is restricted by the perturbative unitarity bound in the minimal left-right model.


Introduction
The QCD θ term is P -and T -odd, and then CP -odd under the CP T invariance.Because it is identical to the total derivative, it never locally affects physics at the classical level (as long as the momentum conservation holds), while its effect occurs only via nonperturbative processes [1][2][3][4][5].It is known that this interaction induces the neutron electric dipole moment (EDM) [6,7].Measurement of the neutron EDM by the nEDM collaboration has set the severe upper bound: |d n | exp < 1.8 × 10 −26 e cm (90% CL) [8].Using the latest lattice result d n = −0.00148(34) θ e fm [9] #1 and assuming that θ is the only source of CP violation, one obtains the upper bound on the angle, | θ| 1.2 × 10 −10 (90% CL) , (1.1) # 1 The first nonzero calculation by using the lattice QCD simulation was achieved in Ref. [10], then the first statistically significant result was obtained in Ref. [9].
where θ is the physical CP -violating angle in the QCD Lagrangian which will be defined explicitly in the next section.Although the experimental bound requires that θ must be around zero, such a CPviolating phase is not restricted in the Standard Model (SM).In fact, the CP -violating phase in the Cabibbo-Kobayashi-Maskawa (CKM) matrix [11,12] is O(1); δ CKM 66 • = 1.2 rad [13] (in the standard parameterization [14,15]).If there is no trick in the full theory, θ 1, or equivalently θ δ CKM = O(1), requires a fine-tuning at O(10 −10 ) level.This is known as the strong CP problem.
The massless up quark could be a solution to the strong CP problem if the observed hadron masses are explained by the nonperturbative effect [16][17][18][19][20].However, some lattice studies ruled out this solution [21].Thus, the strong CP problem would suggest that the SM has to be extended to suppress θ without the fine-tuning.Axion is the simplest solution to the strong CP problem [22][23][24], though it suffers from another fine-tuning in the quantum gravity sector (axion quality problem) [25][26][27].
Alternatively, one may resolve the strong CP problem by extended parity symmetry [28][29][30][31].#2 In such scenarios, the extended parity involves the left-right (LR) gauge symmetry.The parity symmetry forbids the bare θ parameter, while O(1) of δ CKM is allowed.It is known that even though the bare θ parameter is strictly forbidden by the parity symmetry, radiative correction to θ is regenerated since the parity symmetry must be softly broken in nature.Eventually, one has to consider the experimental bound on the model from the neutron EDM measurements in Eq. (1.1) through the radiatively regenerated θ.
The ordinary method for calculating the radiatively generated θ parameter, adopted in many papers, utilizes (1.2) Here, m loop u,d are the up-and down-type quark mass matrices including the radiative corrections.This relation is based on the chiral rotations for the complex quark masses and an anomalous divergence of the axial-vector current, known as the Adler-Bell-Jackiw anomaly [35,36].Or, it is also derived using the path-integral formalism referred to as the Fujikawa method [37].
The ordinary method is simple though it may not be accurate.For instance, within the SM, even if one sets the bare θ parameter to be zero, it is radiatively produced via the CP -violating phase in the CKM matrix.It is shown with the above method in Ref. [38] that the contribution via the radiative corrections to quark masses is of O(G 2 F α 3 s ) (G F is the Fermi coupling constant and α s is the QCD coupling constant).On the other hand, the direct loop calculation of the correction to θ shows that it is derived at four-loop order (O(G 2 F α s )) [39].#3 It is consistent with the fact that the EDMs (and also the chromo-EDMs) of quarks are induced at three-loop order (O(G 2 F α s )) [41].In fact, the ordinary method corresponds to the diagrams where the external gluons are attached in the same #2 Other possibilities are spontaneous CP violation referred to as the Nelson-Barr mechanism [32][33][34].#3 A hadronic long-distance evaluation has confirmed that the radiative θ parameter is produced at order of O(G 2 F αs) [40].
fermion line in loop diagrams contributing to θ.It implies that the leading θ in the SM [39] comes from diagrams with external gluons attached to different fermion lines.
In this paper, we formulate a novel approach to evaluate the radiative corrections to θ through a direct loop-diagrammatic calculation, which should be more robust than the ordinary one.In Ref. [39], the external gluon field is introduced to calculate the correction to θ from the CKM matrix, while details of the technique are not written.We introduce the Fock-Swinger gauge method to directly calculate the radiative corrections to θ under the gluon field-strength background.
As an application, we investigate the radiative θ in the minimal LR symmetric model [30,31], in which the bare and one-loop level θ parameters are strictly forbidden by the LR symmetry.Although the extra heavy quarks whose Yukawa interactions violate CP symmetry are introduced, the CP -violating Yukawa interactions do not contribute to the θ parameters at one-loop level.Furthermore, it is known that two-loop level θ parameter also vanishes, which corresponds to one-loop corrections to the quark masses in the ordinary method [30,42,43].However, the ordinary method is unsettled for the extra quark contributions to θ and indeed does not contain its full radiative corrections, like the SM calculations [39].We first confirm this seminal result by using the proposed method.While new type diagrams contribute to θ at two-loop level, the sum of the diagrams still gives no contribution to θ. Next, we estimate the size of a non-vanishing radiative θ at three-loop level.It will be found that the resultant induced neutron EDM is comparable to the current experimental bound.We will also investigate a relation between the radiative three-loop level θ and the perturbative unitarity bounds of the Yukawa couplings in the minimal LR symmetric model.This paper is organized as follows.In Sec. 2, we discuss methods of direct calculation of the loop diagrams contributing to θ.We show that the operator Schwinger method and the Fock-Schwinger gauge method are applicable, though the latter method has a merit to extend the calculation to the higher-loop diagrams.In Sec. 3, the minimal LR symmetric model is briefly summarized.We also derive the parameterization by the physical parameters based on the seesaw mechanism in the LR symmetric model.We confirm that the two-loop level radiative θ vanishes by using the proposed method in Sec. 4. In Sec. 5, we investigate the numerical size of the non-vanishing radiative θ and compare both experimental (from neutron EDM) and theoretical bounds (from the perturbative unitarity bound).Section 6 is devoted to conclusions and discussion.Details of the loop calculations are given in the Appendix.

Loop-diagrammatic evaluation of QCD θ parameter
In the QCD Lagrangian, imaginary parts of the quark masses and the QCD θ term are P -odd and T -odd interactions that are not restricted from the SU (3) C gauge symmetry, where m q stands for the complex quark masses with m q ≡ |m q | exp(iθ q ), G a µν is the gluon field-strength tensor, Gaµν ≡ 1 2 µνρσ G a ρσ with 0123 = +1, and α s = g 2 s /(4π) is the SU (3) C coupling constant.It is well-known that the axial rotation of quarks turns off the imaginary part of the quark masses and generates an additional QCD θ term, where is a physical θ parameter, if all quarks are massive.This is derived using the path-integral formalism referred to as the Fujikawa method [37] or with the Adler-Bell-Jackiw anomaly [35,36,44,45] The contribution of the quark mass phase to the QCD θ term should be able to be directly evaluated by loop-diagrammatically integrating out quarks, not via the Adler-Bell-Jackiw anomaly or transformation of measure in the path integral (the Fujikawa method).However, it does not generate the QCD θ term if the momenta in the diagrams are conserved.It is because the θ term is equivalent to total-derivative and the total momentum has to be zero.Thus, we have to abandon the momentum conservation or equivalently the translation invariance in order to evaluate the QCD θ term with the loop-diagrammatic calculation.
It can be realized by introducing the gluon field strength background.In this section, we evaluate the QCD θ term with two different methods, 1) the operator Schwinger method and 2) the Fock-Schwinger gauge method.We show that they produce consistent results with Eq. (2.4).#4

Operator Schwinger method
First, we consider the operator Schwinger method.#5 The effective action ∆S induced by the integration of quarks at one-loop level is given by the log-determinant (or trace-log) of the Dirac operator.Now we introduce the complex mass parameters for quarks.In this case, the effective action is given as #4 An alternative way to evaluate the QCD θ term is the CP -odd spurion trick [46].Equation (2.4) can be derived by introducing a spurion, whose vacuum expectation value produces the CP -violating phase of the quark mass, with an external momentum injection via the spurion.#5 See Ref. [47] for the review about the operator Schwinger method.
where P µ = iD µ , D µ ≡ ∂ µ + ig s T a G a µ is the QCD covariant derivative, and P L/R = (1 ∓ γ 5 )/2.In this method, the following basic commutation relation is used, since the gluon field-strength tensor appears from it.Then, the derivative of ∆S over the fermion mass m q for P R is where Since the Levi-Civita tensor appears from the trace of a product of four γ's and γ 5 , the second order of 1  2 g s σ µν G µν leads to the G G term as where Tr(T a T b ) = (1/2)δ ab is used.Here, P 2 is replaced by −∂ 2 , and it is integrated in momentum space.Similarly, d∆S/dm * q leads to the G G term.By Integrating d∆S/dm q with m q and d∆S/dm * q with m * q , we get Gaµν . (2.10) Now we obtain the contribution from a quark with complex mass to the QCD θ term by integrating out the quark, which is consistent with the axial rotation (2.3).

Fock-Schwinger gauge method
In the previous section, we showed that the operator Schwinger method enables us to derive the physical QCD θ term by the diagrammatic evaluation.However, the operator Schwinger method is not suitable to calculate effective operators induced at higher loops because it is the method to obtain an effective action by integrating fermions out with the log-determinant of the Dirac operator.Here alternatively, we introduce the Fock-Schwinger gauge method, which is more applicable to diagrammatic calculation.#6  The Fock-Schwinger gauge is to take such a gauge #6 See Ref. [47] for the review about the Fock-Schwinger gauge method.
< l a t e x i t s h a 1 _ b a s e 6 4 = " d S J C d z J 9 j E S 8 x c 8 a q z P V 2 s P q w s r C + W l J / Z 1 T 9 I d u k t z q P W I l u g Z 1 a m B K h m 9 p w M 6 d J 4 6 o S O d 7 j B 1 o m Q x t + i v 4 e y f A v v 7 0 C M = < / l a t e x i t > Im(m q ) < l a t e x i t s h a 1 _ b a s e 6 4 = " s o B C 5 S a H 5 S J h / q 1 j x e / E Q U < l a t e x i t s h a 1 _ b a s e 6 4 = " C Z X Z 5 i v 5 T M x W 0 7 n c + e w r W d q 2 f q k Feynman amplitude iΠ q X for the loop-diagrammatic evaluation of the θ parameter.
which violates the translation symmetry.Because of breaking the translation symmetry, we can derive perturbatively the QCD θ term as will be shown below.While the Fock-Schwinger gauge violates the translation symmetry, the physical observables do not depend on it.In the below argument, we take gauge dependence parameter x 0 as x 0 = 0 for simplicity.
The gluon field G a µ can be expanded under this gauge around x = 0 and it is given with the gluon field-strength tensor at x = 0, G a µν (0), as [47] Here, the discarded terms are covariant derivatives of the background gluon field-strength tensor, which are irrelevant to the calculation of the QCD θ term.We can systematically evaluate the interaction of the propagating quarks with the background gluon field-strength tensor in this gauge fixing.However, we found that the effective gluon operators such as the QCD θ term cannot be evaluated from the simple quark bubble diagrams.The background gluon fields bring momenta, k in Eq. (2.12), which are taken to be zero in the last step of the calculation due to δ (4) (k).Thus, the quark momentum is not constant due to interaction with the background field and the quark line cannot be closed without violating momentum conservation.
In order to fix this problem, we introduce an auxiliary (dimensionless) background field X, and it is coupled to the CP -odd quark mass terms as We evaluate the leading contribution of Im(m q ) to the QCD θ term in perturbative way, assuming Im(m q ) Re(m q ).The field X is taken to be 1 in the last step of the calculation.#7  The radiative QCD θ term comes from a bubble diagram.The Feynman diagram in Fig. 1 shows the leading contribution, which is realized by integrating the delta function in Eq. (2.12) as where X(k) = d 4 xe ik•x X(x) and p is the loop momentum.We followed the Feynman rules under the Fock-Schwinger gauge, which includes the gluon field G a µ expressed as Eq.(2.12) and the modified CP -odd quark mass term in Eq. (2.13).Until integration of the delta functions, two independent momenta k 1 , k 2 flow into the vertices with the background fieldstrength tensors G a ρµ and G b σν , respectively, and the artificial background field X brings a momentum −k 1 − k 2 (see Fig. 1).After some calculation, we get Here, we take X(0) = 1.
After integrating the quark q out in the full theory, ∆L = Π q X is obtained in the effective action of the gluon.Eventually, one can derive the QCD θ term in the Fock-Schwinger gauge method.This result is consistent with that of the chiral rotation, q → q = exp(− i 2 θ q γ 5 ) q, in Eq. ( 2.3) and also the operator Schwinger method in Eq. (2.10).Hence, we reached a clarification of the equivalence among the Fock-Schwinger gauge method, the operator Schwinger method, and the ordinary chiral rotation, and we noticed that the Fock-Schwinger gauge method is more intuitive than the operator Schwinger method for higher-loop order calculations.
It might be concerned that the diagrammatic evaluations of the light quark contribution to the QCD θ term is not justified from the viewpoint of perturbation, since the loop momentum around the quark mass dominates the integrals in Eqs.(2.9) and (2.15).It might be healthy to evaluate the light quark mass phases above the Λ QCD scale and derive the QCD θ parameter by the chiral rotation.However, since the diagrammatic evaluations are consistent with those of the chiral rotation, we may forget the problem in practical cases.
In this paper, to evaluate the QCD θ term diagrammatically, we will use the Fock-Schwinger gauge method.Note that the auxiliary background field X should be attached to any perturbative interactions, but we suppress them in the following calculations for the sake of clarity.
The indices i and a represent the flavors for the doublet and singlet quarks, respectively.

Model
From this section, we introduce the minimal LR symmetric model that can solve the strong CP problem.The LR symmetry, which is formed by introducing a new SU (2) R gauge symmetry, with spatial parity symmetry is motivated to forbid the QCD θ term at tree level.In particular, we focus on the minimal LR symmetric model, which embeds the Furthermore, a SU (2) R doublet Higgs, H , and three flavors of the up-type and down-type vector-like quarks, U L , U R , D L and D R , have to be introduced.The matter contents are listed in Table 1.
To solve the strong CP problem, the spatial parity symmetry has to be extended to symmetrize the left-handed and right-handed sectors as well as the SU (2) L and SU (2) R gauge bosons, The SU (2) R gauge bosons absorb the Nambu-Goldstone (NG) bosons in the doublet H (ϕ + and ϕ 0 ) to become massive states (W + and Z ).The physical neutral Higgs boson associated with this symmetry breaking is denoted as h .Then, the SU (2) L × U (1) Y gauge symmetry is broken to U (1) EM by the VEV of H 0 , H = (0, v).The W + and Z bosons absorb the NG bosons in the doublet H (ϕ + and ϕ 0 ), and the physical (SM) neutral Higgs boson with this symmetry breaking is denoted as h.#8  The resultant parity violation in nature comes from v = v.Namely, we assume that soft parity breaking terms are contained in the H and H Higgs potentials which lead to v v = 0.These two VEVs can be chosen as real and positive without loss of generality.Since the extended parity is a discrete symmetry, its spontaneous breaking leads to the formation of the domain walls, which dominate the energy density in the Universe.This domain wall problem can be naturally solved by the Planck suppressed higher-dimensional operators, which explicitly violate the parity symmetry [43].
The Yukawa interactions and Dirac mass terms for the vector-like quarks are represented as where i = 1-3 is a flavor index for SU (2) L/R doublets, a = 1-3 is that for the singlets (vector-like quarks), and H( ) = H ( ) * ( 12 = 1).The LR symmetry requires that the Yukawa x ia u/d in the first two terms (in both lines) must be the same complex matrices, and the Dirac mass terms M u and M d must be Hermitian.The Dirac mass terms M a u and M a d in Eq. (3.2) are diagonalized to real and positive eigenvalues by the field redefinitions of the vector-like quarks.#9 The SM quark masses are realized by the seesaw mechanism such as higher dimensional operators induced by integrating out the vector-like quarks.#10  Before discussing the mass matrices in detail in the next section, let us count on the number of physical CP phases in this model.The Yukawa couplings x ia u/d are 3×3 complex matrices.Nine real parameters are removed from the Yukawa matrices by field redefinition of R with a unitary matrix U .Furthermore, phase redefinition of U a L/R and D a L/R removes five phases in the total.A remaining phase rotation corresponds to the baryon number conservation, and it does not change x u and x d .Thus, x u and x d have a total of 22 physical real parameters.We parametrize these 22 parameters as The h-h and Z-Z mixings are induced in the model at tree level, though the mixings are suppressed by v/v and (v/v ) 2 , respectively [30].Since we take v/v → 0 in the calculation of the θ parameter, they are ignored.#9 One can also consider non-Hermitian vector-like quark mass matrices which correspond to soft parity breaking terms [30].However, such contributions produce large quark EDM and radiative θ, and thus they are severely constrained from the EDM bounds [43,49].#10 One can also extend the lepton sector that is insensitive to the QCD θ term.If one considers SU (5)L × SU (5)R grand unification [50][51][52], vector-like neutral leptons are absent and the neutrinos keep massless at the tree level.Interestingly, suitable Dirac neutrino masses are generated from the two-loop radiative corrections [53] with predicting a nonzero ∆N eff [54].Furthermore, an O(10) keV sterile neutrino dark matter with the leptogenesis mechanism can be incorporated [55,56].
with Φ(θ 3 , θ 8 ) ≡ exp(iτ 3 θ 3 ) exp(iτ 8 θ 8 ) , (3.4) where τ 3 and τ 8 are the third and eighth Gell-Mann matrices.Here, xu/d are real diagonal matrices and V Q , V U , and V D are CKM-like unitary matrices which have three rotation angles and one CP -violating phase.It is found that there are seven CP -violating phases in this model (θ u3/u8 , θ d3/d8 , and three phases in V Q/U/D ).When the Dirac masses are assumed to be universal such as M a u = M u and M a d = M d for a = 1-3, the parameters θ u3/u8 , θ d3/d8 , V U/D become unphysical and only V Q remains physical.
In this paper, we assume that v < ∼ M a q and will utilize expansions by v /M a q .This inequality is motivated because the seesaw mechanism may explain the SM fermion mass hierarchy naturally.On the other hand, if v M a q , a copy of the SM fermions has a mass spectrum similar to the SM fermions, which spread over five orders of magnitude.A new naturalness problem might appear in such a model, but the QCD θ term is suppressed by M a q /v since only the CKM phase survives in a limit of M a q → 0. In the following, we will consider the case of v < ∼ M a q .

Parametrization of Yukawa coupling constants
In this section, we show the quark mass matrices and define the mass eigenstates.In the mass matrices the Yukawa coupling constants, x u and x d , appear, and we have to determine them from the observed quark masses and CKM matrix in order to evaluate the radiative θ parameter.We give the parameterization of the Yukawa coupling constants assuming the SM quark masses are given by the seesaw mechanism with v M a q .From Eq. (3.2), the quark mass matrices in the flavor eigenstates are given as where ) are the up-and down-type flavor eigenstates, and v 174.1 GeV.Here, M a u/d are real diagonal, while x ia u/d are complex matrices.It is obvious that arg det M (0) u/d = 0.Then, the 6×6 fermion mass matrices are diagonalized by bi-unitary matrices, V qL and V qR , as with diagonal mass matrices Mq .The mass eigenstates, U P M L/R and D P M L/R for P = 1-6, are given as It is difficult to reconstruct the model parameters from the experimental data in general.Here we assume that v < ∼ M a q and we take leading terms in the expansion of v /M a q for the quark mass eigenvalues.In this expansion, the SM quark masses are given by the following seesaw relation with a 3 × 3 unitary matrix V q and I = 1-3, while the heavy quark masses are given by M a q , for q = u and d.Now let us rewrite Eq. (3.8) as where I 3 is a unit matrix in the three-dimensional space and we ignore the indices of matrices.#11 Thus, the Yukawa matrices x q are given with a unitary matrix U q by According to the previous section, one can remove some unphysical parameters in U q and V q by the field redefinitions.Then, we get where V CKM (≡ V u V † d ) corresponds to the CKM matrix in the SM.Here, V U/D are CKMlike unitary matrices with three mixing angles and one CP -violating phase, though they are different from those in Eq. (3.3).Now we have seven physical CP -violating phases (θ u3/u8 , θ d3/d8 , and three phases in V U/D and V CKM ), which is consistent with our previous counting, and all phases can be O(1) under the extended parity symmetry.
Since we assume that v < ∼ M a q , the 6 × 6 diagonalization matrices in Eq. (3.7) are given of leading terms in the expansion of v /M a q as where x q are given by Eq. (3.11).Here, the diagonal eigenvalue matrices Mq are given as, M P q = diag(m I q , M a q ) for q = u and d . (3.13) #11 A similar parameterization technique, referred to as the Casas-Ibarra parameterization, is applied in the minimal seesaw model [57,58].

Quark EDMs
Before discussing the radiative θ parameter in the minimal LR symmetric model, let us comment on contributions to quark (chromo) EDMs.As long as the vector-like mass matrices are Hermitian, the quark EDMs vanish completely at one-loop level.The W ( )± contribution at one-loop level vanishes trivially since the chirality is conserved in the diagrams.On the other hand, the one-loop quark EDM contributions from neutral Higgses and Z ( ) may have a chirality flip in the diagrams.Nevertheless, they also vanish because the extended parity symmetry restricts the product of two vertices to be strictly real, as shown in Ref. [43].

Confirmation of vanishing QCD θ parameter in two-loop order
The parity symmetry is spontaneously broken by the H ( H ) in the LR symmetric model.Since arg det M (0) u/d = 0 holds, fermion one-loop contributions to the QCD θ term remain zero.However, it is expected that fermion-loop diagrams at a higher than one-loop level would generate it.In this section, we show fermion two-loop contributions to the QCD θ term still vanish.
Integrating out quarks, the following higher-dimensional operators are expected to be generated, with M q as the scale of vector-like quark masses, and the QCD θ term are induced by the spontaneous parity symmetry breaking We evaluate the Wilson coefficients of the operators C n in the following.First, we consider the contributions to the QCD θ term at two-loop level coming from an exchange of the W ± boson.The two-loop fermion bubble diagrams mediated by the W ± boson under the gluon field-strength background conserve chirality in the fermion line, and then it is proportional to where i and j run 1-3 as the flavor index for the SU (2) R doublet, while P and Q run 1-6 for the quark mass eigenstates.Here, a two-loop function f dR and it corresponds to Eq. (4.3) by an exchange of i ↔ j, the contribution is real so that it does not generate the QCD θ term.
< l a t e x i t s h a 1 _ b a s e 6 4 = " E 5 p n G f J C d N g U 7 u l s b y E X 8 y n h P 0 0 r 2 s / / q t y t H w = = < / l a t e x i t > U M < l a t e x i t s h a 1 _ b a s e 6 4 = " A 4 r y 8 w F n M m 6 7 9 + 4 I n j a 1 7 1 d r 9 6 v L 6 c n n t Y f 6 6 Z + k G 3 a Q l 1 F q h N X p M d W q g y g G 9 o U P 6 4 D x x t J M 6 2 T h 1 p p B j r t J f w 3 n x C 2 w 3 0 y w = < / l a t e x i t > ' (0)

< l a t e x i t s h a 1 _ b a s e 6 4 = " H M e B w E M E H g D 1 x w a A f z 4 c u c q W E C A = " >
w F n M m 6 7 9 + 4 I n j a 1 7 1 d r 9 6 v L 6 c n n t Y f 6 6 Z + k G 3 a Q l 1 F q h N X p M d W q g y g G 9 o U P 6 4 D x x t J M 6 2 T h 1 p p B j r t J f w 3 n x C 2 w 3 0 y w = < / l a t e x i t > ' (0)

< l a t e x i t s h a 1 _ b a s e 6 4 = " H M e B w E M E H g D 1 x w a A f z 4 c u c q W E C A = " >
w F n M m 6 7 9 + 4 I n j a 1 7 1 d r 9 6 v L 6 c n n t Y f 6 6 Z + k G 3 a Q l 1 F q h N X p M d W q g y g G 9 o U P 6 4 D x x t J M 6 2 T h 1 p p B j r t J f w 3 n x C 2 w 3 0 y w = < / l a t e x i t > The charged NG boson contributions to the QCD θ term at two-loop level.
The reason why the exchange of the W ± boson does not contribute to the QCD θ term at two-loop level is clear.However, the above discussion is based on the structure of the mixing matrices in the contribution, not on the structure of Lagrangian parameters, such as x u/d and M u/d , and then, it is unclear what is required to generate the QCD θ term in higher-order diagrams.We make it clear by explicit calculation of the loop diagrams in the following.
In the unitary gauge, the lowest dimension operator (n = 1) in Eq. ( 4.1) might come from diagrams which include the longitudinal mode of the W ± boson.It is because the propagator is proportional to k µ k ν /m 2 W (k ν the momentum of the W ± boson), and it could give the lowest order contribution with regard to v 2 .The Yukawa coupling constants x u and x d are multiplied with the Higgs VEVs in the mixing matrices as in Eq. (3.12).
In our calculation, we adopt the R ξ gauge with the Feynman-'t Hooft gauge ξ = 1 (for SU (2) L × SU (2) R × U (1) B−L gauge), in order to avoid the messy calculation in the unitary gauge.The lowest dimension operator (n = 1) in Eq. (4.1) could arise the charged NG boson exchange ϕ ± in this gauge.The charged NG boson is absorbed by W ± boson in the Higgs mechanism, and its mass, m ϕ , is equal to the W ± mass.The charged NG boson interactions are given as 4) Both the left-and right-handed quarks are coupled with the charged NG boson.We will show that the charged NG boson diagrams at two-loop level do not contribute to the QCD θ term.
Three diagrams (dubbed as diagrams A, B and C) in Fig. 2 could give contributions to the θ parameter.By using the Fock-Schwinger gauge method (for SU (3) C gauge) in Sec.2.2, we obtain ) ) where the two-loop functions Ī(1;3) , Ī(3;1) , and I (2;2) are defined in Appendix A. In the above evaluation, we pick up contributions proportional to both x u and x d , which are also proportional to the quark masses in the mass eigenstate propagators.The terms proportional to x u and x * u (or x d and x * d ) is real.Here, we use the dimensional regularization (d = 4 − 2 ) for loop momentum integrals and the partial diagrams produce UV divergence.However, the contributions from diagrams A and B, proportional to Ī (1;3) and Ī (3;1) in Eq. (A.13), respectively, vanish, so that the correction to the θ parameter is finite and scale-independent.The UV divergent parts (1/ terms) in Ī (1;3) and Ī (3;1) cancel out since To derive the above equations we use the following equations, for q = u and d , ( in addition to Eq. (3.6) with ( M P q ) * = M P q .Furthermore, we observed that the following combinations also vanish #12   Im Im Therefore, the second terms of Ī (1;3) and Ī (3;1) in Eq. (A.13) also do not affect the θ parameter.On the other hand, the loop function I (2;2) in Eq. (4.7) is UV finite.Similar to the ϕ ± contribution, the contribution from the SM charged NG boson ϕ ± , absorbed into W ± , is derived by replacing L(R) with R(L), and m 2 ϕ with m 2 ϕ in the above ) and (4.12) can be perturbatively proved by assuming the off-diagonal terms in the quark mass matrices are small.The similar trick is also used around Eq. (4.26).We also checked Eqs.(4.11) and (4.12) numerically.
formulae.Furthermore, (−1) is multiplied since the chiralities of circulating fermions are opposite to diagrams of Fig. 2. It means that, if one sets v = v corresponding to the LR symmetric limit, those two contributions of ϕ ± and ϕ ± cancel each other.The diagrams A and B correspond to the one-loop correction to the fermion mass terms.The two-loop function Ī(3;1) (x 1 ; x 2 ; x 3 ) is expressed as where F 0 (p 2 , x 2 , x 3 ) is a loop function of one-loop diagrams for the fermion mass correction, with Q 2 ≡ 4πµ 2 e −γ E and µ is the renormalization scale.( Ī(1;3) (x 1 ; x 2 ; x 3 ) also has a similar expression, see Appendix A.) Ī(3;1) (x 1 ; x 2 ; x 3 ) has an IR-singular behavior when x 1 x 3 as while small x 2 and x 3 do not lead to IR singularities.This behavior is expected.It is because if a fermion with real mass m f gets a constant radiative correction to the fermion mass m f + δm f , the correction to the θ parameter is given by δθ −Im(δm f )/m f , see Eq. (2.15).#13 However, this evaluation of δθ is justified only when the correction to the fermion mass term is independent of fermion momentum.
The IR-singular behaviors of the SM fermion masses in Ī(3;1) and Ī(1;3) are not physical in δθ| A and δθ| B in Eqs.(4.5) and (4.6), and they can be removed indeed using Eq.(4.10) as ) where A and B run 4-6 as the heavy quark mass eigenstates, see Eq. (3.13).Here, the SM quark masses in the loop function are taken to be zero.On the other hand, the contribution of diagram C is not associated with the correction to the quark masses, and then it is a new type contribution to the θ parameter.It is suppressed by the heavier fermion or ϕ ± masses.By taking the SM quark masses to be zero in the loop function (see Appendix A), it is given as where A and B run 4-6 as the heavy quark mass eigenstates.It is found that the diagram C does not contain any IR-singular behavior, unlike the diagrams A and B.
When v < ∼ M a q , the leading contributions of O(v 2 /(M a q ) 2 ) are given as with a Hermitian matrix A a q , (A a q ) ij ≡ x ia q x * ja q for q = u, d and not sum a index .(4.22)This can be derived from the above formulae by It is found that these radiative corrections to the θ parameter vanish.For example, δθ| A is given as where is the real function, and the Hermitian property of the matrix A a q is used.The same conclusions are applicable to δθ| B and δθ| C at this order.Now we showed that the charged NG boson contribution to the θ parameter at twoloop level vanishes at the leading order of v (n = 1 in Eq. (4.1)).It comes from the fact that the contributions are proportional to the fourth power of x u/d .We have also checked that the contributions of the sixth power of x u/d , corresponding to O(v 4 /(M a q ) 4 ) contributions, also vanish.The contributions are derived from the above formulae with the mass-insertion approximation, ) Each first term is the aforementioned leading contribution, which vanishes (see Eq. (4.26)).The next-to-leading contributions are These loop functions, which come from the mass-insertion approximation, are given in Appendix A. The sequences of masses connected by commas in the loop functions are introduced by the mass insertion.Again, it is found that these contributions are zero.For example, the first term in δθ| A is given as Here, above two real loop functions f and g are symmetric under exchanges of ( , respectively.Then, the above equation vanishes, see Appendix A. The symmetry comes from the mass-insertion approximation.Even if we include the higher-order contributions of x q in the mass-insertion approximation, they still vanish since the loop function is real and symmetric for the exchange of the heavy fermion masses.Now we found that the charged NG boson does not give a contribution to the θ parameter at two-loop level.We also numerically checked this fact by using Eqs.(4.16)- (4.18).The W ± contributions to the θ parameter in the Feynman-'t Hooft gauge at two-loop level vanish.The Yukawa coupling dependence comes from only the mixing matrices, and then the leading contributions, which are proportional to the fourth power of x u/d at most, vanish.The higher-order contributions, coming from mass-insertion approximation, also vanish due to the symmetry of heavy fermion masses in loop functions.
A similar discussion is applicable for the other contribution, such as Z , h , and ϕ 0 at two-loop level.Then, we confirmed the two-loop contribution to the θ parameter vanishes as far as M q > ∼ H H .

Non-vanishing contribution to QCD θ parameter in three-loop order
In the previous section, we confirmed that the QCD θ term is not generated in the twoloop level contribution, i.e., up to the fourth order of the Yukawa interaction x q .We also found that it is valid even if one considers the higher-order contributions of x q by using the mass-insertion approximation.In order to give non-vanishing contributions to the θ parameter, the commutation relation [A a q , A b q ] must be nonzero, see Eq. (4.33).It implies that non-vanishing contribution should be proportional to Im Tr(A a q [A b q , A c q ]) for q, q = u and/or d rather than Im Tr([A b q , A c q ]), and the loop function has to be asymmetric under exchange between (M b q ) 2 and (M c q ) 2 .Thus, the contributions of the following form might  4 from the current neutron EDM measurement, under an assumption of f 13 duu = 1.The numbers in parentheses are those in Fig. 5 under an assumption of f 13 duu = f 23 duu = 1.The parameter sets are restricted by the perturbative bound of the Yukawa coupling x u as 4π or √ 4π. Here, the SM (running) quark masses at µ = 1 TeV.#15 Furthermore, by assuming f 13 duu = 1, the vertical line in the figures stands for the experimental upper bound from the neutron EDM measurement in Eq. (1.1) and the right area is excluded at 90% CL.
Given the fixed values of Φ(θ u3 , θ u8 )V U , one can obtain the eigenvalues of the up-type Yukawa matrix x u .In the PDF analysis, we find that the eigenvalues of the Yukawa matrix x u can be larger than O(1) easily depending on Φ(θ u3 , θ u8 )V U .Therefore, we impose the maximal eigenvalue to be smaller than 4π or √ 4π for the dashed or dotted lines, respectively, as the perturbative unitarity bound.When the maximal eigenvalue exceeds the unitary bound, we discarded these points in the PDFs.In order to display the reduction in statistics as a result of setting the perturbative bound, we do not normalize the dashed and dotted line PDFs by 1, and hence their total integration is less than 1.
The ratios of the excluded parameter regions by the neutron EDM measurement in each of the PDFs are shown in Table 2 under an assumption of f 13 duu = 1.Here, the ratios of the excluded parameter regions are obtained by comparing the areas of the dashed or dotted line PDFs; right area of the vertical line over the total area where the unitarity bound is imposed.We found that some fractions of parameter regions have already been excluded even if the perturbative bound is imposed.In the case of the perturbative bound of 4π to the eigenvalues of x u matrix, about 30-40% of the whole parameter region is already excluded.On the other hand, in the case of √ 4π, the dependence of M 3 u appears obviously.In particular, only 3.65% is excluded in the case of M 1 u = 10 3 M 3 u .These results show that this model is sensitive to the current bound from the neutron EDM experiments and has a possibility to be explored by future improvement of the experiments.
Moreover, in Figs.5a and 5b, we show the PDFs of the absolute value of δθ duu with (b, c) = (1, 3) plus (2, 3) (also included (3, 1) and (3, 2)), with assuming M 1 u = M 2 u which should be a somewhat aggressive parameter choice in light of m u m c .We observed that the PDFs in Fig. 5 are slightly larger than the PDFs in Fig. 4. The ratios of excluded parameter regions by the neutron EDM measurement are shown in Table 2 at the numbers in parentheses.
Above estimation of the leading contribution to the radiative θ parameter is numerically consistent with the latest analysis Eq. (28) of Ref. [49], where the ordinary calculation #15 We take mt = 143 GeV, m b = 2.41 GeV, mc = 528 MeV, and ms = 45.4MeV at µ = 1 TeV [60].method is used.In this paper, we showed for the first time that an effect from the perturbative unitarity bound is important and it reduces the radiative corrections to θ.

Case for the GIM by universal vector-like mass
Next, we investigate a case that all Dirac quark masses are degenerate as M = M a u = M a d .In this case, the GIM-like mechanism occurs and then the CKM matrix becomes a unique source of the CP -violating phase.As a consequence, the radiative θ parameter would be significantly suppressed, and it is expected as the minimum value of the θ parameter in the minimal LR symmetric model.The induced θ parameter should be proportional to the Jarlskog invariant of the CKM matrix J CKM , which is given by Im [61] and J CKM = (3.08 +0.15 −0.13 )×10 −5 [15].If the θ parameter is induced at three-loop level, the contribution would be given as where f CKM is a dimensionless loop function of O (1).It is assumed that in the three-loop diagrams, two scalar bosons are exchanged inside the fermion loop, and the other scalar lines are replaced by v .The above contribution is proportional to v 8 .The perturbativity of the top Yukawa requires M v .Then, δθ minimum is expected as at most 10 −19 or smaller when all vector-like quark masses are degenerate.This size is comparable to or smaller than the CKM phase contribution to the θ parameter in the SM at four-loop level, evaluated in Ref. [39] (though the top quark had been assumed to be lighter than the W boson).This is because the quark mass suppression of the contribution in the SM is milder than Eq.(5.4).In addition, the neutron EDM induced by the CKM phase via longdistance hadronic contributions in the SM [62,63] is much larger than by the contribution to the θ parameter in the above benchmark point.Then, a more precise evaluation of the θ parameter in the benchmark point is too academic and beyond our scope.

Conclusions and discussion
When the QCD axion is absent, one has to solve the strong CP problem by an additional discrete symmetry.The extended parity with the LR gauge symmetry can solve it with generating the SM as the low-energy theory.However, it is known that the radiative θ parameter is induced from the soft symmetry breaking of the parity.In this paper, we first formulated a novel method of direct loop-diagrammatic calculation of the radiative θ parameter by using Fock-Schwinger gauge.This approach should be more robust than the the ordinary calculation method based on the chiral rotations.By using the Fock-Schwinger gauge method, we confirmed a seminal result that two-loop level θ vanishes completely in the minimal LR symmetric model.Furthermore, we estimated the size of the leading contributions to the non-vanishing radiative θ parameter at three-loop level.We derive the parameterization by the physical parameters based on the seesaw mechanism in the LR symmetric model, and we obtained the probability density functions of the radiative θ parameter by varying all free but physical parameters.Here, we also investigated the impact of the perturbative unitarity conditions for the LR symmetric Yukawa matrices.It is found that the resultant θ parameters are partially excluded by the current neutron EDM bound.It implies that this model has a possibility to be explored by future improvement of the experiments.One should note that the minimal LR symmetric model predicts that all hadronic EDMs are dominated by the radiative θ parameter.Therefore, this model can predict a distinctive and non-vanishing correlation between the neutron, proton, nuclei ( 2 H, 3 He), diamagnetic atoms (Hg, Ra), paramagnetic atoms and molecules (YbF, HfF, ThO) EDMs [49,[64][65][66].
The large mass-scale difference between v and v can be explained by a Higgs parity mechanism that predicts λ SM (µ = v ) 0 and v = O(10 10 ) GeV [42,[67][68][69].Since the above estimation of the radiative θ parameter is insensitive to the mass scale itself of the right-handed sector, one could predict the neutron EDM for the case of v = O(10 10 ) GeV.
We also comment on the radiative θ parameter in the spontaneous CP violation (Nelson-Barr) model [32][33][34], which is another model that can explain the strong CP problem by the discrete symmetry.In Ref. [70], the radiative θ parameter has been investigated in detail.It is found that although the reducible θ, which comes from quartic couplings of scalar fields responsible for the spontaneous CP violation with the SM Higgs one, is induced at two-loop level, it can be numerically neglected if the couplings are small enough.On the other hand, the irreducible θ, which is related to the CKM phase, is induced at threeloop level and is safely below the current experimental bounds.The loop-diagrammatic approach we proposed would provide a more robust estimation if possible.
It would be an interesting direction to investigate correlations between the radiative θ parameter and other (flavor) observables in the minimal LR symmetric model: The lepton flavor universality violation in B → D ( * ) lν (R(D ( * ) ) anomaly) (the recent review [71,72]), the Cabibbo angle anomaly (the recent review [73]) and the W -boson mass anomaly [74] could be explained [75][76][77].Furthermore, the investigation of a correlation with electroweak-like baryogenesis in the right-handed sector should be an attractive prospect [78,79].SPRING, Grant Number JPMJSP2125.The author (A.Y.) would like to take this opportunity to thank the "Interdisciplinary Frontier Next-Generation Researcher Program of the Tokai Higher Education and Research System."

A Loop functions
We define the loop functions in this section.The two-loop functions used in this paper are given as where the dimensional regularization is used on d = 4 − 2 dimension and µ is the renormalization scale.The functions can also be derived by derivative or finite difference of I(x 1 ; x 2 ; x 3 ) (≡ I (1;1) (x 1 ; x 2 ; x 3 )) as I (n;m) (x 1 ; x 2 ; x 3 ) = 1 (n − 1)!(m − 1)!
B H b P P 7 C P 7 y n 7 O 5 c o t h 1 E 7 x u p P s S J t r 7 y 4 t v P j v 6 g I q 1 H 9G 3 W q Z k 1 d n M 1 o l d C e W o 8 5 R T D F D 5 8 f H u 9 s b p f z m + w d + w b 9 b 9 k R + 4 Q T x M P v w f u 6 2 H 5 z i p 4 T J R o 3 Y O 5 B g b M M V u P R y N 6 k 2 5 g B b I 6 K F X t P P X h S c J h 1 F j + / l n I d 7 m P X c Z 2 M Y Y 9 s D y L L F S O S w z / N n d i 3 0 P r D c 4 L O q Y T I x L F w x E w P R 4 h I 5 K Q 2 d 4 y + Z M j m 0 K A c 0 3 x 1 J t d k C p xl g j d d / f c F z x q 7 6 5 X q n c p G f a O 0 9 c C 9 7 k W 6 T j f o F m r d p S 1 6 R D V q o I q g l / S K X n v 3 v b b X 9 f r T 1 I W C w 1 y l v 4 b 3 9 B d a z c t C < / l a t e x i t > X < l a t e x i t s h a 1 _ b a s e 6 4 = "

1 < l a t e x i t s h a 1 _ b a s e 6 4 =
w e N 9 a 3 q s U t 9 o 5 9 g / 6 3 7 I h 9 w g m S 4 f f w / a b Y e n O K n h M l G j d g 7 k G B s w p W 4 9 H I X q c 7 m C F s j o o 1 e 0 9 d e D J w m H U WP 7 + W c h 3 u Y d d x n U x g j 2 w P Y s u V I F L A P 8 2 d 2 L f Q + s N z g i 6 o g s j E s X D E T A 9 H i E j k Z D Z 3 j L 7 k y O b Q o B z T f H U m 1 2 Q K n G W C N + 3 / + 4 Jn j e 2 7 N f 9 e b W 1 z r b L x y L 3 u R b p B N + k 2 a t 2 n D X p C d W q i S p d e 0 W s 6 9 B 5 6 3 O t 5 B 9 P U h Z L D X K O / h q d + A c d h y / k = < / l a t e x i t > k " 8 N M S c S t A S s o h e s Q 1 6 b 5 7 L z n W m D w 6 H f w e 7 P e e N s T d V l e U P o C P E l 4 G 7 S G X 2 h X 1 g p + w z + 8 i + s p 8 z u X L L Y d S O s P o T r E h b i 2 9 u 1 3 / 8 F x V h N a p / o 8 7 V r G k f Z z N a J b S n 1 m N O E U z w g 9 e H p / W 1 r X J + j 7 1 j 3 6 D / m J 2 w T z h B P P g e v N 8 U W 0 f n 6 D l T o n E D 5 h 4 2 b W W G b e 4 u z W s z u J 0 y 8 C 7 O l b x w c 3 a 5 M L 4 9 E o + J M f G f + U 3 E j v n A H b u P O P F + S y 8 f P 5 A n T U 6 9 y V k l P 1 y X e T c 7 F S Z K 6 r E T k W 2 u b y f 8 p / c Y 8 q C P k W G n G L r w S h b 8 v w F O w N p 4 v v M t P L E 3 k Z u b S y 9 G F I Q x j j D 6 T m M E 8 F l F k k v c 4 w g l O t U P t U v u k X b W W a m 2 p 5 i 3 + a N r n X 3 7 3 r Q 4 = < / l a t e x i t > D M < l a t e x i t s h a 1 _ b a s e 6 4 = " A 4 r y 8 Q T 2 b W W G b e 4 u z W s z u J 0 y 8 C 7 O l b x w c 3 a 5 M L 4 9 E o + J M f G f + U 3 E j v n A H b u P O P F + S y 8 f P 5 A n T U 6 9 y V k l P 1 y X e T c 7 F S Z K 6 r E T k W 2 u b y f 8 p / c Y 8 q C P k W G n G L r w S h b 8 v w F O w N p 4 v v M t P L E 3 k Z u b S y 9 G F I Q x j j D 6 T m M E 8 F l F k k v c 4 w g l O t U P t U v u k X b W W a m 2 p 5 i 3 + a N r n X 3 7 3 r Q 4 = < / l a t e x i t > D M < l a t e x i t s h a 1 _ b a s e 6 4 = " E 5 p n G f J C d N g U 7 u l s b y E X 8 y n h P 0 0 = " > A A A D A n i c h V I 7 T 9 x A G B x M e I b H Q Z p I a U 4 5 g a h O e x E K i A o p F D S R e B 0 g w e l k + x b O w 4 l w c i 0 v x 7 0 m v K P G I 0 2 5 z N J p a 6 Z f 7 v 7 6 d / / u i y u E Y p 7 5V P Z t Z Y Y 1 7 i 7 N a z O 4 n T L w L s 6 m v 7 + x d z 0 / M D U X D 4 k B c M f + + u B B n 3 I F b / 2 M e z s q 5 b 8 / k C d N T r 3 J W S U / X J d 5 K z s V J k r q s R O S b a x v J / y n d Y W 7 U E X K s N G I X X o n C w w v w G C x + y B c + 5 k d n R 3 O T U + n l 6 M Q 7 v M c I f c Y w i W n M o M g k m 9 j D d + x r u 9 o P 7 U Q 7 b S 7 V W l L N G 9 xr 2 s / / q t y t H w = = < / l a t e x i t > U M < l a t e x i t s h a 1 _ b a s e 6 4 = " E 5 p n G f J C d N g U 7 u l s b y E X 8 y n h P 0 0 = " > A A A D A n i c h V I 7 T 9 x A G B x M e I b H Q Z p I a U 4 5 g a h O e x E K i A o p F D S R e B 0 g w e l k + x b O w 4 l w c i 0 v x 7 0 m v K P G I 0 2 5 z N J p a 6 Z f 7 v 7 6 d / / u i y u E Y p 7 5V P Z t Z Y Y 1 7 i 7 N a z O 4 n T L w L s 6 m v 7 + x d z 0 / M D U X D 4 k B c M f + + u B B n 3 I F b / 2 M e z s q 5 b 8 / k C d N T r 3 J W S U / X J d 5 K z s V J k r q s R O S b a x v J / y n d Y W 7 U E X K s N G I X X o n C w w v w G C x + y B c + 5 k d n R 3 O T U + n l 6 M Q 7 v M c I f c Y w i W n M o M g k m 9 j D d + x r u 9 o P 7 U Q 7 b S 7 V W l L N G 9 xr 2 s / / q t y t H w = = < / l a t e x i t > U M < l a t e x i t s h a 1 _ b a s e 6 4 = " A 4 r y 8 Q T 2 b W W G b e 4 u z W s z u J 0 y 8 C 7 O l b x w c 3 a 5 M L 4 9 E o + J M f G f + U 3 E j v n A H b u P O P F + S y 8 f P 5 A n T U 6 9 y V k l P 1 y X e T c 7 F S Z K 6 r E T k W 2 u b y f 8 p / c Y 8 q C P k W G n G L r w S h b 8 v w F O w N p 4 v v M t P L E 3 k Z u b S y 9 G F I Q x j j D 6 T m M E 8 F l F k k v c 4 w g l O t U P t U v u k X b W W a m 2 p 5 i 3 + a N r n X 3 7 3 r Q 4 = < / l a t e x i t > D M < l a t e x i t s h a 1 _ b a s e 6 4 = " H M e B w E M E H g D 1 x w a A f z 4 c u c q W

# 12
The factors 1/M log M (M : quark mass) in Eqs.(4.11) and (4.12) correspond to the O( ) term in Eq. (2.15) when changing d 4 p to d d p (d = 4 − 2 ).They become O( 0 ) since 1/ comes from the quark self-energy subdiagrams in diagrams A and B. Equations (4.11

Table 1 :
The matter contents and their gauge charges in the minimal LR symmetric model, where U (1) Y

Table 2 :
Ratios of excluded parameter regions for the parameter sets in Fig.