Dynamical Generation of Elementary Fermion Mass: First Lattice Evidence

TL;DR

This paper presents lattice simulation evidence for a non-perturbative mechanism that dynamically generates elementary fermion masses and weak boson masses in gauge-fermion-scalar models, suggesting new physics beyond the Standard Model.

Contribution

First lattice evidence demonstrating a non-perturbative mass generation mechanism for elementary particles in gauge-fermion-scalar models.

Findings

Fermion mass is dynamically generated in the Nambu-Goldstone phase.

Weak boson masses are generated in models with electro-weak interactions.

A non-perturbative mechanism for mass generation is demonstrated from first principles.

Abstract

Using lattice simulations we demonstrate from first principles the existence of a non-perturbative mechanism for elementary particle mass generation in models with gauge fields, fermions and scalars, if an exact invariance forbids power divergent fermion masses and fermionic chiral symmetries broken at UV scale are maximally restored. We show that in the Nambu-Goldstone phase a fermion mass term, unrelated to the Yukawa operator, is dynamically generated. In models with electro-weak interactions weak boson masses are also generated opening new scenarios for beyond the Standard Model physics.

Figures

Tables

Table 1. Table 1: Gauge and scalar bare parameters β 𝛽 \beta , m ϕ 2 superscript subscript 𝑚 italic-ϕ 2 m_{\phi}^{2} , λ 0 subscript 𝜆 0 \lambda_{0} and the corresponding measured value of the renormalized quantities r 0 2 ​ M σ 2 superscript subscript 𝑟 0 2 superscript subscript 𝑀 𝜎 2 r_{0}^{2}M_{\sigma}^{2} and λ R subscript 𝜆 𝑅 \lambda_{R} for the three lattice spacings we have considered. The value of the critical mass b 2 ​ m c ​ r 2 superscript 𝑏 2 superscript subscript 𝑚 𝑐 𝑟 2 b^{2}m_{cr}^{2} which defines the phase transition point is also given.

$b$/fm	$b^2{m}_{cr}^2$	$r_0^2{M}_{\sigma }^2$	${\lambda }_R$	$\beta$	$b^2{m}_{\varphi }^2$	${\lambda }_0$
0.152	-0.5269(12)	1.278(6)	0.437(4)	5.75	-0.5941	0.5807
0.123	-0.5357(11)	1.286(6)	0.441(4)	5.85	-0.5805	0.5917
0.102	-0.5460(10)	1.290(7)	0.444(5)	5.95	-0.5756	0.6022

Table 2. Table 2: Simulation parameters in the Wigner phase. We used 480 gauge and scalar paired configurations at each simulation point. Each of these bosonic configurations is different from any other by either the scalar field or both the scalar and the gauge ones.

$b$/fm	$\left(L^3\times T\right)/b^4$	$\eta$	$b\mu$
0.152	${16}^3\times 32$	-1.1505	0.0180 0.0280 0.0480
($\beta =5.75$)		-1.1898	0.0180 0.0280 0.0480
		-1.3668	0.0180 0.0280 0.0480
0.123	${16}^3\times 40$	-1.0983	0.0224 0.0316 0.0387
($\beta =5.85$)		-1.1375	0.0120 0.0224 0.0316 0.0387
		-1.2944	0.0224 0.0387
0.102	${20}^3\times 48$	-0.9761	0.0186, 0.0321
($\beta =5.95$)		-1.0354	0.0186, 0.0321
		-1.0771	0.0100 0.0186, 0.0321

Table 3. Table 3: The values of η c ​ r subscript 𝜂 𝑐 𝑟 \eta_{cr} and Z P − 1 superscript subscript 𝑍 𝑃 1 Z_{P}^{-1} at our three lattice spacings for ρ = 1.96 𝜌 1.96 \rho=1.96 , the listed values of ρ 𝜌 \rho and λ 0 subscript 𝜆 0 \lambda_{0} and b 2 ​ m ϕ 2 superscript 𝑏 2 superscript subscript 𝑚 italic-ϕ 2 b^{2}m_{\phi}^{2} such that r 0 2 ​ ( m ϕ 2 − m c ​ r 2 ) ≃ 1.22 similar-to-or-equals superscript subscript 𝑟 0 2 superscript subscript 𝑚 italic-ϕ 2 superscript subscript 𝑚 𝑐 𝑟 2 1.22 r_{0}^{2}(m_{\phi}^{2}-m_{cr}^{2})\simeq 1.22 . We also quote the values of r 0 subscript 𝑟 0 r_{0} in lattice units, taken from Eq. (2.6) of Necco and Sommer ( 2002 ) , that we use in the present work.

$\beta$	${\lambda }_0$	${\eta }_{cr}$	$Z_P^{-1}$	$r_0/b$
5.75	0.5807	-1.271(10)	20.6(7)	3.291(16)
5.85	0.5917	-1.207(8)	20.1(7)	4.060(22)
5.95	0.6022	-1.145(6)	20.7(9)	4.912(29)

Table 4. Table 4: Simulation parameters in the NG phase. We used 60 gauge and scalar paired configurations at each simulation point. The spatial extension of the lattice is fixed to L ∼ 2.4 similar-to 𝐿 2.4 L\sim 2.4  fm, while the time extension is T = 6.0 𝑇 6.0 T=6.0  fm for the coarsest lattice spacing and T = 4.9 𝑇 4.9 T=4.9  fm for the two other ones.

$b$/fm	$\eta$	$b\mu$
0.152	-1.2714	0.0050 0.0087, 0.0131, 0.0183, 0.0277
($\beta =5.75$)	-1.2656	0.0131
	-1.2539	0.0183
	-1.2404	0.0131
	-1.231	0.0087, 0.0183
0.123	-1.2105	0.0040, 0.0070, 0.0100, 0.0120
($\beta =5.85$)	-1.2068	0.0040, 0.0070, 0.0100 ,0.0224
	-1.2028	0.0040, 0.0070, 0.0100, 0.0120, 0.0224
	-1.1949	0.0100
	-1.1776	0.0070, 0.0100, 0.0140, 0.0224
	-1.1677	0.0316
0.102	-1.1474	0.0066, 0.0077, 0.0116, 0.0145, 0.0185
($\beta =5.95$)	-1.1449	0.0060, 0.0077, 0.0116, 0.0145
	-1.1215	0.0077
	-1.1134	0.0077, 0.0108

Table 5. Table 5: Simulation parameters at ρ = 2.94 𝜌 2.94 \rho=2.94 , β = 5.85 𝛽 5.85 \beta=5.85 . We use 480 gauge and scalar paired configurations in the Wigner phase and 60 in the NG phase at each simulation point.

phase	$\left(L^3\times T\right)/b^4$	$\eta$	$b\mu$
Wigner	${16}^3\times 32$	-1.765	0.0140, 0.0210, 0.0280
		-1.7160	0.0210
		-1.6670	0.0210
NG	${20}^3\times 40$	-1.8376	0.0100
		-1.8197	0.0050 0.0100 0.0150
		-1.7808	0.0100

Table 6. Table 6: The values of η c ​ r A superscript subscript 𝜂 𝑐 𝑟 𝐴 \eta_{cr}^{A} and η c ​ r B superscript subscript 𝜂 𝑐 𝑟 𝐵 \eta_{cr}^{B} at fixed bare couplings g 0 2 = 6 / 5.95 superscript subscript 𝑔 0 2 6 5.95 g_{0}^{2}=6/5.95 , ρ = 1.96 𝜌 1.96 \rho=1.96 , λ 0 = 0.6022 subscript 𝜆 0 0.6022 \lambda_{0}=0.6022 for the values of μ ϕ 2 ≡ ( m ϕ 2 − m c ​ r 2 ) superscript subscript 𝜇 italic-ϕ 2 superscript subscript 𝑚 italic-ϕ 2 superscript subscript 𝑚 𝑐 𝑟 2 \mu_{\phi}^{2}\equiv(m_{\phi}^{2}-m_{cr}^{2}) given in the 3 rd (4 th ) column in lattice ( r 0 subscript 𝑟 0 r_{0} ) units.

${\eta }_{cr}^A$	${\eta }_{cr}^B$	$b^2{\mu }_{\varphi }^2$	$r_0^2{\mu }_{\varphi }^2$
-1.145(6)	-1.136(6)	0.0504(10)	1.22(3)
-1.138(2)	-1.132(2)	0.0171(10)	0.41(3)

Table 7. Table 7: The values of r A ​ W ​ I A ​ ( B ) superscript subscript 𝑟 𝐴 𝑊 𝐼 𝐴 𝐵 r_{AWI}^{A(B)} , b 2 ​ G P ​ S superscript 𝑏 2 subscript 𝐺 𝑃 𝑆 b^{2}G_{PS} and b 2 ​ M P ​ S 2 superscript 𝑏 2 superscript subscript 𝑀 𝑃 𝑆 2 b^{2}M_{PS}^{2} measured in the Wigner phase at ρ = 1.96 𝜌 1.96 \rho=1.96 and ρ = 2.94 𝜌 2.94 \rho=2.94 . The values of r A ​ W ​ I subscript 𝑟 𝐴 𝑊 𝐼 r_{AWI} in the time plateau region A ​ ( B ) 𝐴 𝐵 A(B) are denoted by r A ​ W ​ I A ​ ( B ) superscript subscript 𝑟 𝐴 𝑊 𝐼 𝐴 𝐵 r_{AWI}^{A(B)} (see Fig.  5 and Sect.  II.1 ). We recall Z P − 1 = [ G P ​ S W ​ r 0 2 ] μ → 0 η → η c ​ r superscript subscript 𝑍 𝑃 1 superscript subscript delimited-[] superscript subscript 𝐺 𝑃 𝑆 𝑊 superscript subscript 𝑟 0 2 → 𝜇 0 → 𝜂 subscript 𝜂 𝑐 𝑟 Z_{P}^{-1}=[G_{PS}^{W}r_{0}^{2}]_{\mu\to 0}^{\eta\to\eta_{cr}} .

$b\mu$	$\eta$	$r_{AWI}^A$	$r_{AWI}^B$	$b^2{G}_{PS}^W$	$b^2{M}_{PS}^2$
$\beta =5.75$,  ${\lambda }_0=0.5807$,  $b^2{m}_{\varphi }^2=-0.4150$,  $\rho =1.96$,  $L/b=16$, $T/b=32$
0.018	-1.1505	-0.1108(60)	-0.0748(44)	2.085(58)	0.14672(86)
0.028	-1.1505	-0.1063(64)	-0.0631(45)	2.188(56)	0.21952(87)
0.048	-1.1505	-0.0998(79)	-0.0410(60)	2.376(53)	0.35619(90)
0.018	-1.1898	-0.0718(59)	-0.0355(44)	2.095(59)	0.14710(81)
0.028	-1.1898	-0.0675(64)	-0.0233(45)	2.196(56)	0.21999(82)
0.048	-1.1898	-0.0609(79)	-0.0015(61)	2.382(53)	0.35671(87)
0.018	-1.3668	0.1017(61)	0.1389(44)	2.065(58)	0.14808(85)
0.028	-1.3668	0.1063(66)	0.1537(49)	2.176(56)	0.22072(87)
0.048	-1.3668	0.1137(81)	0.1757(70)	2.373(53)	0.35730(96)
$\beta =5.85$, ${\lambda }_0=0.5917$, $b^2{m}_{\varphi }^2=-0.4615$ $\rho =1.96$,  $L/b=16$, $T/b=40$
0.0224	-1.0983	-0.1090(65)	-0.0613(48)	1.482(34)	0.1652(11)
0.0316	-1.0983	-0.1078(75)	-0.0509(57)	1.597(32)	0.2280(11)
0.0387	-1.0983	-0.1063(82)	-0.0432(63)	1.684(31)	0.2755(11)
0.012	-1.1375	-0.0693(56)	-0.0372(40)	1.364(38)	0.0924(09)
0.0224	-1.1375	-0.0703(66)	-0.0238(48)	1.491(34)	0.1658(11)
0.0316	-1.1375	-0.0690(75)	-0.0132(57)	1.605(32)	0.2287(11)
0.0387	-1.1375	-0.0674(82)	-0.0055(63)	1.692(32)	0.2762(11)
0.0224	-1.2944	0.0849(74)	0.1266(52)	1.483(35)	0.1664(12)
0.0387	-1.2944	0.0879(86)	0.1456(65)	1.690(34)	0.2766(13)
$\beta =5.95$,  ${\lambda }_0=0.6022$,  $b^2{m}_{\varphi }^2=-0.4956$,  $\rho =1.96$,  $L/b=20$, $T/b=48$
0.0186	-0.9761	-0.1681(41)	-0.1377(37)	1.027(37)	0.11941(74)
0.0321	-0.9761	-0.1671(45)	-0.1255(35)	1.181(36)	0.20216(82)
0.0186	-1.0354	-0.1094(37)	-0.0798(36)	1.044(36)	0.12033(68)
0.0321	-1.0354	-0.1087(42)	-0.0676(34)	1.194(35)	0.20340(77)
0.01	-1.0771	-0.0675(35)	-0.0478(36)	0.955(37)	0.06673(65)
0.0186	-1.0771	-0.0680(36)	-0.0390(36)	1.051(36)	0.12072(66)
0.0321	-1.0771	-0.0677(41)	-0.0268(35)	1.200(35)	0.20396(75)
$\beta =5.95$,  ${\lambda }_0=0.6022$,  $b^2{m}_{\varphi }^2=-0.5289$,  $\rho =1.96$,  $L/b=20$, $T/b=48$
0.0186	-0.9761	-0.1572(18)	-0.1392(16)	1.005(033)	0.11813(91)
0.0321	-0.9761	-0.1561(24)	-0.1302(17)	1.170(033)	0.20087(99)
0.0186	-1.0354	-0.0996(18)	-0.0811(15)	1.034(034)	0.11957(76)
0.0321	-1.0354	-0.0983(24)	-0.0722(17)	1.191(034)	0.20276(85)
0.01	-1.0771	-0.0592(14)	-0.0465(16)	0.946(035)	0.06619(65)
0.0186	-1.0771	-0.0596(19)	-0.0401(15)	1.046(035)	0.12019(71)
0.0321	-1.0771	-0.0581(24)	-0.0316(19)	1.199(034)	0.20360(82)
$\beta =5.85$, ${\lambda }_0=0.5917$, $b^2{m}_{\varphi }^2=-0.4615$, $\rho =2.94$,  $L/b=16$, $T/b=40$
0.021	-1.667	-0.1527(100)	-0.0861(73)	1.587(36)	0.16168(100)
0.021	-1.716	-0.1073(100)	-0.0417(74)	1.604(37)	0.16264(97)
0.014	-1.765	-0.0634(100)	-0.0160(70)	1.518(40)	0.11239(88)
0.021	-1.765	-0.0622(110)	0.0026(77)	1.615(37)	0.16601(80)
0.028	-1.765	-0.0598(100)	0.0168(89)	1.709(36)	0.21274(100)

Table 8. Table 8: The values of b 2 ​ M P ​ S 2 superscript 𝑏 2 superscript subscript 𝑀 𝑃 𝑆 2 b^{2}M_{PS}^{2} , 2 ​ b ​ m A ​ W ​ I 2 𝑏 subscript 𝑚 𝐴 𝑊 𝐼 2bm_{AWI} and Z V subscript 𝑍 𝑉 Z_{V} measured in the NG phase at ρ = 1.96 𝜌 1.96 \rho=1.96 and ρ = 2.94 𝜌 2.94 \rho=2.94 .

$b\mu$	$\eta$	$b^2{M}_{PS}^2$	$2b{m}_{AWI}$	$Z_V$
$\beta =5.75$,  ${\lambda }_0=0.5807$,  $b^2{m}_{\varphi }^2=-0.5941$,  $\rho =1.96$,   $L/b=16$, $T/b=40$
0.005	-1.2715	0.0882(22)	0.01738(41)	0.9904(240)
0.0087	-1.2715	0.1214(15)	0.02206(34)	0.9752(93)
0.0131	-1.2715	0.1561(13)	0.02611(32)	0.9764(56)
0.0183	-1.2715	0.1945(12)	0.02989(32)	0.9791(43)
0.0277	-1.2715	0.2600(12)	0.03541(33)	0.9834(33)
0.0131	-1.2657	0.1489(12)	0.02359(31)	0.9735(52)
0.0183	-1.2540	0.1757(10)	0.02212(29)	0.9730(35)
0.0131	-1.2404	0.1244(09)	0.01357(26)	0.9630(38)
0.0131	-1.2277	0.1161(08)	0.00894(24)	0.9592(33)
0.0183	-1.2277	0.1557(08)	0.00894(24)	0.9658(26)
$\beta =5.85$,  ${\lambda }_0=0.5917$,  $b^2{m}_{\varphi }^2=-0.5805$  $\rho =1.96$,   $L/b=20$, $T/b=40$
0.004	-1.2106	0.0851(20)	0.01894(44)	0.9420(190)
0.007	-1.2106	0.1032(17)	0.02099(39)	0.9518(100)
0.01	-1.2106	0.1220(16)	0.02271(39)	0.9577(70)
0.012	-1.2106	0.1348(15)	0.02371(39)	0.9609(57)
0.004	-1.2069	0.0785(19)	0.01717(40)	0.9422(170)
0.007	-1.2069	0.0974(17)	0.01930(37)	0.9517(95)
0.01	-1.2069	0.1169(15)	0.02104(37)	0.9577(65)
0.0224	-1.2067	0.1981(13)	0.02609(39)	0.9701(30)
0.004	-1.2028	0.0717(18)	0.01529(36)	0.9420(150)
0.007	-1.2028	0.0914(16)	0.01749(35)	0.9517(87)
0.01	-1.2028	0.1115(14)	0.01924(36)	0.9577(61)
0.012	-1.2028	0.1250(14)	0.02024(36)	0.9607(50)
0.0224	-1.2028	0.1943(13)	0.02432(38)	0.9697(29)
0.01	-1.1950	0.1022(13)	0.01591(33)	0.9578(52)
0.007	-1.1777	0.0636(10)	0.00756(25)	0.9538(56)
0.01	-1.1777	0.0861(10)	0.00896(28)	0.9580(39)
0.014	-1.1777	0.1153(10)	0.01048(30)	0.9620(28)
0.0224	-1.1777	0.1746(10)	0.01309(32)	0.9674(20)
0.0316	-1.1677	0.2325(10)	0.01103(33)	0.9706(15)
$\beta =5.95$,  ${\lambda }_0=0.6022$,  $b^2{m}_{\varphi }^2=-0.5756$,  $\rho =1.96$,  $L/b=24$, $T/b=48$
0.006	-1.1475	0.0615(11)	0.01380(33)	0.9725(110)
0.0077	-1.1475	0.0705(11)	0.01444(32)	0.9729(84)
0.0116	-1.1475	0.0922(11)	0.01565(31)	0.9741(54)
0.0145	-1.1475	0.1089(11)	0.01640(30)	0.9749(43)
0.0185	-1.1475	0.1325(11)	0.01732(30)	0.9757(34)
0.006	-1.1450	0.0589(11)	0.01279(32)	0.9713(110)
0.0077	-1.1450	0.0681(11)	0.01344(31)	0.9720(80)
0.0116	-1.1450	0.0902(11)	0.01463(30)	0.9735(51)
0.0145	-1.1450	0.1071(11)	0.01538(30)	0.9744(41)
0.0077	-1.1216	0.0522(09)	0.00462(23)	0.9671(47)
0.0108	-1.1134	0.0691(09)	0.00224(21)	0.9669(34)
0.0077	-1.1134	0.0496(08)	0.00166(21)	0.9642(44)
$\beta =5.85$, ${\lambda }_0=0.5917$, $b^2{m}_{\varphi }^2=-0.5805$, $\rho =2.94$,   $L/b=20$, $T/b=40$
0.01	-1.8376	0.2044(30)	0.04219(65)	0.922(17)
0.005	-1.8197	0.1512(33)	0.03254(69)	0.847(35)
0.01	-1.8197	0.1785(28)	0.03560(59)	0.916(16)
0.015	-1.8197	0.2085(25)	0.03869(55)	0.942(10)
0.01	-1.7808	0.1292(22)	0.02232(48)	0.915(12)

Table 9. Table 9: The values of b 2 ​ M P ​ S 2 superscript 𝑏 2 superscript subscript 𝑀 𝑃 𝑆 2 b^{2}M_{PS}^{2} , 2 ​ b ​ m A ​ W ​ I 2 𝑏 subscript 𝑚 𝐴 𝑊 𝐼 2bm_{AWI} and Z V ~ subscript 𝑍 ~ 𝑉 Z_{\tilde{V}} measured at ρ = 0 𝜌 0 \rho=0 and η = 0 𝜂 0 \eta=0 . Boundary conditions were the same as in the study of the model with non-zero ρ 𝜌 \rho and η 𝜂 \eta .

$b\mu$	$\eta$	$b^2{M}_{PS}^2$	$2b{m}_{AWI}$	$1-Z_{\tilde{V}}$
$\beta =5.95$,  ${\lambda }_0=0.6022$,  $b^2{m}_{\varphi }^2=-0.5756$,  $\rho =0$,  $L/b=24$, $T/b=48$
0.0040	0	0.02729(39)	-18(6)$\times {10}^{-12}$	1(13)$\times {10}^{-12}$
0.0050	0	0.03352(41)	-2(8)$\times {10}^{-12}$	15(14)$\times {10}^{-12}$
0.0060	0	0.03967(42)	-3(6)$\times {10}^{-12}$	10(12)$\times {10}^{-12}$
0.0077	0	0.05005(45)	2(9)$\times {10}^{-12}$	-1(11)$\times {10}^{-12}$
0.0090	0	0.05794(47)	-6(7)$\times {10}^{-10}$	-7(7)$\times {10}^{-9}$
0.0100	0	0.06399(48)	-4(7)$\times {10}^{-10}$	18(13)$\times {10}^{-12}$
0.0116	0	0.07363(49)	-6(8)$\times {10}^{-12}$	14(10)$\times {10}^{-12}$

Equations

$$\displaystyle{{S}_{{L}}\!=\!b^{4}\!\sum_{x}\!\Big{\{}{\cal L}_{k}^{g}[U_{.}]\!+\!{\cal L}_{k}^{s}(\Phi)\!+\!{\cal V}(\Phi)}{+\bar{Q}{D_{L}[U_{.},\Phi]}Q\Big{\}}}\,,\;$$

$$\displaystyle{\cal L}_{k}^{g}[U_{.}]\,\,{\mbox{: SU($3$) plaquette action}}\,,$$

$${\mathcal{L}}_k^s(\Phi )+\mathcal{V}(\Phi )=$$

$$=\frac{1}{2}tr[{\Phi }^{†}(-{\partial }_{\mu }^{\ast }{\partial }_{\mu })\Phi ]+\frac{m_{\varphi }^2}{2}tr[{\Phi }^{†}\Phi ]+\frac{{\lambda }_0}{4}(tr[{\Phi }^{†}\Phi ]{)}^2,$$

$$(D_L[U_{.},\Phi ]Q)(x)={\gamma }_{\mu }{\tilde{\nabla }}_{\mu }Q(x)+\eta F(x)Q(x)+$$

$$-b^2\rho \frac{1}{2}F(x){\tilde{\nabla }}_{\mu }{\tilde{\nabla }}_{\mu }Q(x)+$$

$$\displaystyle-b^{2}\rho\frac{1}{4}\Big{[}(\partial_{\mu}F)(x)U_{\mu}(x)\widetilde{\nabla}_{\mu}Q(x+\hat{\mu})+$$

$$\displaystyle+(\partial_{\mu}^{*}F)(x)U_{\mu}^{\dagger}(x-\hat{\mu})\widetilde{\nabla}_{\mu}Q(x-\hat{\mu})\Big{]}\,.$$

$${\tilde{\nabla }}_{\mu }f(x)\equiv \frac{1}{2}({\nabla }_{\mu }+{\nabla }_{\mu }^{\ast })f(x),$$

$$b{\nabla }_{\mu }f(x)\equiv U_{\mu }(x)f(x+\hat{\mu })-f(x),$$

$$b{\nabla }_{\mu }^{\ast }f(x)\equiv f(x)-U_{\mu }^{†}(x-\hat{\mu })f(x-\hat{\mu }).$$

$$q_{\alpha ,a}^B(y)=\sum _{\xi }\overline{U}(2y,2y+\xi )[{\mathcal{A}}_{\xi }{]}_{\alpha ,a}\frac{(1-b{\xi }_{\mu }{\tilde{\nabla }}_{\mu })}{8}{\chi }^B(2y+\xi ),$$

$$S_L^F=\sum _{y,B}{\bar{q}}^B(y)\left\{\sum _{\mu }({\gamma }_{\mu }\otimes 11)D_{\mu }+\eta {\mathcal{F}}^B(y)\right\}{q}^B(y)+\text{mbox}O(b^2),$$

$$P^i(x)=\overline{Q}(x){\gamma }_5\frac{{\tau }^i}{2}Q(x),x_{\mu }=2y_{\mu }+{\xi }_{\mu }$$

$$\frac{1}{16}\sum _{\xi }{P}^i(2y+\xi )=\sum _{B=1}^4{s}_B{\overline{q}}^B(y)({\gamma }_5\otimes t_5)\frac{{\tau }^i}{2}{q}^B(y)+\text{mbox}O(b^2).$$

$${\tilde{A}}_{\mu }^{i,BW}(x)=\frac{1}{2}\overline{Q}(x-\hat{\mu }){\gamma }_{\mu }{\gamma }_5\frac{{\tau }^i}{2}{U}_{\mu }(x-\hat{\mu })Q(x)+$$

$$+\frac{1}{2}\overline{Q}(x){\gamma }_{\mu }{\gamma }_5\frac{{\tau }^i}{2}{U}_{\mu }^{†}(x-\hat{\mu })Q(x-\hat{\mu }),x=2y+\xi$$

$$\frac{1}{16}\sum _{\xi }{\tilde{A}}_{\mu }^{i,BW}(2y+\xi )=\sum _{B=1}^4{s}_B{\overline{q}}^B(y)({\gamma }_{\mu }{\gamma }_5\otimes t_5)\frac{{\tau }^i}{2}{q}^B(y)+\text{mbox}O(b^2).$$

$$Q(x)\to Q^{\prime }(x)=e^{-ix\cdot {\pi }_H}{M}_{H}Q(x)$$

$$\overline{Q}(x)\to {\overline{Q}}^{\prime }(x)=\overline{Q}(x)M_H^{†}{e}^{ix\cdot {\pi }_H},$$

$$S_L^{tm}=b^4\mu \sum _x\bar{Q}(x)i{\gamma }_5{\tau }^{3}Q(x).$$

$$C^{11}(x_0)=b^3\sum _{\mathbf{x}}⟨P^1(y_{\mathrm{source}})P^1(x)⟩{|}_{y_{\mathrm{source}}=0},$$

$$m_{AWI}^{\rm gen.\,off-diag}=\frac{\sum_{\bf x}\partial_{0}\langle\bar{Q}\gamma_{0}\gamma_{5}\frac{\tau^{1}}{2}Q^{\prime}(x)\bar{Q}^{\prime}\gamma_{5}\frac{\tau^{1}}{2}Q(0)\rangle}{2\sum_{\bf x}\langle\bar{Q}\gamma_{5}\frac{\tau^{1}}{2}Q^{\prime}(x)\bar{Q}^{\prime}\gamma_{5}\frac{\tau^{1}}{2}Q(0)\rangle}\Big{|}_{\eta_{cr}}$$

$$M_{\sigma }^2{r}_0^2=1.284\text{mbox}and{\lambda }_R\equiv \frac{M_{\sigma }^2}{2v_R^2}=0.441.$$

$$Z_{\varphi}^{1/2}=\Big{[}M_{\sigma}^{2}\Big{(}b^{4}\sum_{y}\langle\varphi_{0}(y)\varphi_{0}(0)\rangle\,-\,V_{4}v^{2}\Big{)}\Big{]}^{-1}\,,$$

$$m_{AWI}=\frac{\partial_{0}^{FW}\!\sum_{\bf x}\langle\tilde{A}_{0}^{1,BW}(x)P^{1}(0)\rangle}{2\sum_{\bf x}\langle P^{1}(x)P^{1}(0)\rangle}\Big{|}_{\eta_{cr}}\,.$$

$$\displaystyle Z_{\tilde{V}}\partial_{0}^{FW}\!\sum_{\bf x}\langle\tilde{V}_{0}^{2,BW}(x)P^{1}(0)\rangle\Big{|}_{\eta_{cr}}=$$

$$\displaystyle=2\mu\sum_{\bf x}\langle P^{1}(x)P^{1}(0)\rangle\Big{|}_{\eta_{cr}}+{\mbox{O}}(b^{2})\,.$$

$${\partial }_{\nu }{V}_{\nu }^2(x)=2\mu P^1(x),$$

$$r_0{M}_{PS}$$

$$2r_0{m}_{AWI}{Z}_{\tilde{V}}{Z}_P^{-1}$$

$$L_{\mathrm{toy}}^{Sym}=L_4+\sum _{n=1}^{\infty }{b}^{2n}{L}_{4+2n}.$$

$$Q(x)\to iQ(x^{PT}),$$

Full text

1. I Lattice formulation

2. I.1 Fermionic regularization

3. I.2 Renormalizability and other technicalities

4. I.3 Ignoring fermionic disconnected diagrams

5. I.4 Gauge and scalar sector renormalization

6. II Data analysis details

7. II.1 Analysis of section “Lattice study and results”

8. II.2 Check of analysis stability and systematics

9. II.3 Check and irrelevance of finite size effects

10. II.4 Studies at larger values of $\rho$ and $v$

11. II.5 Summary of numerical results

12. III $m_{\phi}^{2}$ independence of $\eta_{cr}$

13. IV Lattice artifacts ${\rm O}(b^{2n})$

14. V Further consistency checks

15. V.1 The critical model at $\rho=\eta=0$

16. V.2 $M_{PS}$ in the Wigner phase of the critical model

17. VI Towards a realistic model

I Lattice formulation

Lattice simulations of models with gauge, fermion and scalar fields are rather uncommon in the literature. See, however, Ref. De and Shigemitsu (1988) for a quenched study where, unlike here, the scalar field belongs to the fundamental representation of the SU(2) gauge group but no $\tilde{\chi}_{L}\times\tilde{\chi}_{R}$ symmetry violating terms are included.

For our purposes it is crucial to adopt a lattice formulation preserving the exact $\chi_{L}\times\chi_{R}$ symmetry, while including terms breaking the purely fermionic $\tilde{\chi}_{L}\times\tilde{\chi}_{R}$ invariance. We have thus chosen the lattice regularized action of the model (L1) to be of the form

[TABLE]

where the naïve lattice discretization of fermions was used. We set $\Phi=\varphi_{0}1\!\!1+i\varphi_{j}\tau^{j}$ and denote by the symbol $U_{.}$ the set of gauge links. Here and in the following by “(Ln)”, with n integer, we shall refer to the equation (n) of the Letter, also called “main text”. Setting $F(x)\equiv[\varphi_{0}1\!\!1+i\gamma_{5}\tau^{j}\varphi_{j}](x)$ the lattice Dirac operator reads

[TABLE]

Gauge covariant lattice derivatives are defined via

[TABLE]

The partial derivatives $\partial_{\mu}$, $\partial_{\mu}^{*}$ are the forward, backward derivatives $\nabla_{\mu}$, $\nabla_{\mu}^{*}$, respectively, with $U_{\mu}\!=\!1\!\!1$. We use periodic boundary conditions for the boson (gauge and scalar) fields, while for fermions we impose antiperiodic boundary conditions in time and periodic ones in three-dimensional space.

The lattice action (1) describes 2 (isospin) $\times$ 16 (doublers) fermionic flavours, even in the $b\rightarrow 0$ limit. In fact, the fermion Lagrangian term proportional to $\rho$, which is a lattice regularization of the ${\cal L}_{W}$ term of Eq. (L2), being the product of $b^{2}$ times a (dimension six) fermion bilinear operator with two derivatives, does not remove the doublers in the continuum limit. However, for the goals of the present quenched study, the presence of doublers in the valence sector, which does not alter the $\beta$-function and hence preserves the asymptotic freedom of the model, is harmless. With respect to the unquenched case, the quenched approximation is expected to entail only some systematic error (typically of order 10–20% based on Lattice QCD experience) in the value of the coefficient $c_{1}$ of the NP fermion mass term (see Eq. (L13)) that is possibly generated in the NG phase EL. For analogous unquenched studies a lattice Dirac operator $D_{L}[U_{.},\Phi]$ with no doublers would of course be needed, which could e.g. be built by using the domain-wall Kaplan (1992) or overlap Neuberger (1998) fermion kernel for two flavours with the addition of suitable Yukawa and ${\cal L}_{W}$-like terms.

I.1 Fermionic regularization

Following closely Refs. Kluberg-Stern et al. (1983), Patel and Sharpe (1993) and Luo (1997) for staggered fermions, one can analyze the fermion flavour content of the action (1). We first rewrite its fermionic sector in terms of the field $\chi(x)={\cal A}_{x}^{-1}Q(x)$ with ${\cal A}_{x}=\gamma_{1}^{x_{1}}\gamma_{2}^{x_{2}}\gamma_{3}^{x_{3}}\gamma_{4}^{x_{4}}$, then we pass to the so-called flavour basis

[TABLE]

where $\alpha$, $a$ and $B$, all taking values in $\{1,2,3,4\}$, denote Dirac, taste and replica indices. While the $\chi^{B}(x)$ fields are defined on the fine lattice $x_{\mu}=2y_{\mu}+\xi_{\mu}$, $\xi_{\mu}=0,1$ with $\overline{U}(2y,2y+\xi)$ denoting the average of link products along the shortest paths from $2y$ to $2y+\xi$, the flavour basis fields $q^{B}_{\alpha,a}(y)$ live on the coarse lattice with coordinate $y_{\mu}$ (here we have set $b=1$). In terms of these fields close to the continuum limit the lattice fermionic action reads

[TABLE]

where for each $B$ value the matrix ${\cal F}^{B}(y)=\varphi_{0}(2y)(1\!\!1\otimes 1\!\!1)+s_{B}i\tau^{i}\varphi_{i}(2y)(\gamma_{5}\otimes t_{5})$, $s_{1,2}=1=-s_{3,4}$, has a “Dirac $\otimes$ taste” structure with $t_{\mu}=\gamma_{\mu}^{*}$ taste matrices. Unlike in the original $Q(x)$ basis, in the flavour basis $q(y)$ defined via Eq. (5) space-time and flavour degrees of freedom are completely disentangled from each other. As the action (6) is diagonal in taste and replica indices up to $O(b^{2})$, it obviously describes 32 fermion species, namely 4 replicas of the 4 tastes of the isospin doublet $Q^{T}=(u,d)$.

The quark bilinears expressed in the $Q$ basis, once summed over one $2^{4}$ hypercube of the fine lattice, have well defined quantum numbers in the classical continuum limit, as one checks by rewriting them in the flavour $q(y)$ basis. For example the isotriplet pseudoscalar density

[TABLE]

summed over the hypercube ($\xi_{\mu}=0,1$) yields in flavour basis the corresponding isotriplet, taste diagonal density

[TABLE]

Similarly for the lattice backward isotriplet axial current

[TABLE]

summing over the hypercube one finds in the flavour basis the corresponding isotriplet, taste diagonal current,

[TABLE]

I.2 Renormalizability and other technicalities

Beyond tree level, loop effects do not generate $d\leq 4$ terms besides the operators that are already present in the lattice action (1). This follows from the symmetries (gauge, $H(4)$, $\chi_{L}\times\chi_{R}$, C, P, T) of the model (1) which were discussed in Ref. Frezzotti and Rossi (2015) and are preserved by the lattice regularization. The latter is chosen so as to respect the “spectrum doubling symmetry” Montvay and Muenster (1997) of naïve fermions, implying that the action (1) is invariant under the transformations

[TABLE]

where $H$ is an ordered set of four-vector indices $H\equiv\{\mu_{1},...,\mu_{h}\},\,(\mu_{1}<\mu_{2}<...<\mu_{h})$. For $1\leq h\leq 4$ there are 16 four-vectors $\pi_{H}$, with $\pi_{H,\mu}=\pi$ if $\mu\in H$ or $\pi_{H,\mu}=0$ otherwise and 16 matrices $M_{H}\equiv(i\gamma_{5}\gamma_{\mu_{1}})...(i\gamma_{5}\gamma_{\mu_{h}})$. The use of symmetric (gauge covariant) lattice derivatives $\widetilde{\nabla}_{\mu}$ acting on $Q$ in the ${\cal L}_{W}$ term is dictated by the need of preserving the spectrum doubling symmetry. The latter, being an exact symmetry of $S_{L}$, is an invariance of the lattice effective action $\Gamma_{L}[U_{.},\Phi,Q]$. Hence all fermion species, including doublers, enter in a symmetric way in the effective action, in particular in the non-irrelevant local terms $\sum_{y}\sum_{B}\overline{q}^{B}(y)(\gamma_{\mu}\otimes 1\!\!1)D_{\mu}q^{B}(y)$ and $(\eta-\bar{\eta})\sum_{y}\sum_{B}\overline{q}^{B}(y){\cal F}^{B}(y)q^{B}(y)$ (in the flavour basis).

These nice renormalization properties imply in particular that $\eta_{cr}$, as determined from the vanishing of $r_{AWI}$ (Eq. (L15)), is unique for all fermion species.

As the model (1) is studied in the quenched approximation, problems with exceptional configurations of gauge and scalar fields supporting very small eigenvalues of $D_{L}$ plague, as expected, the Monte Carlo sampling of fermionic observables. In order to avoid such problems and have a robust IR cutoff in the fermionic propagators we include in the lattice action (1) an additional “twisted mass” term Frezzotti et al. (2001); Frezzotti and Rossi (2004)

[TABLE]

This term does not spoil power counting renormalizability, it only breaks in a soft way the $\chi_{L}\times\chi_{R}$ symmetry and can thus be safely removed at the end of the calculations by taking (even at fixed lattice spacing) the limit $\mu\to 0$.

In quenched lattice QCD there is overwhelming numerical evidence (see e.g. Refs. Butler et al. (1994); Bernard et al. (1998); Aoki et al. (2000); Blum et al. (2004)) that at large enough Euclidean time separations the correlators (vacuum expectation values of local fields) with fixed spatial momentum admit a spectral decomposition in terms of intermediate states with definite quantum numbers. The products of matrix elements entering this spectral representation (in particular those involving excited states) may at worst take values that are not consistent with time-reflection positivity, because the quenched model, owing to the omission of virtual quark contributions, is not unitary. We shall here assume that the same holds true in our quenched model (1). The validity of this assumption about the lowest-lying states of the spectrum is confirmed by strong numerical evidence for all the correlators we actually evaluated – as we routinely checked in our data analyses. As a typical example, we show i) in Fig. 1 the PS-meson two–points correlator in the NG phase of the critical model for $\beta=5.95$ and two different twisted quark mass values (for more details on simulation parameters see Table 4), together with its best fit to the exponential time dependence that is expected from one-particle state dominance, and ii) in Fig. 2 the corresponding effective mass plot, from the plateau of which the PS-meson mass is eventually estimated.

I.3 Ignoring fermionic disconnected diagrams

In our quenched study of the lattice model (1) the calculation of the correlators appearing in Eq. (L11) and Eqs. (L16)-(L17) was carried out by consistently neglecting all the fermionic disconnected Wick contractions. The latter arise in correlators with isospin non-singlet operators because of the occurrence of isospin changing fermion-antifermion-$\Phi$ vertices and are statistically quite noisy. Neglecting fermionic disconnected diagrams is, however, justified in our case as it affects neither the results for the mixing coefficient $\bar{\eta}$ (which determines $\eta_{cr}$, see Eq. (L5)) nor the pseudoscalar isotriplet meson mass, $M_{PS}$. It entails only an immaterial (and numerically small) modification of the renormalized matrix element ratio $\frac{Z_{\tilde{A}}}{Z_{P}}m_{AWI}$, see Eq. (L11).

The arguments behind the statements above are given in this section. Following Veneziano (1976), we call “one-boundary” diagrams the fermionic connected diagrams in which valence fermions localized at $y_{\,\rm source}=0$ are contracted with the valence quarks sitting at $x\neq y_{\,\rm source}$, and “two-boundaries” diagrams those in which valence quarks and antiquarks are contracted with each other separately at $y_{\,\rm source}=0$ and $x\neq y_{\,\rm source}$. The valence fermion structure of the two-boundaries Wick contractions that were neglected in evaluating all the correlators considered in this paper is shown in the diagrams of Fig. 3 .

Concerning the evaluation of $\bar{\eta}-\eta$ the key observation is that the bare SDE of $\tilde{\chi}_{L}$ and/or $\tilde{\chi}_{R}$ transformations (see Eqs. (3.15)–(3.16) of Ref. Frezzotti and Rossi (2015)), as well as their renormalized counterparts, see e.g. Eq. (L5) for the $\tilde{\chi}_{L}$ SDE, can be split in two independent identities according to the number of diagram boundaries ensuing from the way valence quarks are contracted. The reason is that the one-boundary and the two-boundaries diagrams have different functional dependencies on $x_{0}$, which must be separately matched between the l.h.s. and the r.h.s. of the SDE. In fact, for large enough time separations, the one-boundary Wick contractions are dominated (as found in many numerical studies) by the lowest–lying isospin non-singlet intermediate state, while the two-boundaries diagrams are expected to be saturated by different states involving an odd number of $\Phi$-particles and possibly gluons. The quantity $\bar{\eta}-\eta$, being an operator mixing coefficient, is a factor common to both one-boundary and two-boundaries diagrams. In other words, $\bar{\eta}$, and thus $\eta_{cr}$, can be in principle computed equally well by restricting attention to either the one-boundary or the two-boundaries Wick contractions relevant for the SDE of $\tilde{\chi}_{L,R}$ transformations. We also recall that such SDE make perfectly sense in the quenched approximation, because it can be defined as a renormalizable Lagrangian theory with ghost fields à la Morel–Bernard–Golterman (see Refs. Morel (1987); Bernard and Golterman (1992)).

The mass $M_{PS}$ of the isotriplet pseudoscalar meson, as obtained, say, from the correlator (here $x=({\bf x},x_{0})$)

[TABLE]

can also be evaluated by just restricting attention to the one-boundary Wick contractions. To see that two-boundaries diagrams do not affect the value of $M_{PS}$, let us examine them in terms of their possible cuts, see Fig. 3. One can either cut through one of the two fermion loops stemming from the Wick contractions at $y_{\,\rm source}=0$ and at $x\neq y_{\,\rm source}$, or one can cut in the middle of the two-boundaries diagram (no fermionic lines are cut in this case). In the first case one exposes terms decreasing exponentially in Euclidean time with $M_{PS}$, but one also recognizes from the cut diagram topology that these terms merely contribute to the renormalization of either $P^{1}(x)$ (as in left panel of Fig. 3) or $P^{1}(y_{\,\rm source})$ without affecting the value of $M_{PS}$. In the other case (right panel of Fig. 3) one obtains terms that, owing to the quenched approximation, are expected to be dominated by a set of intermediate states that have nothing to do with a pseudoscalar meson, rather they involve an odd number of scalar $\Phi$-particles plus gluons.

Unlike $\bar{\eta}$ and $M_{PS}$, the value of bare correlator ratio $m_{AWI}$ (Eq. (L11)) gets instead modified if two-boundaries diagrams are neglected. The reason is that two-boundaries diagrams contribute, as we just saw discussing the case of $M_{PS}$, to the renormalization of the isotriplet pseudoscalar quark density and the axial current four-divergence, as well as to mixing of these quark bilinears with other operators involving scalar fields but no fermions. In the quenched lattice ${\cal L}_{\rm toy}$ model of interest to us, however, the properly renormalized ratio, $m_{AWI}^{R}\equiv\frac{Z_{\tilde{A}}}{Z_{P}}m_{AWI}$, is obtained (at $\eta=\eta_{cr}$) if one consistently ignores all two-boundaries diagrams in the evaluation of the renormalization constants $Z_{\tilde{A}}$ and $Z_{P}$, which is what we have done in our numerical study.

A proof of the above statement can be given as follows. Let us consider in the quenched approximation the model ${\cal L}_{\rm toy}$ (see Eq. (1)) and its extension to two generations of fermion doublets, i.e. an analogous model with fermion sector action $S_{L}^{F}=\overline{Q}D_{L}Q+\overline{Q}^{\prime}D_{L}Q^{\prime}$, with $Q=(u,d)$ and $Q^{\prime}=(u^{\prime},d^{\prime})$ entering with identical $\tilde{\chi}$ violating parameters, $\rho$ and $\eta$. Owing to fermion quenching, the latter model will admit an EL description analogous to the one of ${\cal L}_{\rm toy}$ (see Eq. (L13)), with the fermion sector simply replicated for the $Q$- and the $Q^{\prime}$- generations. In particular the $\tilde{\chi}$ violating NP terms in the EL will occur multiplied by coefficients $c_{1}$, $c_{2}$, …that are generation-independent and have the same value as for ${\cal L}_{\rm toy}$. In the two generations lattice model introduced above the non-vanishing of $c_{1}$ is equivalent to the non-vanishing of the renormalized correlator ratio $(\frac{Z_{\tilde{A}}}{Z_{P}}m_{AWI})|^{\rm gen.\,off-diag}$, where

[TABLE]

and $\bar{Q}\gamma_{\mu}\gamma_{5}\frac{\tau^{1}}{2}Q^{\prime}$ is the current of a $\tilde{\chi}$-axial generation off-diagonal transformation. But this renormalized correlator ratio of the two generations model, to which no two-boundaries diagrams contribute owing to the generation off-diagonal nature of the axial current, is just equal to the ratio $\frac{Z_{\tilde{A}}}{Z_{P}}m_{AWI}$ that is computed in the one generation model ${\cal L}_{\rm toy}$ by neglecting two-boundaries diagrams both in $m_{AWI}$ and in $\frac{Z_{\tilde{A}}}{Z_{P}}$. The latter ratio (at $\eta=\eta_{cr}$) is hence a renormalized and $\tilde{\chi}$ covariant quantity that is non-vanishing if and only if the coefficient $c_{1}$ in the EL expression (L13) is non-zero.

I.4 Gauge and scalar sector renormalization

Owing to the quenched approximation, gauge and scalar configurations can be independently generated. Similarly, action parameter renormalization can be carried out independently for gauge and scalar fields. The gauge coupling is renormalized by keeping the Sommer length scale $r_{0}$ Guagnelli et al. (1998); Necco and Sommer (2002) fixed in physical units. In the following for mere orientation we take its value from calculations of the static quark potential in QCD and set $r_{0}=0.5$ fm $\sim 1/(394~{}{\rm MeV})$.

The choice of the phase in which the theory lives is made by taking $m_{\phi}^{2}$ either larger (Wigner) or smaller (NG phase) than the critical value $m_{cr}^{2}$ at which the scalar susceptibility diverges (in practice it has a sharp peak on a large enough space-time volume). The values of the critical mass are listed in Tab. 1.

In our study of $m_{AWI}$ and $M_{PS}$ in the NG phase the bare scalar parameters in (1) at each value of the lattice spacing are determined by the renormalization conditions

[TABLE]

Moreover the renormalized scalar vev, $v_{R}=Z_{\varphi}^{1/2}v$, where $Z_{\varphi}$ denotes the scalar field ($\varphi$) renormalization constant, is fixed in physical units by the condition $v_{R}^{2}r_{0}^{2}\!=\!1.458$. The $\varphi$ field renormalization constant is computed from

[TABLE]

where $V_{4}\!=\!L^{3}\!\times\!T$ is the lattice four-volume and the mass $M_{\sigma}$ of the one-$\varphi_{0}$ state is extracted from the time decay of the correlator $b^{3}\sum_{\bf x}\langle\varphi_{0}(x_{0},{\bf x})\varphi_{0}(0,{\bf 0})\rangle$. Once $M_{\sigma}$, $Z_{\varphi}$ and the bare scalar vev $v$ are known, $\lambda_{R}$ is easily evaluated.

In a numerical simulation on a finite lattice the scalar vev is always zero, even if $m_{\phi}^{2}<m_{cr}^{2}$. Hence following Refs. Bulava et al. (2013); Goeckeler and Leutwyler (1991) the technique of “fixing” the axial part of the global $\chi_{L}\times\chi_{R}$ symmetry was exploited in order to get $\langle\Phi\rangle=v1\!\!1\ \neq 0$ in the NG phase. As this technique yields the correct NG phase results up to finite size effects, by simulating on lattices of different spatial size $L$ we checked that for the observables discussed in the section “Lattice study and results” the finite-$L$ systematic errors on the quoted values at $L\simeq 2.4$ fm are well below the uncertainties shown in Figs. 2 and 3 of the main text, or those quoted in Sect. II.2. Indeed, as detailed in Sect. II.3, simulations at three $L$-values in the range $2.0~{}{\rm fm}<L<3.4~{}{\rm fm}$ show that the deviations of $M_{PS}$ and $m_{AWI}$ evaluated at $L\simeq 2.4$ fm from their infinite volume counterparts can be conservatively estimated to be around 1-2% and 5-7%, respectively, and are hence marginally significant as compared to the corresponding purely statistical uncertainties of about 1.5% and 3%.

The bare parameters used in the NG-phase simulations are listed in Tab. 1. Concerning the simulations in the Wigner phase aimed at determining $\eta_{cr}$, we have worked at the same lattice resolution (i.e. $\beta$) and the same bare quartic scalar coupling ($\lambda_{0}$) as employed in the NG phase. For each lattice resolution we have chosen the value of the bare scalar mass parameter $b^{2}m_{\phi}^{2}$ so as to obtain $(m_{\phi}^{2}-m_{cr}^{2})r_{0}^{2}=1.22$ within errors. The latter condition entails a small and uniform O($b^{2}(m_{\phi}^{2}-m_{cr}^{2})$) contribution to the lattice artifacts in the resulting $\eta_{cr}$ values at the three lattice spacings. These UV cutoff artifacts are then propagated in the analysis of NG observables (see Sect. II) and safely removed by extrapolating to the continuum limit.

A comment about the triviality of the scalar quartic coupling is in order here. Triviality implies that in order to really take the limit $\Lambda_{UV}\sim b^{-1}\to\infty$ while keeping fixed the aforementioned low energy renormalization conditions, among which the one on $\lambda_{R}$ (see Eq. (15)), the bare quartic coupling $\lambda_{0}$ should be made to diverge logarithmically with increasing $\Lambda_{UV}$. In practice however we shall estimate the limit $b^{-1}\to\infty$ of renormalized observables by considering their values at a few values of $1/b^{2}$ which differ from each other by a factor of about two and correspond to values of $\lambda_{0}$ and $g_{0}^{2}$ that change only by a few percent in the range $\lambda_{0}\lesssim 0.6$ and $g_{0}^{2}\gtrsim 1.0$ . Under these circumstances it is a very plausible assumption that scaling violations behave effectively as polynomials in $b^{2}$ (about the absence of odd-powers of $b$ see Sect. IV), because terms O($b^{2}\log^{k}(b^{2})$) – with any integer $k$ – can not be resolved from genuine O($b^{2}$) terms within numerical errors and in the limited range of $1/b^{2}$-values that we explore. This assumption that we make in our continuum limit analyses i) makes practically harmless the problem of triviality and ii) is fully analogous to the assumption that are made (and well checked) when taking the continuum limit in all the numerical studies of lattice QCD, and more recently lattice QCD+QED.

II Data analysis details

II.1 Analysis of section “Lattice study and results”

In the Wigner phase we compute $r_{AWI}$ from the ratio (L15) for the values of the bare parameters $\eta$ and $\mu$ listed in Tab. 2 and for a value of $\tau\equiv y_{0}-x_{0}=0.6$ fm chosen conveniently (see below) and kept fixed for the three lattice spacings. We fit the dependence of $r_{AWI}$ on $\eta$ and $\mu$ through the linear ansatz $r_{AWI}=K_{0}+K_{1}\eta+K_{2}\mu$, which turns out to describe very well our data. We thus identify $\eta_{cr}$ with $-K_{0}/K_{1}$. The resulting values of $\eta_{cr}$ for the three lattice spacings corresponding to $\beta=5.75,5.85,5.95$ are given in Tab. 3. There we also give the values of the Sommer scale $r_{0}$ and of the inverse of the pseudoscalar density renormalization constant $Z_{P}$, which is defined in a hadronic scheme in the Wigner phase as $Z_{P}^{-1}=r_{0}^{2}\langle 0|P^{1}|{\rm PS~{}meson}\rangle^{W}$ (see main text, Sect. “Lattice study and results” ). The data for $\langle 0|P^{1}|{\rm PS~{}meson}\rangle^{W}$ at fixed lattice spacing are within statistical errors well linear in $\mu$ and in $\eta$ and were extrapolated to $\mu=0$ and $\eta=\eta_{cr}$ in order to obtain the quoted results for $Z_{P}^{-1}$.

A remark is in order about the physics of the model in the Wigner phase, i.e. for $m_{\phi}^{2}-m_{cr}^{2}>0$. Owing to the presence in the lattice action of the twisted mass term $\propto\mu$ (see Eq. (12)), which is needed to give a robust IR cutoff to all observables, here, even if $\Phi$ takes zero vacuum expectation value, strong interactions induce anyway, just like in QCD, the spontaneous breaking of the $\tilde{\chi}$ symmetries. At $\eta=\eta_{cr}$ the latter get restored up to negligible O($b^{2}$), the scalar field decouple from quark and gluons and (in large volume) one observes pseudoscalar mesons with squared mass $M_{PS}^{2}\sim\mu$ – see the discussion in Sect. V.2. In this respect the name of Wigner phase may be slightly abusive. The lattice spatial size $L$ is large enough that the pseudoscalar meson mass at all our “Wigner phase” simulation points satisfies $M_{PS}(\eta,\mu)L\geq 5.3$.

Another technical remark concerns our particular choice of the lattice estimator for $\eta_{cr}$. By definition this special value of $\eta$ is the one at which the correlator ratio $r_{AWI}$ (see Eq. (L15)) vanishes for time separations $x_{0}$ and $y_{0}-x_{0}$ both sufficiently away from zero, up to O($b^{2}$) lattice artifacts. Due to the latter, in practice one obtains different lattice estimators of $\eta_{cr}$ for different ranges of $x_{0}$ (while we keep $y_{0}-x_{0}\equiv\tau=0.6$ fm fixed), because different combinations of intermediate states contribute to three-point correlators entering $r_{AWI}$. In principle any choice of the $x_{0}$–range is legitimate, however in order to keep lattice artifacts under control we make a choice of $x_{0}$-range (the region A discussed in next paragraph) where only intermediate states with low mass contribute to the correlators. By considering a second reasonable choice of $x_{0}$–range (the region B below) we also explicitly check that the difference between the two resulting lattice estimators of $\eta_{cr}$ vanishes in the continuum limit, as expected. The selection of these two $x_{0}$ ranges is suggested by the fact that the profile of $r_{AWI}$ as a function of $x_{0}$ (Fig. 5) typically exhibits two plateau regions.

The $r_{AWI}$ values shown in Fig. 1 of the main text and the $\eta_{cr}$ results of Tab. 3 all refer to the time plateau region at $x_{0}$ smaller than $T/2$ that we call region A. In this case the $x_{0}$ range is chosen in such a way that the three-point correlators in Eqs. (L16)-(L17) are dominated by the PS-meson state at times between 0 and $x_{0}$, by the one-$\Phi$-particle state at times between $x_{0}$ and $y_{0}=x_{0}+\tau$ and by the vacuum at times larger than $y_{0}$. The time plateau region of $r_{AWI}$ called region B corresponds instead to $x_{0}$ larger than $T/2$. In this situation the same three-point correlators are dominated by the vacuum at times between 0 and $x_{0}$, by a PS-meson plus one-$\Phi$-particle state at times between $x_{0}$ and $y_{0}$ and by one-PS-meson state at times larger than $y_{0}$. We check in Fig. 5 that the difference in the values of $\eta_{cr}$ obtained by evaluating $r_{AWI}$ from the two plateau regions A and B vanishes as an O($b^{2}$) quantity in the continuum limit.

In the NG phase we compute $M_{PS}$ and $m_{AWI}$ for the values of $\eta$ and $\mu$ listed in Tab. 4. Again strong interactions are responsible for the spontaneous breaking of $\tilde{\chi}$ symmetries (restored at $\eta_{cr}$), with the vacuum being polarized by the combined effect of the twisted mass term ($\propto\mu$) and the NP $\tilde{\chi}$ breaking term(s) (e.g. $\propto c_{1}\Lambda_{S}$) in the effective Lagrangian (see Eq. (L13)). The lattice spatial size $L$, which is kept fixed around $2.4$ fm (if we assume $r_{0}=0.5$ fm for the Sommer scale $r_{0}$) at the three lattice spacings, is large enough that the pseudoscalar meson mass $M_{PS}$ always satisfies $M_{PS}(\eta,\mu)L\geq 4.7$.

The pseudoscalar (PS) meson mass $M_{PS}$ is straightforwardly evaluated from the two-point correlator $\sum_{\bf x}\langle P^{1}(x)P^{1}(0)\rangle$, which at all our simulation points exhibits a clear exponential dependence on the time separation $x_{0}$, with ground PS-state dominance for large $x_{0}$ (typically $\gtrsim 1$ fm). In terms of the lattice fermion operators (7)-(9) $m_{AWI}$ reads (see Eq. (L11))

[TABLE]

We extrapolate these observables to $\eta_{cr}$ and $\mu=0$ with the ansatz $M_{PS}^{2}=Y_{0}+Y_{1}\eta+Y_{2}\mu+Y_{3}\eta^{2}+Y_{4}\mu^{2}+Y_{5}\eta\mu\;$ and $m_{AWI}=Y_{0}+Y_{1}\eta+Y_{2}\mu+Y_{3}\mu^{2}\;$.

In this way our analysis assumes no special value for $c_{1}\Lambda_{S}$. In particular the chosen fit ansatz for the $\mu$-dependence of our observables and the variants described in Sect. II.2 are suitable to describe the data even in the case of $c_{1}\Lambda_{S}=0$, when $m_{AWI}$ should be (owing to lattice artifacts) an analytic function of $\mu$ that at $\eta=\eta_{cr}$ vanishes smoothly as $\mu\to 0$, and in the same limit the relation $M_{PS}^{2}\sim\mu$ is expected to hold – up to small anomalous “chiral log” corrections in the quenched approximation Bernard and Golterman (1992); Sharpe (1992). If instead $c_{1}\Lambda_{S}\neq 0$, then the low energy effective theory of the critical model must involve both massive (even in the $\mu\to 0$ limit) pseudoscalar meson and the Goldstone modes of the scalar field. In this situation at $\eta=\eta_{cr}$ our observables will depend on the effective quark mass $\propto\sqrt{(c_{1}\Lambda_{S})^{2}+\zeta^{2}\mu^{2}}$, which in turn admits a smooth Taylor expansion around the $\mu=0$ point. Here $\zeta$ denotes a relative normalization factor which is immaterial for our discussion.

The definition of a renormalized counterpart of the matrix element ratio that we have suggestively denoted by $m_{AWI}$ (see Eq. (L11)) requires computing both $Z_{P}^{-1}$, which, as discussed in the section “Lattice study and results”, is provided by an hadronic matrix element evaluated for $(\eta,\mu)\to(\eta_{cr},0)$ in the Wigner phase, and $Z_{\tilde{A}}=Z_{\tilde{J}}|_{\eta_{cr}}=Z_{\tilde{V}}$. Since this last renormalization factor is independent of the twisted mass parameter $\mu$ appearing in the action term (12), it can be accurately computed by imposing, at $\eta=\eta_{cr}$, the validity of the $\tilde{\chi}$ SDE involving the bare lattice vector current $\tilde{V}^{2,BW}=\tilde{J}^{2,BW}_{L}+\tilde{J}^{2,BW}_{R}$, namely

[TABLE]

Indeed from the exact WTI of $\chi_{L}\times\chi_{R}$ symmetry that is satisfied by the lattice $\chi$-vector current $V_{\nu}^{2}$,

[TABLE]

we learn that the combination $2\mu P^{1}$ needs neither to be renormalized nor rescaled by any factor. In practice the estimates of $Z_{\tilde{V}}$ that we obtained in the NG phase by imposing Eq. (18) at the simulation points in Tab. 4 exhibit only a smooth dependence on $\eta$ and $\mu$, which allowed us to extrapolate them to $\eta\to\eta_{cr}$ and $\mu\to 0$ by assuming the polynomial ansatz $a+b\eta+c\mu+d\mu^{2}$.

II.2 Check of analysis stability and systematics

In order to check the stability of the results presented in the section “Lattice study and results” of the main text and the underlying analysis of Sect. II.1 we have studied the impact of modifying in turn

the plateau choice for $r_{AWI}$, which was taken in the time region either A or B of Fig. 5; 2. 2.

the fit function for $r_{AWI}$ in the Wigner phase, which was chosen to be either of the form $K_{0}+K_{1}\eta+K_{2}\mu$ or extended by adding the non-linear term $K_{4}\mu^{2}$ (data are very well linear in $\eta$). In addition the mean value over the results from the two fits above was considered; 3. 3.

the choice of the time plateau for $M_{PS}$ and $m_{AWI}$, taken to be either $[1.4,2.2]$ fm or $[1.1,2.4]$ fm; 4. 4.

the fit functions for $m_{AWI}$ and $M_{PS}^{2}$ in the NG phase, which were both either truncated to $Y_{0}+Y_{1}\eta+Y_{2}\mu$ or extended by adding in turn terms non-linear in $\eta$ and/or $\mu$, such as $Y_{4}\mu^{2}$, $Y_{3}\eta^{2}+Y_{4}\mu^{2}$ or $Y_{3}\eta^{2}+Y_{4}\mu^{2}+Y_{5}\eta\mu$.

Our care in checking the determination of $\eta_{cr}$ at each lattice spacing is due to the fact that its uncertainty, once propagated onto $m_{AWI}$ and $M_{PS}^{2}$, turns out to be the main source of error in the final results (20)–(21). For each combination of the above choices the properly renormalized results of the fits to $(\eta,\mu)=(\eta_{cr},0)$, i.e. $r_{0}M_{PS}$ and $2r_{0}m_{AWI}Z_{\tilde{V}}Z_{P}^{-1}$, were extrapolated linearly in $b^{2}/r_{0}^{2}$ to the continuum limit. To estimate the impact of possible O($b^{4}$) and higher order lattice artifacts, we have also considered a continuum extrapolation excluding the coarsest lattice spacing, but imposing that the two choices of taking either $\eta^{A}_{cr}$ or $\eta^{B}_{cr}$ yield the same continuum value. For a typical choice of the time plateau and the fit function in $\eta$ and $\mu$ for $M_{PS}r_{0}$, in Fig. 6 we show the continuum limit extrapolations obtained by working at $\eta^{A}_{cr}$ or at $\eta^{B}_{cr}$ and using data at three lattice spacings, as well as the combined continuum extrapolation that results from using data on the two finer lattices at both $\eta^{A}_{cr}$ and $\eta^{B}_{cr}$. The analogous comparison for the renormalized ratio $2r_{0}m_{AWI}Z_{\tilde{V}}Z_{P}^{-1}$ is presented in Fig. 7 .

The combination of the variants discussed above produced 72 different analyses, none of which gave a result for $m_{AWI}$ or $M_{PS}$ that is compatible with zero (i.e. with the no-mechanism scenario). The median over all the analyses, excluding only a few ones for $M_{PS}^{2}$ that yield a $\chi^{2}$ value per degree of freedom larger than 2, gives

[TABLE]

where the first error is statistical and the second systematic. These two error terms were estimated by attributing equal weight to all analyses, following the method of Ref. Carrasco et al. (2014) (see there the text around Eq. (28)).

Further strong evidence for a non-zero fermion mass in the NG phase EL comes from the following “reductio ad absurdum” argument. Let us assume that no fermion mass term is generated. If this were the case, at $\eta=\eta_{cr}$ one should have $M_{PS}={\mbox{O}}(b^{2})$, and hence $M_{PS}^{2}$ would approach a vanishing continuum limit with only O($b^{4}$) artifacts.

In Fig. 8 we plot $r_{0}^{2}M_{PS}^{2}$ as a function of $(b/r_{0})^{4}$ (the original data are those at $\eta=\eta_{cr}^{A}$ already shown in Fig. 3 of the main text, and a subset of those in Fig. 6). On a $(b/r_{0})^{4}$ scale the lattice results for $r_{0}^{2}M_{PS}^{2}$ lie very close to the continuum limit, which appears to be non-zero by more than five standard deviations, whatever ansatz is taken for the continuum extrapolation. Thus the hypothesis of “no mechanism” is not supported by lattice data.

II.3 Check and irrelevance of finite size effects

Our observables in the NG phase were evaluated on lattices of linear size $L=2.4$ fm. In order to check the impact of finite-size effects (FSE), for one lattice resolution ($\beta=5.85$, $b\simeq 0.123$ fm) we carried out simulations in the NG phase at precisely the same bare parameters as specified in Table 4 for two additional lattice volumes, namely $L/b=16\leftrightarrow L\simeq 2.0$ fm and $L/b=28\leftrightarrow L\simeq 3.4$ fm. The lattice time extension was kept fixed at $T/b=40$, i.e. $T\simeq 4.8$ fm, and we restricted attention to $\eta=-1.2068$, practically coinciding with $\eta_{cr}$, and two twisted mass values, $b\mu=(0.0040,0.0070)$.

As the results at $L=2.4$ fm indicate a clear evidence in favour of the dynamical generation of a NP terms $\propto c_{1}\Lambda_{S}$ in the effective Lagrangian (L13) of our model, we expect that, owing to this and analogous NP effective action terms, the Goldstone modes of the scalar field $\Phi$ are still coupled to the fermion-gauge sector in spite of working at $\eta_{cr}$. Thus the low energy effective theory should include both these Goldstone modes and the massive PS mesons, with interactions giving rise to FSE due to around-the-world propagation of Goldstone modes and PS mesons. While the latter will be exponentially suppressed at large $L$ (at $L=2.4$ fm for $b\mu=0.0040$ we have $M_{PS}L\simeq 5.8$), the former kind of FSE are expected to decrease only as $1/L^{2}$ and will hence be dominating. For this reason in Figs. 9 and 10 we plot as a function of $(b/L)^{2}$ our results for $bM_{PS}$ and $bm_{AWI}$, together with their linear extrapolation to infinite $L$. The use of the simple extrapolation ansatz $K_{0}+K_{1}L^{-2}$ is justified a posteriori by the relative smallness of the observed FSE, well below 3% for $M_{PS}$ and 10% for $m_{AWI}$ in the whole $L^{-2}$–range. The error bars on the simulated data points are merely statistical, plus those entailed by the extrapolation in $L^{-2}$ for the infinite volume estimates.

Although given these errors the statistical significance for FSE is not impressive (especially in the case of $M_{PS}$), still we do observe a uniform trend for both $M_{PS}$ and $m_{AWI}$, which is consistent with a slight linear decrease as a function of $L^{-2}$. As anticipated in Sect. I.4, the deviation of $M_{PS}$ and $m_{AWI}$ evaluated at $L\simeq 2.4$ fm from their infinite volume counterparts can be conservatively estimated to be around 1-2% and 5-7%, respectively.

FSE corrections of this size are (in particular for the matrix element ratio named $m_{AWI}$, see Eq. (L17)) a bit larger than those one should observe (see e.g. Ref. Colangelo et al. (2005)) in a QCD-like theory with scalars completely decoupled, but seem to be of a reasonable magnitude for a model with Goldstone modes non-perturbatively coupled to the fermions owing to a non-zero $c_{1}$. The latter can be roughly estimated to be of order 0.1. In fact on the basis of the $M_{PS}$ results quoted in Eq. (20) $c_{1}$ should describe a NP fermion mass close to one third of the strange quark mass (if we calibrate our model with $r_{0}$ from QCD).

Anyway the observed finite-$L$ corrections on the results for $M_{PS}$ and $m_{AWI}$ at $L\simeq 2.4$ fm are well below the uncertainties shown in Figs. 2 and 3 of the main text, or those quoted in Eqs. (20)–(21).

A more detailed understanding of the FSE on specific hadron masses and operator matrix elements is beyond the scope of the present work. It is however important to remark that the evidence for a NP terms $\propto c_{1}\Lambda_{S}$ in the effective Lagrangian of our critical model can be perfectly well established in finite volume, provided the system volume ($L^{3}T$) is large enough to allow for vacuum orientation under the phenomenon of (dynamical) spontaneous breaking of the $\tilde{\chi}$ flavour symmetries.

In our case we can check the large volume condition by looking at the value of $M_{PS}L$ for $L=2.4$ fm, or equivalently $L/r_{0}=4.8$, and $M_{PS}\sim 0.93/r_{0}$, as obtained in the limits $\mu\to 0$ and $b^{2}\to 0$, see Eq. (20). In this way we find $M_{PS}L\sim 4.5$, a confortably large value which guarantees that the vacuum of the system under the (strong interaction driven) spontaneous breaking of the $\tilde{\chi}$ flavour symmetries gets actually polarized by the dominating effective action terms breaking explicitly these symmetries. For this reason, in the results of Eqs. (20)–(21), which establish the mass generation mechanism at $L=2.4$ fm, we can safely ignore the finite size corrections, which also turn out to be much smaller than other errors.

II.4 Studies at larger values of $\rho$ and $v$

Once a fermion mass term breaking $\tilde{\chi}_{L}\times\tilde{\chi}_{R}$ symmetry is found in the continuum limit of the critical theory in its NG phase, as conjectured in Ref. Frezzotti and Rossi (2015), it is important to check the behaviour of this mass as a function of $\rho$ and the scalar field expectation value $v$.

Since as $\rho\to 0$ also $\eta_{cr}$ vanishes Frezzotti and Rossi (2015), at $\rho=0$ the $\tilde{\chi}_{L}\times\tilde{\chi}_{R}$ transformations become exact symmetries of the lattice model (1) and $m_{AWI}$ vanishes, as well as the infinite volume value of $M_{PS}$. One thus expects that $m_{AWI}$ and $M_{PS}^{2}$ are increasing functions of $\rho$, and arguments were given in Ref. Frezzotti and Rossi (2015) according to which for small $\rho$ they should increase at least as O($\rho^{2}$). The main results of this work were produced at $\rho=1.96$ and $r_{0}v_{R}=1.207(3)$, for lattice spacings where the combination $\rho bv$, which roughly speaking plays here a role analogous to the one of the $r$-parameter for Wilson fermions as far as the magnitude of the $\tilde{\chi}_{L}\times\tilde{\chi}_{R}$ breaking is concerned, ranges from 0.5 to 0.7.

In order to check the $\rho$ dependence of $m_{AWI}$ and $M_{PS}^{2}$ we have thus studied the case of a 1.5 times larger $\rho$ value, namely $\rho=2.94$, for the intermediate lattice spacing ($\beta=5.85$, $a\sim 0.123$ fm) and for the same scalar sector parameters as specified in Sect. I.4.

The list of simulations we performed at $\rho=2.94$ is given in Tab. 5. We limited ourselves to a few simulations since due to the larger $\rho$ value the CPUtime effort for inverting the Dirac matrix for similar volumes and twisted mass values increased by a factor 5 to 10 with respect to the case of $\rho=1.96$.

In the Wigner phase, using the data analysis procedure of Sect. II.1, we found that changing $\rho$ from $1.96$ to $2.94$ leads to a substantial change of $\eta_{cr}$, from $-1.207(8)$ to $-1.838(13)$. Notice that the increase of $\eta_{cr}$ is almost perfectly linear in $\rho$, in line with the theoretical remark Frezzotti and Rossi (2015) that $\eta_{cr}$ is an O($\rho$) quantity, and an odd function of $\rho$.

In the NG phase we observed a significant increase of $2r_{0}m_{AWI}Z_{\tilde{A}}Z_{P}^{-1}$ from $1.15(0.21)$ to $2.80(0.36)$ (see Fig. 11) as well as a similar increase of $M_{PS}^{2}r_{0}^{2}$. These findings are in nice agreement with theoretical expectations from Ref. Frezzotti and Rossi (2015). For their interpretation in the context of a revised concept of universality in more realistic models based on the NP mass mechanism under discussion, see the remarks at the end of the section “Lattice study and results” and in Sect. VI.

Due to CPU time limitations we leave for a future study the investigation of the behaviour of $m_{AWI}$ and $M_{PS}^{2}$ in the critical theory as a function of the scalar field expectation value $v$, which is motivated by the remarks in the main text (section “Nambu–Goldstone phase and NP anomaly”) on the dependence of $c_{1}$ on $\hat{\mu}_{\phi}^{2}$.

II.5 Summary of numerical results

For completeness we summarize in the Tables 7 (for the Wigner phase) and 8 (for the NG phase) all the simulation data we have obtained, which were used to produce the results and the figures presented in the paper. In the Table 9 we also provide the data from the additional simulations at $\eta=\rho=0$ that were carried out in view of the consistency checks discussed in the section V.1.

III $m_{\phi}^{2}$ independence of $\eta_{cr}$

The critical value of the bare Yukawa coupling, $\eta_{cr}$, is defined as the value of $\eta$ that solves the equation $\eta-\bar{\eta}(\eta;g_{0}^{2},\rho,\lambda_{0})=0$. The dimensionless coefficient function $\bar{\eta}(\eta;g_{0}^{2},\rho,\lambda_{0})$ characterizes the UV power divergent mixing of the $d=6$ operator ${O}_{6}^{Li}$ with the $d=4$ density $\tilde{D}_{L}^{i}$. Here $\tilde{D}_{L}^{i}$ and ${O}_{6}^{Li}$ stand for the operators that are obtained from the variation of the ${\cal L}_{Y}$ and ${\cal L}_{W}$ terms in Eq. (1), respectively, under an isotriplet (index $i$) $\tilde{\chi}_{L}$ transformation. The expression of $\tilde{D}_{L}^{i}$ is given in Eq. (L6), while the detailed form of ${O}_{6}^{Li}$ can be found in Sect. III D of Ref. Frezzotti and Rossi (2015).

In what follows we explain why the dimensionless mixing coefficient $\bar{\eta}$, and consequently also $\eta_{cr}$, is independent of the scalar squared mass parameter $m_{\phi}^{2}$, or more precisely of its subtracted counterpart $\mu_{\phi}^{2}\equiv m_{\phi}^{2}-m_{cr}^{2}$ (we neglect here a dimensionless renormalization factor, as it is irrelevant for our arguments). This property is a consequence of renormalizability, locality and straightforward dimensional considerations. We start by observing that in a correlator of renormalized fields with one insertion of ${O}_{6}^{Li}$ the leading UV divergency is given by quadratically divergent Feynman diagrams and is insensitive to the shift $\mu_{\phi}^{2}\to\mu_{\phi}^{2}+\delta m^{2}$. Indeed, if we imagine expanding in $\delta m^{2}$ these diagrams, we see that for dimensional reasons the terms in the expansion beyond the first one will be necessarily suppressed by at least one relative power of $b^{2}\delta m^{2}$. Hence they do not contribute to $\bar{\eta}$.

At the NP level one might worry about the possibility that $\mu_{\phi}^{2}$ enters $\bar{\eta}$ scaled by $\Lambda_{S}^{2}$. We now prove that the presence of either positive or negative powers of the ratio $\mu_{\phi}^{2}/\Lambda_{s}^{2}$ in $\bar{\eta}$ would violate the general principles of renormalizability in Quantum Field Theories, which are assumed to hold in our ${\cal L}_{toy}$ model.

1. Positive mass powers - Terms of the kind $\mu_{\phi}^{2}/\Lambda_{s}^{2}$, i.e. terms displaying inverse powers of the RGI scale of the theory, cannot appear in a renormalizable field theory. Apart from the fact that this kind of dependence would completely ruin the Feynman diagrams power counting, from the argument given above and renormalizability to the NP level in $g_{0}^{2}$ it follows that corrections $\sim\mu_{\phi}^{2}/\Lambda_{s}^{2}$ to $\bar{\eta}$ can never occur.

In fact, one explicitly checks that any correlator of renormalized fields with one insertion of ${O}_{6}^{Li}$ can be expanded in powers of $\eta$, $\rho$ and $\lambda_{0}$ yielding a series where each term corresponds to a diagram consisting only of renormalized effective vertices with external legs and/or scalar lines sewn with a number of scalar field propagators. In all terms of this series the $\mu_{\phi}^{2}$ dependence comes only from the scalar propagators, since the renormalized effective vertices, which by construction here resum diagrams with gluon and fermion (but no scalar) internal lines, are independent of $\mu_{\phi}^{2}$, up to immaterial O($b^{2}$) artifacts. Hence, under the shift $\mu_{\phi}^{2}\to\mu_{\phi}^{2}+\delta m^{2}$, a correction to the correlator of the kind $\delta m^{2}/\Lambda_{s}^{2}$ does not occur. Indeed, were such a correction present, it could only arise as the product of a relative O($b^{2}\delta m^{2}$) correction from the scalar propagators times an O($b^{-2}\Lambda_{S}^{-2}$) contribution from the $\mu_{\phi}^{2}$-independent renormalized effective vertices. But, as the renormalized effective vertices are by definition UV finite, this is manifestly impossible, thereby ruling out terms of the kind $\mu_{\phi}^{2}/\Lambda_{s}^{2}$ in the correlator.

2. Negative mass powers - Terms of the kind $\Lambda_{s}^{2}/\mu_{\phi}^{2}$, i.e. terms displaying inverse powers of mass parameters cannot show up in mixing coefficients. In fact, IR divergencies, revealed in the limit of vanishingly small mass, can only appear in matrix elements and not in the expression of the renormalized operators in terms of the bare ones. This is the consequence of the fact that one can determine the mixing coefficients by only using correlators at non-exceptional momenta that do not develop IR divergencies in the massless limit.

Finally, as a check of consistency, we repeated the numerical computation of $\eta_{cr}$ at $\beta=5.95$ reducing the subtracted squared scalar mass, $r_{0}^{2}\mu_{\phi}^{2}\equiv r_{0}^{2}(m_{\phi}^{2}-m_{cr}^{2})$, from $1.22$ to $0.41$, see Table 6. For this smaller value of the subtracted scalar mass we worked at the same bare parameters and lattices as employed for the larger one, see the block of Tab. 2 corresponding to $\beta=5.95$. Adopting the same analysis procedure as discussed in Sect. II.1, we find that the change in the estimated value of $\eta_{cr}$ observed upon reducing $\mu_{\phi}^{2}$ by a factor of three is not significant within our statistical errors. Moreover the change in the central $\eta_{cr}$ value obtained with the two different values of $\mu_{\phi}^{2}$ in Table 6 but the same choice of time plateau for $r_{AWI}$ is at a given $\mu_{\phi}^{2}$ as small as the central value of $\eta_{cr}^{A}-\eta_{cr}^{B}$, which is certainly an O($b^{2}$) lattice artifact in the computation of $\eta_{cr}$. These findings altogether are nicely consistent with $\eta_{cr}=\eta_{cr}(g_{0}^{2},\rho,\lambda_{0})|_{g_{0}^{2}=6/5.95,\rho=1.96,\lambda_{0}=0.6022}$ being independent of the value of $\mu_{\phi}^{2}$ up to O($b^{2}$) lattice artifacts.

IV Lattice artifacts ${\rm O}(b^{2n})$

We show in this section that in the lattice model (1) only discretization errors of order ${\rm O}(b^{2n})$ ($n$ positive integer) can show up in the continuum limit of renormalized quantities. This property is a consequence of the exact invariance of the lattice action under the global $\chi_{L}\times\chi_{R}$ transformations, the form of which is given by Eqs. (L3)-(L4), with $\Omega_{L}$ and $\Omega_{R}$ independent SU(2) matrices.

The proof is based on the observation that for symmetry reasons the Symanzik local effective Lagrangian Symanzik (1983a, b); Luscher et al. (1996) describing the lattice artifacts of the model (1), with the subtracted scalar squared mass $\mu_{\phi}^{2}\equiv m_{\phi}^{2}-m_{cr}^{2}$ set to an O($b^{0}$) value, cannot involve operator terms with odd engineering mass dimension, thereby taking the form

[TABLE]

In fact, by properly combining parity ($P$), time reversal ($T$), a finite $\chi_{L}\times\chi_{R}$ transformation and a finite vector U(1) transformation, which all correspond to exact lattice symmetries, one can define a lattice transformation which, in analogy with the Wilson Lattice QCD case discussed in Ref. Frezzotti and Rossi (2004), we call ${\cal D}_{d}$. The action of ${\cal D}_{d}$ on the lattice fields reads

[TABLE]

where the SU(3) gauge link transformation implies for the gauge potentials $A_{\mu}(x)\to-A_{\mu}(x^{PT})$, $\mu=0,1,2,3$. One can also check that the gauge covariant derivatives transform so that e.g. $\widetilde{\nabla}_{\mu}Q(x)\to-i\widetilde{\nabla}_{\mu}Q(x^{PT})$, while for the plaquette $U_{\mu\neq\nu}(x)$ we have $U_{\mu\neq\nu}(x)\to U_{\mu\neq\nu}^{\dagger}(x^{PT})$.

We thus see that all the fields entering the basic lattice action (1) under the operation ${\cal D}_{d}$ get transformed in their counterparts localized at the $PT$-reflected coordinate $x^{PT}=-x$ times a phase $(-1)^{d_{F}}$, where $d_{F}$ is the engineering mass dimension of the field $F$, namely $d_{A}=d_{\Phi}=-1$ for gluons and scalars, $d_{Q}=d_{\bar{Q}}=(-1)^{3/2}=i$ for fermion fields. By the arguments of Ref. Luscher et al. (1996), the invariance of the lattice action (1) under ${\cal D}_{d}$ implies an analogous invariance of the formal Symanzik effective Lagrangian $L_{\;\rm toy}^{Sym}$ under the continuum counterpart of ${\cal D}_{d}$, which in turn rules out from it all the local terms with odd engineering mass dimension, leading to Eq. (22).

By using again the exact ${\cal D}_{d}$ invariance and the standard arguments of Ref. Luscher et al. (1996) about the lattice artifact corrections to local operators, it is also clear that the discretization errors on the matrix elements of a local operator $O(x)$ taken between states at momenta parametrically much smaller than $b^{-1}$ can be fully described in terms of the Symanzik effective Lagrangian $L_{\;\rm toy}^{Sym}$ and the effective local operator $O_{0}(x)+\sum_{n=1}^{\infty}b^{2n}O_{2n}$.

The same kind of arguments can be straightforwardly extended to prove that the power divergent mixings of a local operator can only involve even powers of $b^{-1}$.

V Further consistency checks

After discussing in Sects. II to IV several technical aspects concerning the numerical evidence for the NP mass generation mechanism under investigation, here we present evidence that no such a phenomenon occurs if the fermionic chiral transformations $\tilde{\chi}_{L}\times\tilde{\chi}_{R}$ (see Eqs. (L3)–(L4)) are exact in the UV regulated lattice model, i.e. for the case of $\rho=0$ and $\eta=\eta_{cr}(g_{0}^{2},0,\lambda_{0})=0$. Since the scalar $\Phi$ field is then completely decoupled, this case actually corresponds to quenched lattice QCD with naïve fermions.

We also show that in the Wigner phase of the critical model corresponding to $\rho=1.96$ the pseudoscalar (PS) meson squared mass in the combined chiral and continuum limits vanishes within errors, as expected from the arguments given in the main text (Sect. 2) and in Ref. Frezzotti and Rossi (2015).

V.1 The critical model at $\rho=\eta=0$

In this section we present a study of the model obtained by taking $\eta=0$, $\rho=0$ with the purpose of showing that in the absence of chirally breaking terms the dynamical mass generation phenomenon identified in this work does not take place. Indeed, in this situation we get quenched lattice QCD with the scalars completely decoupled from (naïve) fermions and gluons. We thus expect both $M_{PS}$ and $m_{AWI}$ to vanish in the infinite volume limit. As in the case of the toy-model (1) discussed before, reliable simulations at $\eta=0$, $\rho=0$ require introducing a twisted mass term, like Eq. (12), that can be safely removed by taking the limit $\mu\to 0$ at fixed lattice spacing.

For the sake of this check we will limit our investigation to the finest lattice resolution at our disposal, namely $\beta=5.95$, corresponding to $b\simeq 0.102(6)$ fm, and consider seven (small) values of the twisted quark mass $\mu$. Details about simulation parameters (chosen as in the study of the model with non-zero $\tilde{\chi}$ symmetry breaking) and numerical results are given in Table 9.

In Figs. 12 and 13 we compare the values of $2r_{0}m_{AWI}$ and $r_{0}^{2}M_{PS}^{2}$, respectively, at $\rho=\eta=0$ with the values obtained for the same lattice resolution ($\beta=5.95$) at $\rho=1.96$ and $\eta=\eta_{cr}$ in the NG phase (Tables 4 and 8).

First of all we see that in the case $\eta=0$, $\rho=0$ $m_{AWI}$ vanishes for all values of $b\mu$. This is a somewhat trivial result, however, because actually the twisted mass term (12) does not contribute to the SDE associated with the $\tilde{\chi}_{L}\times\tilde{\chi}_{R}$ transformations (L11) for $i=1$, from which $m_{AWI}$ is extracted at non-zero $\mu$–values (see Eq. (17)). More interesting is the case of $M_{PS}^{2}$. The value of $r_{0}^{2}M_{PS}^{2}$ at $\rho=\eta=0$ extrapolated to $\mu=0$ is $\sim 20$ times smaller than the value we found at $\eta=\eta_{cr}$, $\rho=1.96$.

While in the latter case the data from simulations at several $\eta$- and $\mu$-values at $\beta=5.95$ (Table 4) were exploited and fit according to the ansatz $M_{PS}^{2}=Y_{0}+Y_{1}\eta+Y_{2}\mu+Y_{3}\eta^{2}+Y_{4}\mu^{2}+Y_{5}\eta\mu\,$, for the case of $\eta=\rho=0$ the extrapolation to $\mu=0$ yields $r_{0}^{2}M_{PS}^{2}|_{\mu=0;\beta=5.95}=0.041(11)$ and is performed by fitting our data (with $M_{PS}L\geq 4$) to a third order polynomial in the twisted mass parameter $\mu$. Following several authors Butler et al. (1994); Bernard et al. (1998); Blum et al. (2004) who studied the meson spectrum in the quenched approximation with similar lattice sizes and statistics, we choose for $M_{PS}^{2}$ as a function of the quark mass $\mu$ a polynomial fit ansatz that reproduces the LO behaviour expected in Chiral Perturbation Theory and effectively describes higher order chiral terms and lattice artifacts proportional to $b^{2}\mu^{2}\Lambda_{S}$ or $b^{2}\mu^{3}$. We do not try to disentangle quenched and/or standard “chiral log” effects, a task that is beyond the scope of our work and typically requires higher statistics than we have. For instance, using an effective statistics about 30 times larger (due to more gauge configurations and bigger lattices) than here the authors of Ref. Aoki et al. (2000) could estimate only with a 50–60% relative uncertainty the parameter $\delta\sim 0.1$ of the anomalous quenched “chiral log” Bernard and Golterman (1992); Sharpe (1992) behaviour $M_{PS}^{2}\sim(2B\mu)^{1/(1+\delta)}$.

Although lattice symmetries here imply that in infinite volume one should get $M_{PS}^{2}\to 0$ as $\mu\to 0$, the marginally significant deviation from zero that we observe in $r_{0}^{2}M_{PS}^{2}|_{\mu=0;\beta=5.95}=0.041(11)$ is likely to arise from a combination of statistical uncertainties and small residual finite size and/or quenched “chiral log” effects. Our result appears anyway to be well in line with what was obtained (with larger errors on a bit smaller lattices) in previous studies of quenched twisted mass lattice QCD, see e.g. Refs. Jansen et al. (2005a, b), which for $L=1.6$ fm at fixed (very similar) lattice spacing find $r_{0}^{2}M_{PS}^{2}|_{\mu=0;\beta=6}=0.014(58)$.

Most importantly for the goals of the present work, if one assumes a systematic absolute uncertainty of $0.04$ in the extrapolation of $r_{0}^{2}M_{PS}^{2}$ to $\mu=0$ the resulting error on $r_{0}M_{PS}$ amounts to only $\sim 2.5\%$ of our final result for $r_{0}M_{PS}$ in the NG phase, see Eq. (20), and is thus four times below the relative systematic error quoted there. The present study of $m_{AWI}$ and $M_{PS}^{2}$ in the model at $\rho=\eta=0$ provides hence a nice consistency check of our results in the main text and Sect. II.

The findings above are also in agreement with the numerical evidence presented in Sect. II.4 that, when the non-zero parameter $\rho$ that controls the explicit breaking of the $\tilde{\chi}$ symmetries is increased by a factor 1.5, the NP fermion mass that is dynamically generated in the NG-phase of the critical model increases by a factor $\sim 2.4(3)$, i.e. nearly like $\rho^{2}$, as expected from arguments in Frezzotti and Rossi (2015).

It is finally worth comparing the model at $\eta=\rho=0$ with the model at $\rho=1.96$ and $\eta=\eta_{cr}$ in the Wigner phase. These two choices of parameters are expected to describe the same model in the continuum limit, but at fixed lattice spacing they will differ by lattice artifacts. Indeed Fig. 14 shows that at $\beta=5.95$ the extrapolated values of $r_{0}^{2}M_{PS}^{2}$ at $\mu=0$ do not significantly differ, while the small misalignment of the blue ($\eta=0$, $\rho=0$) and red curves ($\eta=\eta_{cr}$, $\rho=1.96$) should be ascribed to lattice artifacts and/or small (differences in) finite size effects.

V.2 $M_{PS}$ in the Wigner phase of the critical model

Here we show evidence that in the Wigner phase of the quenched critical model at $\rho=1.96$ the quantity $r_{0}^{2}M_{PS}^{2}$ approaches zero within errors in the combined continuum and $\mu\to 0$ limits. Our final error on the extrapolated value of $r_{0}^{2}M_{PS}^{2}$ is not much smaller than it was in the quenched twisted mass lattice QCD study of Refs. Jansen et al. (2005a, b), because in the present work the simulations in the Wigner phase are rather expensive due to $\Phi$-field induced fluctuations of the fermion propagators and have been thus limited to only three or four values of $\mu$ (times three values of $\eta$) per each lattice spacing. This was enough for our main goal of determining $\eta_{cr}$. Although $r_{0}^{2}M_{PS}^{2}$ should vanish as $\mu\to 0$ even at finite lattice spacing, combining data at three different lattice resolutions effectively increases the statistical information and allows for a better control of the $\mu$-dependence in the available PS-meson mass range.

In this consistency check we employ PS-meson mass data coming from our Wigner phase simulations (satisfying $M_{PS}L\geq 5.3$), which for each lattice spacing are extrapolated to $\eta=\eta_{cr}$ linearly in $(\eta-\eta_{cr})^{2}$. Even if performing an extrapolation linear in $(\eta-\eta_{cr})$ would yield equivalent numerical results within our present errors, from a theoretical viewpoint in the phase with zero vacuum expectation value of $\Phi$ the two–points PS–meson correlator and all derived quantities are expected to admit a Taylor expansion in powers of $(\eta-\eta_{cr})$ with only terms of even integer degree.

In Fig. 15 we plot the critical model lattice data for $r_{0}^{2}M_{PS}^{2}$ at the three lattice spacings and several values of the twisted mass as a function of $r_{0}\mu_{R}=r_{0}\mu Z_{P}^{-1}$, where $\mu_{R}$ is the renormalized twisted mass and $Z_{P}^{-1}=G_{PS}^{W}r_{0}^{2}=\langle 0|P^{1}|{\rm PS~{}meson}\rangle^{W}$ (see main text, Sect. “Lattice study and results”). The $r_{0}^{2}M_{PS}^{2}$ data are then fit to the ansatz $K_{0}+K_{1}\mu_{R}r_{0}+\tilde{K}_{1}\mu_{R}b^{2}r_{0}^{-1}+K_{2}\mu_{R}^{2}r_{0}^{2}$, while the errors on $bM_{PS}$, $r_{0}/b$ and $Z_{P}^{-1}$ at each $\beta$ are propagated taking correlations into account through a jackknife analysis.

The resulting best fit describes well the data and is shown for each lattice resolution (continuous curves) and in the continuum limit (dotted curve), where at $\mu_{R}=0$ we obtain $r_{0}^{2}M_{PS}^{2}=0.032(51)$. The deviation from zero is statistically not significant, but does not exclude small finite size effects and/or quenched “chiral log” corrections (see also the discussion in Sect. V.1).

VI Towards a realistic model

We extend here the discussion started in the section “Conclusions and outlook” of the main text about how one can possibly build particle physics models based on the NP mass generation mechanism discussed and identified in this paper. In order to get a phenomenologically sensible theory one has to extend the model (L1) by introducing, besides EW interactions, also a “Tera-strong” force and a set of “Tera-fermions” (we take this name from Ref. Glashow (2005)) that communicate with ordinary matter (quarks or leptons) via strong and/or EW interactions Frezzotti and Rossi (2018). The “Tera-strong” force should be given by a new non-Abelian gauge interaction that becomes strong at the scale $\Lambda_{T}\gg\Lambda_{QCD}$, where $\Lambda_{T}$ can be roughly estimated to be in the few TeV range Frezzotti and Rossi (2015, 2018).

Models of this type without tree-level flavour changing neutral currents can be easily formulated, if SU(2)${}_{L}\times$U(1)Y invariant EW interaction terms are introduced in the same way as in the SM, e.g. by coupling SU(2)L gauge bosons only to left-handed Tera-fermions arranged in doublets of definite hypercharge. By then taking all right-handed Tera-fermions to be weak isospin singlets with given hypercharge the usual relation $Q_{em}=T_{3}+Y/2$ between electric charge $Q_{em}$, third weak isospin component $T_{3}$ and hypercharge $Y$ is always fulfilled. Owing to the sizeable NP mass of the “Tera-fermions”, the lightest Tera-mesons are expected to be rather heavy, i.e. with a mass around 2-3 times $\Lambda_{T}$, implying that these extended models have a good chance of passing the EW precision tests, in particular those concerning the $S$-parameter Peskin and Takeuchi (1990, 1992). By treating the chiral SU(2)${}_{L}\times$U(1)Y interactions in perturbation theory and carrying out unquenched lattice simulation of models including QCD and “Tera- strong” gauge forces, it is in principle possible to relate with controllable errors the NP mass of the top quark to the mass of the lightest “Tera-mesons”, which should be seen as resonances in accelerator experiments at sufficiently high energy.

Lattice methods can also be useful to check and possibly establish the existence in the $WW+ZZ$ channel of one composite Higgs state with mass $M_{h}\sim 125$ GeV, which is expected to form owing to the interactions of weak gauge bosons with “Tera-fermions” that are subjected to EW and “Tera-strong” interactions. For a first step in this direction see Ref. Romero-Lopez et al. (2018). If the existence of such a composite Higgs state is confirmed, then perturbative unitarity of the (longitudinal) $WW$ scattering amplitude at momentum scales well below $\Lambda_{T}$, for which the Higgs state plays a crucial role Lee et al. (1977), would make very plausible that the Higgs couplings are close to the present experimental findings and deviate from their SM value only by small O($p^{2}/\Lambda_{T}^{2}$) effects, with $p^{2}$ standing for $M_{h}^{2}$ or the typical momentum scale of the process where the Higgs couplings are measured.

A striking and interesting feature of our mass generation mechanism is the observed dependence of low energy physical quantities on the parameters of the $d\geq 6$ $\tilde{\chi}$ breaking Lagrangian terms. To all orders in the perturbative expansion in the gauge coupling these terms would represent irrelevant details of the UV completion of the critical model. But at the NP level, in the NG phase of the critical theory, the $d\geq 6$ $\tilde{\chi}$ breaking operators appearing in the UV regulated classical Lagrangian produce relevant O($\Lambda_{UV}^{0}$) effects, because the symmetries of our model allow for the occurrence in the EL of RGI finite NP terms, such as those with coefficients $c_{1}$ and $c_{2}$ in Eq. (L13). As first pointed out in Ref. Frezzotti and Rossi (2015) and discussed in the section “Nambu–Goldstone phase and NP anomaly”, this dynamical phenomenon is due to a subtle interplay between the residual O($b^{2}$) terms violating $\tilde{\chi}_{L}\times\tilde{\chi}_{R}$ and the strong interactions, which are responsible for spontaneous breaking of the approximate $\tilde{\chi}_{L}\times\tilde{\chi}_{R}$ symmetries. Of course no such a phenomenon would be observed if the critical model were defined in a UV regularization preserving exactly $\tilde{\chi}_{L}\times\tilde{\chi}_{R}$ invariance, e.g. at $\rho=0$ and $\eta=\eta_{cr}|_{\rho=0}=0$. This remark, which is indeed consistent with our numerical finding that $c_{1}$ increases upon increasing $\rho$, see Sect. II C, also implies that non-Abelian gauge models including fermions and scalars in general fall into different universality classes depending on whether the UV regularization breaks or possibly preserves the purely fermionic chiral symmetries.

In other words, the “NP anomaly” to the recovery of broken fermionic chiral symmetries that we demonstrated in our work appears to require a non-trivial revision of the concept of universality, as compared to its standard formulation which can be proven to all orders of perturbation theory (only). This important point certainly deserves a more thorough investigation. At the moment it seems to us conceivable that the kind of universality that is realized in the context of maximal $\tilde{\chi}$ restoration carries over also to the NP mass terms in the EL Frezzotti and Rossi (2018). In this new NP framework universality should thus mean that a set of UV complete models endowed with different $\tilde{\chi}$ breaking Lagrangian terms lead to the same low energy observables, once the conditions of maximal $\tilde{\chi}$ restoration are fulfilled at correspondingly different critical values of the relevant bare parameters. The low energy physics should thus be controlled only by the gauge couplings, the number of fermions and their transformation properties under the gauge groups, plus possibly a finite number of dimensionless coefficients parametrizing the relative strength of the $\tilde{\chi}$ breaking terms for the different fermion generations and having a magnitude limited by the conditions of maximal $\tilde{\chi}_{L}$ restoration.

We finally note that a beyond the SM model supporting the mass generation mechanism we have discussed and endowed with a reasonable content of (standard and Tera) matter fields, if a suitable, non standard choice of hypercharges for Tera-fermions is made, can yield unification of EW and strong (and perhaps even Tera-strong) gauge couplings Frezzotti et al. (2016) at an energy scale around $10^{18}$ GeV and to a precision level that is as good as for the unification achieved in the Minimal SUSY Model.