The chiral condensate from renormalization group optimized perturbation
aa r X i v : . [ h e p - ph ] S e p The chiral condensate from renormalization group optimized perturbation
Jean-Lo¨ıc Kneur and Andr´e Neveu Laboratoire Charles Coulomb (L2C), UMR 5221 CNRS-Universit´e de Montpellier, Montpellier, France
Our recently developed variant of variationnally optimized perturbation (OPT), in particularconsistently incorporating renormalization group properties (RGOPT), is adapted to the calculationof the QCD spectral density of the Dirac operator and the related chiral quark condensate h ¯ qq i inthe chiral limit, for n f = 2 and n f = 3 massless quarks. The results of successive sequencesof approximations at two-, three-, and four-loop orders of this modified perturbation, exhibit aremarkable stability. We obtain h ¯ qq i / n f =2 (2 GeV) = − (0 . − . , and h ¯ qq i / n f =3 (2 GeV) = − (0 . − . where the range spanned by the first and second numbers (respectively four-and three-loop order results) defines our theoretical error, and ¯Λ n f is the basic QCD scale in the MS -scheme. We obtain a moderate suppression of the chiral condensate when going from n f = 2 to n f = 3. We compare these results with some other recent determinations from other nonperturbativemethods (mainly lattice and spectral sum rules). PACS numbers: 12.38.Aw, 12.38.Lg, 12.38.Cy
I. INTRODUCTION
The chiral quark condensate h ¯ qq i plays a central rˆole in QCD nonperturbative dynamics, being together withthe pion decay constant the other principal (lowest dimensional) order parameter of spontaneous chiral symmetrybreaking, SU ( n f ) L × SU ( n f ) R → SU ( n f ) V for n f massless quarks, the physically relevant cases being n f = 2 or 3.It is considered a nonperturbative quantity by excellence, in the sense that it is trivially vanishing at any finite orderof (ordinary) perturbative QCD in the chiral (massless quark) limit, as we recall in more details below.There is a long history of its determination from various models and analytic methods, the best known being theGell-Mann–Oakes-Renner (GMOR) relation [1], relating the quark condensate to the pion mass, the decay constant F π , and the light current quark masses m u,d , typically for the degenerate two-flavor case: F π m π = − ( m u + m d ) h ¯ uu i + O ( m q ) . (1.1)Nowadays the light current quark masses can be precisely extracted from lattice simulations [2] or spectral sumrules (see e.g. [3]), giving from using the GMOR relation above an indirect precise determination of the conden-sate. However, as indicated relation (1.1) is valid upon neglecting possible higher order terms O ( m q ). Indeed, theGMOR relation entails explicit chiral symmetry breaking from current quark masses. Since the condensate is morea property of the QCD vacuum in the strict chiral limit, it is highly desirable to obtain a possible ‘first principle’determination of this dynamical quantity in the strict chiral limit, to disentangle from quark current mass effects.An early analytic determination was in the framework of the Nambu–Jona-Lasinio (NJL) model [4] and its variousextensions as a low-energy effective model of QCD, valid at least for the gross features of chiral symmetry breakingproperties. In the NJL model, the condensate is evaluated in the strict chiral limit, or taking into account explicitbreaking from small current quark masses, in the leading (large- N ) approximation [5] as a function of the physicalcutoff and other parameters of the model to be fitted from data. There are also related other attempts based onanalytic methods, like Schwinger-Dyson equations [6, 7] typically. Phenomenological values of the condensate canalso be extracted [3, 8] indirectly from data using spectral QCD sum rule methods [9], where the quark condensateand other higher-dimensional condensates enter as nonperturbative parameters of the operator product expansion ininverse powers of momenta.More recently ‘ab initio’ lattice calculations have determined the quark condensate by several independent ap-proaches, some actually related to the GMOR relation, or some more direct determinations with different methods [10](see [2] for a review of various recent lattice determinations). In particular after early pioneering work [11], therehas been more recently a renewed intense interest in computing on the lattice the spectral density of the Diracoperator [12, 13], directly related to the quark condensate through the Banks-Casher relation [14]. However, whilemany recent lattice results are statistically very precise, lattice determinations rely in the end on extrapolations to thechiral limit, often using chiral perturbation theory [15] input for that purpose. Earlier general results on the spectraldensity in the nonperturbative low eigenvalue range had been obtained in [16], and attempts to use chiral perturbationtheory information was elaborated e.g. in [17]. The link between different definitions of the quark condensate, inparticular between the spectral density and the condensate appearing in the operator product expansion (OPE), hasbeen carefully discussed in [7].On more phenomenological grounds there are long-standing questions also on the dependence of the quark conden-sate on the number of flavors ( n f = 2 or n f = 3 for the physically relevant cases). In particular it had been advocatedthat the GMOR relation may receive substantial corrections, and some authors indeed found significant suppressionof the three-flavor case with respect to the two-flavor case [18], which may be attributed to the relatively large explicitchiral breaking from the not-so-small mass of the strange quark, roughly of order m s < ∼ ¯Λ QCD /
3. There are alsosome hints that chiral perturbation theory might not converge very well for n f = 3 [18, 19], for which value latticesimulations also show large discrepancies between different collaborations and methods (with some results with, andsome other results without relative suppression of the n f = 3 condensate [2]).Our recently developed renormalization group optimized perturbation (RGOPT) method [20–22] appears particu-larly adapted to estimate this quantity, since it gives a nonperturbative sequence of approximations starting from apurely perturbative expression. This allows by construction an analytic ‘first principle’ calculation of this quantityand gives a nontrivial result by construction in the chiral limit, in contrast with the standard perturbative one (seealso [23, 24] for earlier attempts in that direction). Morever this also gives us a very simple analytic handle on theexploration of the chiral limit for arbitrary number of flavors (2 or 3 in practice), which in our framework is simplycontained in the known flavor dependence of the first few perturbative coefficients, while this appears more difficultat present both for chiral perturbation theory and lattice simulations.The paper is organized as follows. In section II we shortly recall basics of the spectral density and its Banks-Casherconnection with the condensate. In section III we recall the main OPT method and our RGOPT version incorporatingconsistent RG properties. We adapt the RGOPT to the spectral density case in section III.C. Section IV is a digressionwhere we first consider the spectral density RGOPT calculation in the Gross-Neveu model, where it can be comparedwith the exact result for the fermion condensate in the large- N limit, known from standard methods. Section V dealswith the actual computation of the optimized spectral density in QCD at the three presently available orders (two,three and four loops) of the variationally modified perturbation. Detailed numerical results are presented as well assome comparison with other recent determinations of the quark condensate. Finally section VI is a conclusion. II. SPECTRAL DENSITY AND THE QUARK CONDENSATE
We shall just recall in this section some rather well-known features of the spectral density and its connection withthe chiral condensate, known as the Banks-Casher relation [14] (see also e.g. [16] for more details), to be exploitedbelow. We thus start from the (Euclidean) Dirac operator which formally has eigenvalues λ n and eigenvectors u n ,i /D u n ( x ) = λ n u n ( x ); /D ≡ /∂ + g /A , (2.1)where /D is the covariant derivative operator and A the gluon field. Except for zero modes, the eigenvectors comein pairs { u n ( x ); γ u n ( x ) } , with respective eigenvalues { λ n ; − λ n } that depend on A . In the discrete case ( i.e. on alattice with finite volume V ), by definition the spectral density is given by ρ ( λ ) ≡ V h X n δ ( λ − λ [ A ] n ) i , (2.2)where δ ( x ) here is the Dirac distribution and h· · ·i designates averaging over the gauge field configurations, hi = Z [d A ] N Y i =1 det(i /D + m ) . (2.3)The quark condensate is given by 1 V Z V d x h ¯ q ( x ) q ( x ) i = − mV X λ n > λ n + m . (2.4)Now when V goes to infinity the operator spectrum becomes dense, so that h ¯ qq i = − m Z ∞ d λ ρ ( λ ) λ + m , (2.5)where ρ ( λ ) is the spectral density.The Banks-Casher relation is the m → m → h ¯ qq i = − πρ (0) , (2.6)if the spectral density at the origin can be known. This is an intrinsically nonperturbative quantity, vanishing to allorders of ordinary perturbation, just as the left-hand side of this last equation. Now taking into account that fornon-zero fermion masses m , h ¯ qq i = h ¯ qq i ( m ) ≡ − Σ( m ) we have from the defining relations (2.2), (2.5) the followinginteresting tautology, ρ ( λ ) = 12 π [Σ(i λ + ǫ ) − Σ(i λ − ǫ )] | ǫ → , (2.7) i.e. ρ ( λ ) is determined by the discontinuities of Σ( m ) across the imaginary axis. The interest of this relation isthat, when the quark mass is nonzero, h ¯ qq i has a purely perturbative series expansion, known to three-loop orderat present, and its discontinuities are simply given by those coming from the perturbative, purely logarithmic massdependence. Therefore it makes sense to calculate a perturbative spectral density using the above relation. Usually itwill be of little use, since taking its λ → δ expansion (see below), is precisely the analytic handle giving a nontrivial result for λ →
0, just as itgives a nontrivial result for m → ρ ( λ ) up tofour loops, using the logarithmic discontinuities involved in Eq.(2.7). Then we perform a variational transformation(so-called δ -expansion, defined in next section) on the perturbative series, and solve appropriate optimization OPTand renormalization group (RG) equations (to be precisely defined in next section) to derive a nontrivial, optimizedvalue of the relation (2.6). III. OPTIMIZED AND RG OPTIMIZED PERTURBATION (RGOPT)A. Standard OPT
The key feature of the optimized perturbation (OPT) method (appearing in the literature under many names andvariations [25]) is to introduce an extra parameter 0 < δ <
1, interpolating between L free and L int for any Lagrangian,in such a way that the mass parameter m is traded for an arbitrary trial parameter. This is perturbatively equivalentto taking any standard perturbative expansions in the coupling g ( µ ), after renormalization in some given scheme ( e.g. M S -scheme with arbitrary scale m u ), reexpanded in powers of δ after substituting, m → m (1 − δ ) a , g → δ g . (3.1)Such a procedure is consistent with renormalizability [23, 24, 26] and gauge invariance [24], whenever the latter isrelevant, provided of course that the above redefinition of the coupling is performed consistently for all interactionterms and counterterms appropriate for renormalizability and gauge invariance, as is the case for QCD. In (3.1) weintroduced an extra parameter a , to reflect a certain freedom in the interpolation form, which will be crucial toimpose compelling RG constraints, as discussed below and in our previous work [21, 22]. Applying (3.1) to a givenperturbative expansion for a physical quantity P ( m, λ ), reexpanded in δ at order k , and taking afterwards the δ → massless theory, leaves a remnant m dependence at any finite δ k order. The arbitrarymass parameter m is then most conveniently fixed by an optimization (OPT) prescription, ∂∂ m P ( k ) ( m, g, δ = 1) | m ≡ ˜ m = 0 , (3.2)which generally determines a nontrivial optimized mass ˜ m ( g ), having a nonperturbative g dependence, realizingdimensional transmutation. (More precisely, for asymptotically free theories, the optimized mass is automaticallyof the order of the basic scale Λ ∼ µ e − / ( b g ) , in contrast with the original vanishing mass). In simpler ( D = 1)models the procedure may be seen as a particular case of “order-dependent mapping”[27], which has been proven[28]to converge exponentially fast for the D = 1 Φ oscillator energy levels. For higher dimensional D > δ is more involved, and no rigorous convergence proof exists, although OPT wasshown to partially damp the factorially divergent (infrared renormalons) perturbative behavior at large orders [29].Nevertheless the OPT can give rather successful approximations for nonperturbative quantities beyond mean-fieldapproximations in a large variety of models [25, 30, 31], including studies of phase transitions at finite temperaturesand densities [33, 34]. B. Renormalization group optimized perturbation (RGOPT)
In most previous OPT applications [25], the linear δ expansion is used, namely assuming a = 1 in Eq. (3.1) mainlyfor simplicity and economy of parameters. However a well-known drawback of this conventional OPT approach is that,beyond lowest order, Eq. (3.2) generally gives more and more solutions at increasing orders, many being complex-valued, as a result of exactly solving algebraic equations in g and/or m . This problem is typically encountered first attwo-loop order. In general, without some insight on the nonperturbative behavior of the solutions, it can be difficultto select the right one, and unphysical nonreal solutions at higher orders are embarrassing. As it turns out, RGconsistency considerations provide a compelling way out, as developed in our more recent approach [20–22], whichdiffers crucially from the more conventional OPT based on the linear δ expansion in two main respects. First, itintroduces a straightforward combination of OPT and RG properties, by requiring the ( δ -modified) expansion tosatisfy, in addition to the OPT Eq. (3.2), a standard RG equation, µ dd µ (cid:16) P ( k ) ( m, g, δ = 1) (cid:17) = 0 , (3.3)where the RG operator is defined as usual , µ dd µ = µ ∂∂µ + β ( g ) ∂∂g − γ m ( g ) m ∂∂m . (3.4)Note, once combined with Eq. (3.2), the RG equation takes a reduced simple form, corresponding to a massless theory, (cid:20) µ ∂∂µ + β ( g ) ∂∂g (cid:21) P ( k ) ( m, g, δ = 1) = 0 . (3.5)Thus Eqs. (3.5) and (3.2) together completely fix optimized m ≡ ˜ m and g ≡ ˜ g values.Remark indeed that RG invariance is in general spoiled after the rather drastic modification from (3.1), reshufflinginteraction and free terms from the original perturbative expansion. This feature has seldom been considered andappreciated in previous applications of the OPT based on the linear δ -expansion method to renormalizable theories.Thus RG invariance has to be restored in some manner, accordingly Eq. (3.5) gives an additional non-trivial constraint.Intuitively, just as the stationary point OPT solutions from Eq. (3.2) are expected to give sensible approximations,at successive orders, to the actually massless theory, one similarly expects that combining the OPT with the RGsolutions of (3.5) should further give a sensible sequence of best approximations to the exactly scale invariant allorder result. (N.B.: An earlier way of reconciling the δ expansion with RG properties was used in [23, 24, 26]:Schematically it amounted to resumming the δ -expansion to all orders, which can be done in practice only for thepure RG dependence up to two-loop order. These resummations came as rather complicated integral representations,rendering difficult generalizations to higher orders, other physical quantities or other models of interest. In contrastthe purely perturbative procedure of imposing Eqs. (3.2), (3.5) is a considerable shortcut, straightforward to applyto any model, being based solely on purely perturbative expansions.)Yet applying Eq. (3.2), (3.5) without further insight still gives multiple solutions at increasing orders. So weproposed [21, 22] a compelling selection criterion, by retaining only the branch solution(s) g ( m ) (or equivalently m ( g )) continuously matching the standard perturbative (asymptotically free (AF) RG behavior in the QCD case) forvanishing coupling, namely, ˜ g ( µ ≫ ˜ m ) ∼ (2 b ln µ ˜ m ) − + O ((ln µ ˜ m ) − ) . (3.6)Now the crucial observation is that requiring at least one of the solutions of Eq. (3.5) to satisfy (3.6) implies astrong necessary condition on the basic interpolation (3.1), fixing the exponent a uniquely in terms of the universal(scheme-independent) first order RG coefficients [21, 22]: a ≡ γ b , (3.7)which is the second important difference of the present RGOPT with respect to the standard OPT . For the criticalvalue (3.7), Eq. (3.5) is in fact exactly satisfied at the lowest δ order, therefore giving no further constraint. Athigher δ orders, (3.7) implies that one at least of both the RG and OPT solutions fulfills Eq. (3.6), and solutions withthis behavior are essentially unique (although not necessarily) at a given perturbative order. Moreover, taking (3.7)drastically improves the convergence of the method: More precisely the known nonperturbative result of generic pureRG-resummed expressions are obtained exactly from RGOPT at the very first δ -order, while the convergence of theconventional OPT with a = 1 is not clear or very slow, if any (see Sec. III.C of ref. [22] for details).The criterion (3.6) can easily be generalized to any model, even nonasymptotically free ones, by similarly selectingthose optimized solutions that simply match the standard perturbative behavior for small coupling values. Thusclearly the resulting unique critical value like in (3.7) is valid for any model with its appropriate RG coefficients. For For QCD our normalization is β ( g ) ≡ dg/d ln µ = − b g − b g + · · · , γ m ( g ) = γ g + γ g + · · · , where g ≡ πα S . The b i , γ i coefficients up to four loops are given in [35]. The important role of the anomalous dimension γ / (2 b ) appeared also in our earlier constructions resumming the RG dependence ofthe δ -expansion [23, 24, 26, 29], although it had not been recognized at that time as a crucial RG consistency property of lowest δ orders. the QCD spectral density, as we will see below the equivalent of the criteria (3.6) indeed selects a unique solution ata given order for both the RG and OPT equations, at least up to the four-loop order (as was also the case for thepion decay constant [22]).Incidentally, a connection of the exponent a with RG anomalous dimensions/critical exponents had also beenestablished previously in a different context, in the D = 3 Φ model for the Bose-Einstein condensate (BEC) criticaltemperature shift, by two independent approaches [30, 31], where for this model it also leads to real OPT solutions [31].Indeed in ref. [30, 32] it was convincingly argued, based on critical behavior considerations, that the OPT can onlyconverge if an appropriate Wegner critical exponent is used in the interpolation (3.1), which appears quite similar toour criterion (3.7). Note however that (3.7) is identified exactly from the known first-order RG coefficients, thus validfor any model, while in [32] the analoguous exponent was determined more approximately by looking for a plateau inthe variational parameter dependence. In any case, from these examples, it is established that it is necessary for theOPT method to give useful results to have in general an exponent a in (3.1) differing from 1 in a well-defined way.Coming back to the present OPT and RG Eqs. (3.2), (3.5), beyond lowest orders AF-compatible solutions withbehavior (3.6) are however not necessarily real in general. A rather simple way out is to further exploit the RGfreedom, considering a perturbative renormalization scheme change to attempt to recover RGOPT solutions bothAF-compatible and real [22]. In the present case of the spectral density, this extra complication is not even necessary,at least up to the highest order here studied (four loops): We shall find also that the unique AF-compatible RG andOPT solution remains real. C. RGOPT for the spectral density
Formally a generic pertubative expansion for the condensate reads typically h ¯ qq i pert = m X p g p p X k =0 f pk ln p − k ( mµ ) , (3.8)where the f pk coefficients are determined by RG properties from lowest orders for k < p . According to Eq. (2.7),calculating the (perturbative) spectral density formally involves calculating all logarithmic discontinuities. This isconveniently given by taking in any typical perturbative expansion like (3.8), all nonlogarithmic terms to zero, thosetrivially having no discontinuities, while replacing all powers of logarithms, using m → i | λ | etc., asln n (cid:18) mµ (cid:19) = 12 n ln n (cid:18) m µ (cid:19) → n π (cid:20)(cid:18) | λ | µ + i π (cid:19) n − (cid:18) | λ | µ − i π (cid:19) n (cid:21) , (3.9)leading to the following simple substitution rules for the first few termsln (cid:18) mµ (cid:19) → /
2; ln (cid:18) mµ (cid:19) → ln | λ | µ ; ln (cid:18) mµ (cid:19) →
32 ln | λ | µ − π , (3.10)and so on (note the appearance of nonlogarithmic ∼ π terms starting at order ln m ). This gives a perturbativeexpression of the spectral density of the generic form ρ pert ( | λ | , g ) = | λ | X p ≥ g p p X k =0 f SDpk ln p − k ( | λ | µ ) , (3.11)where the determination of the coefficients f SDpk follows from the above relations (3.9).To obtain the RG equation for ρ ( g, λ ), we use the defining integral representation of the spectral density in (2.5) andthe basic algebraic identity ∂∂m mλ + m = − ∂∂λ λλ + m . (3.12)Throwing away surface terms in partial integrations as usual in the spirit of dimensional regularization, one thus findsthat ρ ( λ ) actually obeys the same RG equation as h ¯ qq i , with ∂m replaced by ∂λ , (cid:20) µ ∂∂µ + β ( g ) ∂∂g − γ m ( g ) λ ∂∂λ − γ m ( g ) (cid:21) ρ ( g, λ ) = 0 . (3.13)One can next proceed to the modification of the resulting perturbative series ρ ( λ, g ) as implied by the δ -expansion,now, from Eq. (3.12) clearly applied not on the original mass but on the spectral value λ : λ → λ (1 − δ ) a g → δ g . (3.14)Optimizing perturbation theory means that the derivative with respect to m of ∞ X n =0 ( − n n ! m n (cid:18) ∂∂m (cid:19) n h ¯ qq i , (3.15)being formally zero , one should obtain a good approximation for the value at m = 0 of h ¯ qq i at finite order by settingto zero the derivative of a finite number of terms of this series, see Eq. (3.2). Using Eq. (3.12), this mass optimizationon h ¯ qq i thus translates into an optimization of the spectral density with respect to λ , ∂ρ ( k ) ( λ, g ) ∂λ = 0 , (3.16)at successive δ k order. IV. LESSONS FROM THE GROSS-NEVEU MODEL
The fermion condensate can also be defined from the spectral density for the D = 1 + 1 O ( N ) Gross-Neveu (GN)model [36]. This will give us a very useful guidance for the more elaborate QCD case below, and we also set up someformulas actually generically valid for both the GN model and QCD.We start from the known expression of the vacuum energy evaluated in the large- N limit [26] after all necessarymass, coupling, and vacuum energy (additive) renormalizations, in terms of the explicit mass m and mass gap M ( m, g ), E GN = − (cid:18) N π (cid:19) (cid:18) M ( m, g ) + 2 mg M ( m, g ) (cid:19) , (4.1)where the mass gap is defined in compact form as M ( m, g ) = m (cid:18) g ln Mµ (cid:19) − , (4.2)where m ≡ m ( µ ) and g ≡ g ( µ ) are the renormalized mass and coupling in the M S scheme (after convenienly rescalingthe original coupling defined by (1 / g GN ( ¯ΨΨ) as g GN N/π ≡ g ). The fermion condensate is formally given by thederivative with respect to m of the vacuum energy, giving simply after some algebra h ¯ ψψ i GN ( m, g ) ≡ − (cid:18) N π (cid:19) M ( m, g ) g . (4.3)This means that up to a trivial overall factor, the fermion condensate is directly related to the mass gap, as intuitivelyexpected also from mean-field arguments. Eq. (4.3) has a well-defined nontrivial perturbative expansion to arbitraryorder, − (cid:18) πN (cid:19) h ¯ ψψ i GN ( m, g ) ≡ −h ¯ΨΨ i GN ≃ m (cid:18) g − L m + gL m ( L m + 1) − L m ( L m + 2)(2 L m + 1) g + O ( g ) (cid:19) , (4.4)where L m ≡ ln m/µ (and for convenience we redefined the condensate by a trivial rescaling). From properties of theimplicit M ( m ) defined by (4.2) and its reciprocal function m ( M ), one can establish [26] that M ( m ) → Λ ≡ µe − /g for m →
0, which translates here into the simple relation −h ¯ΨΨ i GN ( m →
0) = Λ g , (4.5) We adopt in the following the notation λ ≡ | λ | since it is necessarily positive. For simplicity, we have set a to one in this equation. which provides a consistent bridge between the massive and massless case. But deriving this needs the knowledge ofthe all-order expression (4.2), only known exactly in the large- N limit.Now alternatively, performing the substitution (3.1), expanding at order δ p , setting δ to one, and optimizing theresulting expression with Eqs. (3.2) and (3.5) gives the exact result (4.5) at any order in δ p , at optimized couplingand mass values, ˜ g = 1; ˜ L m = − , (4.6)just as in the mass gap case [20]. This is not very surprising, in view of the rather trivial relation (4.5) between themass gap and the condensate in the large- N limit.However here we shall try to obtain this large- N result in an indirect way, using the spectral density, the aim beingevidently to test on the exactly known result (4.5) the possibility of calculating a spectral density and of estimatingits reliability from the first few perturbative orders, in order to subsequently apply the same procedure in the morechallenging QCD case.Thus from (4.4) and using m → i | λ | and relations (3.9) the perturbative expression of the spectral density, forinstance restricted up to four-loop order g , follows as ρ pertGN ( λ, g ) ≃ λ (cid:18) −
12 + g (cid:18) L λ + 12 (cid:19) − g (cid:0) L λ + 20 L λ + 4 − π (cid:1) + g (cid:0) L λ + 156 L λ + (108 − π ) L λ + 12 − π (cid:1) + O ( g ) (cid:19) , (4.7)where now L λ ≡ ln λ/µ . We can then proceed by applying (3.14), expanded to order δ k , then taking δ →
1, and finallyapplying on the resulting expression the OPT (3.16) and RG (3.13) Eqs. In practice we shall proceed to relativelyhigh perturbative orders, since the exact expression (4.2) may be formally expanded to arbitrary order.Thus with those obvious replacements we can proceed first with the δ -expansion (3.14) operating on λ and g , andtaking the relevant value of the exponent for the large- N case, a ≡ γ / (2 b ) = 1 in (3.1). At first nontrivial δ orderwe simply obtain ρ δ ( δ → λ, g ) = λg (cid:18)
12 + L λ (cid:19) , (4.8)from which the OPT (3.16) and RG (3.13) give the unique solution,˜ L λ (cid:18) ≡ ln | λ | µ (cid:19) = −
32 ; ˜ g = 12 , (4.9)which plugged back into (4.8) gives the final result, using also the Banks-Casher relation (2.6),˜ g h ¯ΨΨ i GN ( m → / Λ ≡ − π ˜ gρ (˜ λ, ˜ g ) / Λ = − π √ e ≃ − . , (4.10)which is to be compared to the exact result g h ¯ΨΨ i GN ( m → / Λ = − M athematica r [37]).The results up to order g (beyond which numerics become really tedious) are given in Table I (where we give forconvenience the scale invariant condensate ˜ g h ¯ΨΨ i values). Actually, starting at order g , there are more than onereal solution. We show here those solutions which are unambiguously determined to be the closest to the standardAF perturbative behavior, L λ ≃ − /g + O (1) for g →
0, by analogy with the exact value − ˜ g ˜ L m ≡ g ≃ .
4, while the optimal spectral parameter ˜ L λ is less stable with wider variations at successive orders.These approximate results at successive orders when optimizing the spectral transform are clearly in contrastwith those obtained from optimizing directly the original perturbative expansion, where as explained above the AF Taking a = 1 as appropriate for the large- N GN model.
TABLE I. RGOPT − g h ¯ΨΨ i GN Λ and corresponding optimized coupling ˜ g and (logarithm of) optimized eigenvalue ln ˜ λ/µ resultsat successive orders in δ . δ k order − g h ¯ΨΨ i GN Λ ˜ g ln ˜ λµ . − .
984 0 . − . .
897 0 . − . .
081 0 . − . .
924 0 . − . .
877 0 . − . .
980 0 . − . .
903 0 . − . .
065 0 . − . .
934 0 . − . .
894 0 . − . .
978 0 . − . .
972 0 . − . .
013 0 . − . .
989 0 . − . .
027 0 . − . .
978 0 . − . .
013 0 . − . compatible solution always gives the exact result for the condensate, with the corresponding optimal coupling andmass values (4.6) obtained already at first and all successive orders. Here we see in Table I that the second orderis very close, less than 2%, from the exact result, but at higher orders the solutions exhibit a rather slow empiricalconvergence, with an oscillating behavior towards the exact large- N limit result.This slow convergence can be essentially traced to the effect of numerous factors π p , p = 1 , ...n/ g n , g n +1 for n ≥ N resummed mass in (4.2). Moreover, starting at order g these π p terms come with largercoefficients relative to the other non-logarithmic coefficients, originating from the ln m coefficients in the originalperturbation (4.2), (4.7) which are roughly all of the same order O (1). More precisely inspecting (4.7) one can seethat at order g the − π / / g the relevant terms to be compared are the two last ones: 12 /
24 and − (13 / π , sothe latter is an order of magnitude larger than the former. This explains the very good result at order g and also thedegraded result at the next g order. A similar behavior is observed at higher orders, until at sufficiently high orderthe original perturbative coefficients of ln m also start to grow fastly, such that a balance with the π p contributionscan again occur, with the (slow) convergence observed. Indeed the spectral parameter optimization (3.16) tends todamp these relative π p contributions: For instance the combined contributions of all L λ -dependent terms for theoptimal value ˜ L λ ∼ − .
28 at order g in Table I almost cancel the large π term. But since the latter terms aregrowing with the order, it is not surprising that the optimal ˜ L λ values are not much stable in Table I.Interestingly however, if instead of taking the exact results at a given order, we consider a well-defined approxima-tion, by keeping, at growing δ n order, only π p terms with a fixed maximal power, then there is a higher δ n order atwhich one recovers the simple exact RGOPT solution : ˜ g = 1 , ˜ L λ = − , h ¯ΨΨ i GN = − Λ. For instance, keeping only π and π terms (the latter appearing first at order g ) and increasing the δ k order, the exact solution is recovered atorder δ . This cancellation mechanism thus indicates that the “maximal convergence” properties of RGOPT, specificof the simpler GN model mass gap expression in (4.2) [20], are not completely lost within the perturbative spectraldensity, but rather hidden, being obstructed by the more involved perturbative coefficients. Those remarks are to bekept in mind when comparing with the QCD case below. This nice property may appear somewhat artificial but it can be explained more rigorously: As observed in [20] the mass gap M ( m, g )at a given (sufficiently high) order δ n , has flat optima roughly at order ∼ n/ i.e. its n/ m vanishes).Now for the spectral density, π p terms appearing at order g p arise from (2.7) as the coefficient of the (logarithmic) derivatives ∂/∂ ln m of order p : So if discarding terms of higher power π q , q > p , there is necessarily a fixed higher order at which the cancellation of all π p terms occur. V. DETERMINATION OF THE QCD QUARK CONDENSATEA. Perturbative three-loop quark condensate x x xx x
FIG. 1. Samples of standard perturbative QCD contributions to the chiral condensate up to three-loop order. The cross denotesa mass insertion.
The perturbative expansion of the QCD quark condensate for a nonzero quark mass can be calculated systematicallyfrom the related vacuum energy graphs. A few representative Feynman graphs contributions at successive orders upto three-loop order are illustrated in Fig. 1 (there are evidently a few more three-loop contributions not shown here).Note that the one-loop order is O (1) = O ( g ). The two-loop contributions were computed long ago [38] and thethree-loop ones in [40]. The three-loop order result in the M S -scheme reads explicitly m h ¯ qq i QCD ( m, g ) = 32 π m (cid:18) − L m + gπ ( L m − L m + 512 ) + ( g π ) q ( λ, n f ) (cid:19) , (5.1)where m ≡ m ( µ ) and g ≡ πα S ( µ ) are the running mass and coupling in the M S scheme, and the three-loop coefficientreads [40] q ( λ, n f ) = 127 (cid:0) − a −
32 ln z + 504 z + (672 z − n f + 528 z (cid:1) + (52 n f − z ) L m −
329 (5 n f − L m + 329 (2 n f − L m , (5.2)where z i ≡ ζ ( i ) and a = Li (1 /
2) .The calculation in dimensional regularization of (5.1) actually still contains divergent terms needing extra subtrac-tion after mass and coupling renormalizations in the
M S scheme. The correct procedure to obtain a RG invariantfinite expression when subtracting those divergences consistently is well known in the standard renormalization ofcomposite operators with mixing [39]. We can define [24] the needed subtraction as a perturbative series,sub( g, m ) ≡ m g X i ≥ s i g i , (5.3)with coefficients determined order by order by µ dd µ sub( g, m ) ≡ Remnant( g, m ) = µ dd µ [ m h ¯ qq i (pert) | finite ] , (5.4)where the remnant part is obtained by applying the RG operator Eq. (3.4) to the finite expression (5.1), not separatelyRG invariant. Eq. (5.3) does not contain any ln m/µ terms and necessarily starts with a s /g term to be consistentwith RG invariance properties. To obtain RG invariance at order g k , fixing s k in (5.3), one needs knowledge of thecoefficient of the ln m term (equivalently the coefficient of 1 /ǫ in dimensional regularization) at order g k +1 . Concerningthe condensate the extra contribution to the RG equation (5.4) is the so-called anomalous dimension of the QCD(quark) vacuum energy, entering the renormalization procedure of the m h ¯ qq i operator due to mixing with m × The originally calculated expressions in refs. [40] are given for arbitrary N c colors, n h massive quarks and n l massless quarks enteringat three-loop order. In our context m is the ( n f -degenerate) light quark mass and its precise mass dependence is what is relevant forthe optimization procedure, so one should trace properly the full n l , n h dependence, and take n l = 0 , n h ≡ n f with n f = 2 (3) for the SU (2) (resp. SU (3)) case. s i coefficients can be expressed in terms of RG coefficientsand other terms using RG properties. In compact form for completeness in our normalization they read , s = 14 π ( b − γ ) ,s = − π + 14 b − γ b − γ ,s = − (112077 + 24519 n f + 2101 n f + 576(15 − n f ) z )1152 π ( −
81 + 2 n f )(15 + 2 n f ) ,s ≃ . − a ( −
81 + 2 n f )( − .
615 + n f )(15 + 2 n f ) + 2 .
47 + 1 . n f + 0 . n f − . n f + 3 . − n f ( −
57 + 2 n f )( −
81 + 2 n f )(15 + 2 n f ) . (5.5) B. Perturbative spectral density
One crucial advantage of using the spectral density with the Banks-Casher relation (2.6) is that it gives a directaccess to the QCD condensate in the chiral limit, unlike the original direct RG invariant expression m h ¯ qq i in (5.1),where the condensate is being screened by the mass for m →
0. Indeed, note that a direct RGOPT optimisation inthe GN model of the corresponding expression m h ¯ΨΨ i gives an exactly vanishing result consistently at any orders [22](even though the optimized GN mass is clearly nonvanishing, ˜ m = Λ).Taking thus the logarithmic discontinuities according to Eqs. (2.7),(3.9) gives us the perturbative spectral densityup to three-loop order, − ρ MSQCD ( λ ) = 32 π λ (cid:18) −
12 + gπ (cid:18) L λ − (cid:19) + ( g π ) q SD ( λ, n f ) (cid:19) , (5.6)where now L λ ≡ ln( | λ | /µ ) and q SD ( λ, n f ) = 12 (52 n f − z ) −
329 (5 n f − L λ + 329 (2 n f − L λ − π . (5.7)Note that the π in the last term arises from the discontinuities of ln ( m ) according to (3.9). Now remark that noneof the nonlogarithmic contributions in the original perturbative expression (5.1) contribute to the spectral density.Thus similarly all subtraction terms in (5.5), which are necessary for RG invariance of the original expression, do notcontribute as well, thus making the final expression to be optimized relatively simpler.This point is worth elaborating in some detail. Just as for the GN model, instead of using the spectral density wecould apply the RGOPT method more directly to the original perturbative expression of the condensate, Eq. (5.1),including in this case the subtractions (5.3), (5.5) required by RG invariance (and also removing an overall factor m from Eq. (5.1) to define a nontrivial h ¯ qq i in the chiral limit). When this is done, one finds rather unstable values forthe optimized mass, coupling, and resulting condensate, showing no clear empirical convergence pattern at successiveorders, at least at the presently available (three-loop) order. Furthermore these results tend to give a wrong sign(positive, or ambiguous) condensate. More precisely, at the first nontrivial δ (one loop) order, firstly there is nocommon nontrivial RG and OPT Eqs. solution. Considering then the OPT or RG Eqs. alone, both give a positivecondensate, of roughly the right order of magnitude: h ¯ qq i ( δ ) = √ e/ (2 π )Λ MS ≃ . MS from the OPT, and a verysimilar value is obtained from the RG. Next at the δ (two-loop) order, the (unique) AF-compatible branch solutionof the combined RG and OPT Eqs. (3.5), (3.2) gives complex-valued optimized coupling, mass and condensate,with a negative real part for the condensate but a much larger imaginary part: h ¯ qq i ( δ ) ≃ ( − . ± . i )Λ MS ,a result clearly ambiguous. These calculations are also not very stable upon different truncations of perturbativehigher order terms in the RG equation (3.5). Furthermore, attempting to recover real AF-compatible solutions bya perturbative renormalization scheme change, both at orders δ and δ , happens to give no solutions (unlike for the The expression of s here approximated to 10 − relative uncertainty, largely sufficient for our purpose, needs the knowledge of thefour-loop ln m/µ coefficient, or alternatively the four-loop vacuum energy anomalous dimension. The latter, not explicitly available inthe published literature so far, has been kindly provided to us by K. Chetyrkin and P. Maier [41] from a related work. N approximation) is precisely given by the very same one-loop first graph of Fig. 1, up to trivial overallfactors, while genuine QCD contributions only enter at the next orders with gluon and further quark loop dressing.In dimensional regularization, at lowest orders in the coupling, one finds that the extra subtractions (5.3), (5.5) havea sign opposite to the sign of the similar terms in the GN model, and that this is due to the fact that the poleof Γ(1 − D/
2) in the perturbative calculation of the condensate changes sign when going from dimension D = 2(corresponding to the logarithmic divergence in the GN model, and quadratic divergence in the D = 4 NJL model)to dimension D = 4 (corresponding to QCD). Indeed, as is well-known most of the phenomenological successes ofthe NJL model rely strongly on the physical cutoff interpretation of the quadratically divergent mass gap in fourdimensions .Hence, the RGOPT appears bound to lead to a wrong sign and/or unstable results, if applied directly to the QCDperturbative expression of the quark condensate evaluated in dimensional regularization and related M S scheme atlow orders. We do not know whether higher orders would cure this problem by ultimately stabilizing the result,but this appears rather unlikely since the first few orders are likely to remain dominant in our approach. (Indeeda general property of the optimized perturbation is that the optimized coupling ˜ g turns to be reasonably small, sothat the first few orders dominate). We did not have this problem in our previous works [20–22] in which we weredealing with the pole mass and the pion decay constant, which are only logarithmically divergent quantities. Yet oneshould not hastily conclude that the OPT or RGOPT approaches are bound to fail in any situation where quadraticdivergergences would be present in a cutoff regularization . Rather, the above problems stress that in a given modelit is crucial to choose carefully the basic entity to be perturbatively modified and optimized within the RGOPTframework. (This is analoguous to the traditional variationnal Rayleigh-Ritz method in quantum mechanics, wherethe trial wave functions should often be appropriately chosen to obtain a sensible result). This is why for QCDone must use instead the spectral density, which in our framework we anyway derived from the very same originalperturbative condensate expression (5.1), but which at the same time formally gets rid of the influence of quadraticdivergences. Indeed, only the infrared part λ → e.g. relying on the GMOR relation.We thus proceed with the actual RGOPT calculations for the spectral density at successive perturbative orders.First remark that, since there is no logarithmic L λ contribution in the spectral density at one-loop order (the one-loopln m contribution in (3.10) only giving the constant 1 / λ = 0 optimized solution of (3.16) atone-loop order. Thus we should start applying our method at the next two-loop order. C. Two-loop O ( δ ) results Let us perform step by step the RGOPT optimization by restricting first (5.6) at the first nontrivial two-loop order.Concerning the δ -expansion given by (3.14), it is crucial [21, 22] to take the right value of the exponent a , determinedby the lowest order anomalous mass dimension, which makes the δ -modified series compatible with RG propertiesand matching asymptotic freedom (AF), as we recalled in some details in Sec. III.B above. In the case of the large- N limit of the GN model, one has simply a = γ / (2 b ) = 1. Actually, since m h ¯ qq i is RG invariant to all orders ratherthan h ¯ qq i , it is easily derived that the correct value to be used for h ¯ qq i , and thus for the related spectral density from(2.5), is a = 43 ( γ b ) . (5.8) The NJL model may be formulated in dimensional regularization to some extent, in dimensions 2 < D <
4, but with rather oddproperties, see e.g. [42] In fact the (standard) OPT has been applied in the framework of the effective NJL model with a cut-off at next-to-leading δ order [34],giving sensible results beyond the large- N approximation, and consistent with important basic properties like the Goldstone theorem,the GMOR relation, etc. δ the modified series reads, − ρ δ QCD = 32 π λ (cid:18) gπ ( L λ −
512 ) (cid:19) , (5.9)and the OPT (3.16) and RG (3.13) equations have a unique solution, using also (2.6), given in the first line of Table II.For a simpler first illustration we actually used the RG Eq. (3.13) at the very first order with the one-loop coefficient b , in order to get simple analytic solutions. Therefore we obtain h ¯ qq i / ( n f = 2)( µ ≃ . ) ≃ − .
96 ¯Λ , a fairlydecent value given this lowest nontrivial order. At two-loop perturbative order expression (5.9) does not dependexplicitly of the number of flavors n f , but a n f dependence enters into the optimized results indirectly from theRG Eq. (3.13) involving b ( n f ), also entering ¯Λ n f . The corresponding optimized coupling is ˜ α S ≡ ˜ g/ (4 π ) ≃ .
83, amoderately large value very similar to the optimal coupling values obtained, at first nontrivial RGOPT order, whenconsidering the pion decay constant F π in ref. [22]. Of course the precise number obtained for the condensate dependson the precise definition of the ¯Λ reference scale, which is generally perturbative and a matter of convention. To getthe numbers in the first lines of Table II we have used the simpler one-loop form, ¯Λ = µe − / (2 b g ) , consistently withthe one-loop RG Eq. used. When comparing below more precisely with other phenomenological determinations of thecondensate, we will use a more precise perturbative definition of ¯Λ, at four-loop order, in agreement with most otherpresent determinations. Remark also that the condensate being scale dependent, our RGOPT optimization fixes alsoa scale, consistently with a defining convention for ¯Λ, as indicated in Table II.For n f = 3 at order δ one finds similarly the optimized results given in the first line of Table III, indeed very close tothe n f = 2 results above.Now, since our basic expression originated from an exact two-loop calculation, it is a priori more sensible to applythe RG Eq. (3.13) at the same two-loop order, in order to capture as much as possible higher-order effects. Doingthis we obtain the results given in the second lines of Tables II, III for n f = 2 ,
3. Those results are therefore tobe considered more consistent at two-loop order. One can already observe the substantial decrease of the optimalcoupling α S to a more perturbative value, and the correspondingly higher optimal scale µ , with respect to the resultsusing pure one-loop RG equation.For completeness and later use we also give in Table II, III the corresponding values of the RG invariant condensate h ¯ qq i RGI , perturbatively defined in our normalization as h ¯ qq i RGI = h ¯ qq i ( µ ) (2 b g ) γ b (cid:18) γ b − γ b b ) g + O ( g ) (cid:19) , (5.10)where higher order terms not shown here are easily derived from integrating exp [ R dgγ m ( g ) /β ( g )] thus know pertur-batively to four-loop g order since only depending on the RG function coefficients [35] b i , γ i . The factor multiplyingthe scale-dependent condensate h ¯ qq i ( µ ) in (5.10) is obviously the inverse of the one defining similarly a scale-invariantmass, given explicitly to four-loop order in the literature (see e.g. [43]). We will calculate the RG invariant condensateat successive orders in Tables II, III using (5.10) consistently at the same perturbative order as the RG order usedin Eq. (3.13), and taking g ≡ ˜ g = 4 π ˜ α S , the corresponding optimized values obtained at each order. (Alternativelywe could optimize directly the expression (5.10) instead of h ¯ qq i ( µ ) (with the last term ∝ γ m ( g ) in the RG Eq. (3.13)removed): This gives the same optimized solutions, as expected since it is formally completely equivalent, up to tinynumerical differences due to perturbative reexpansions, of less than 10 − relative to the numbers given in Tables II,III.) TABLE II. Main optimized results at successive orders for n f = 2, for the optimized spectral parameter ˜ λ , the optimizedcoupling ˜ α S , and resulting optimized condensate. We also give the RG invariant condensate h ¯ qq i / RGI calculated at the consistentperturbative order from (5.10). ¯Λ is conventionally normalized everywhere by Eq. (5.11), except in the very first line wherethe one-loop expression ¯Λ ≡ µ e − / (2 b g ) is rather used. δ k , RG order ln ˜ λµ ˜ α S −h ¯ qq i / ¯Λ (˜ µ ) ˜ µ ¯Λ −h ¯ qq i / RGI ¯Λ δ , RG 1-loop − π ≃ .
83 0 .
962 2 . . δ , RG 2-loop − .
45 0 .
480 0 .
822 2 . . δ , RG 2-loop − .
686 0 .
483 0 .
792 2 .
797 0 . δ , RG 3-loop − .
703 0 .
430 0 .
794 3 .
104 0 . δ , RG 3-loop − .
838 0 .
405 0 .
793 3 .
306 0 . δ , RG 4-loop − .
820 0 .
391 0 .
796 3 .
446 0 . TABLE III. Same as Table II for n f = 3 δ k order ln ˜ λµ ˜ α S −h ¯ qq i / ¯Λ (˜ µ ) ˜ µ ¯Λ −h ¯ qq i / RGI ¯Λ δ , RG 1-loop − π ≃ .
82 0 .
965 2 .
35 0 . δ , RG 2-loop − .
56 0 .
474 0 .
799 3 .
06 0 . δ , RG 2-loop − .
766 0 .
493 0 .
776 2 .
942 0 . δ , RG 3-loop − .
788 0 .
444 0 .
780 3 .
273 0 . δ , RG 3-loop − .
967 0 .
414 0 .
769 3 .
540 0 . δ , RG 4-loop − .
958 0 .
400 0 .
773 3 .
700 0 . D. Three-loop O ( δ ) results At three-loop g order the n f dependence enters explicitly within the perturbative expression of the spectral density,see Fig 1 and the last g coefficient in Eq. (5.6), (5.7). This is interesting also in view of other results on the variationof the condensate value with the number of flavors [2, 18].We find a unique real AF-compatible optimized solution. More precisely, at this three-loop order there are actuallytwo real optimized solutions for ˜ L λ , ˜ α S , but the selection of the right physical solution is unambiguous since only oneis clearly compatible with AF behavior for g = 4 πα S →
0, ln(˜ λ/µ ) ≃ − d k / (2 b g ) + O (1) with d k = O (1), both forthe RG and OPT equations. In contrast the other real solution has for g → L Λ = ln ˜ λ/µ >
0, which also means incompatibility with perturbativity, since we expect µ ≫ ˜ λ just likethe perturbative range being µ ≫ ˜ m ∼ ¯Λ for the original expansion with mass dependence. Explicitly we obtain for n f = 2 , α S decreased by almost a factor two with respect to the lowest nontrivial orderresult above, which indicates that the resulting series is much more perturbative. But ˜ α S almost does not changeas compared to the more consistent two-loop results. Similarly ˜ α S further slightly decreases and the optimal scale ˜ µ increases when going from two to three loops in the RG Eq. (3.13).As a matter of numerical detail, to obtain the results in Tables II, III we took a convention of the QCD scale ¯Λbased on a perturbative four-loop expression [44]:¯Λ − loopn f ( g ) ≡ µ e − b g ( b g ) − b b exp (cid:20) − g b · (cid:18) ( b b − b b ) + ( b b − b b b + b b ) g (cid:19)(cid:21) . (5.11)This is convenient and important to compare precisely below with most recent other determinations, using the samefour-loop perturbative order conventions for ¯Λ. Actually when performing the RG Eq. (3.13) at order δ k it would bemore natural to adopt a ¯Λ convention at the consistent ( k + 1)-loop order ¯Λ k +1 , given from (5.11) by taking b = 0(and b = 0) respectively at three loops (two loops). But this only affects an overall normalization of the final result,as ¯Λ itself is not involved in the actual optimization process when performing Eq. (3.2) and (3.13). Besides, startingat three-loop order the differences obtained from such different conventions are minor. (The ¯Λ convention also affectsthe precise value of the optimal scales ˜ µ/ ¯Λ in Tables II, III from which we shall start evolution to higher scale tocompare with other determinations in the literature, see below). Strictly speaking, the different values of h ¯ qq i ( µ )obtained in Tables II, III e.g. at three-loop order cannot be directly compared, being obtained at different scales ˜ µ .Thus we also give the scale-invariant condensate values h ¯ qq i RGI which can be more appropriately compared.Notice also that in spite of the more than 10% change in the optimal coupling ˜ α S when taking two- or three-loopRG, the final physical value of the condensate only varied by 0 . : This also reflects a strong stability. Moreover Due to the nontrivial relations between ln m and ln λ , Eqs. (3.9), (3.10), the 1 /g coefficient of the correct AF branch L λ ( g ) for g → − / (2 b ), like it is for ln m/µ (Eq. (3.6)), essentially determined by the leading logarithms g k ln k ( m/µ ) behavior. Butthis 1 /g coefficient has the correct AF sign, being at order δ k : − d k / (2 b ) with d k > d k → k increases. This result appears so stable partly due to the 1 / h ¯ qq i , the corresponding variation is ≃ . h ¯ qq i / / ¯Λ changed by about 20% with respect to the crude two-loop order result (first lines of Tables),but much less when compared to the more consistent two-loop order result (second lines). This shows a posteriori thatstability appears at the first nontrivial two-loop result, with an already quite realistic value. This stability, suggestingthat we remain within the domain of validity of perturbation theory, is an important requirement for the usefulness ofour method. A similar behavior was observed when optimizing the pion decay constant in [22]. Remark also that theoptimized coupling values ˜ α S at successive orders happen to be rather close to those for which the scale invariancefactor in (5.10) multiplying h ¯ qq i ( µ ) would be exactly one (which for n f = 2 happens for α S ≃ . , . , . α S is close to a (variational) “fixed-point” scale-invariant behavior. Had we found optimized ˜ α S values a factor 2-3 smaller, or larger, we would obtain no valuableresults beyond ordinary perturbation in the first case, or much more unstable results in the second case. However westress that the optimal coupling ˜ α S or optimal mass ˜ m do not have really a universal physical interpretation sincethe precise ˜ α S and ˜ m values depend on the physical quantity being optimized. For instance when optimizing F π in ref. [22] at a given perturbative order, the corresponding ˜ α S values were pretty much close to the present ones,but nevertheless slightly different. The physically meaningful result is obtained when replacing ˜ α S and ˜ λ within thequantity being optimized, like here the condensate.Comparing Tables II, III it also clearly appears that the ratio of the quark condensate to ¯Λ has a moderatedependence on the number of flavors n f , although there is a definite trend that h ¯ qq i / n f =3 is smaller by about 2 − h ¯ qq i / n f =2 , in units of ¯Λ n f , at the same perturbative orders. The smallness of this difference wasquite expected, due to the n f dependence only appearing at three-loop order and the overall stability of the modifiedperturbation. However from various different estimations, including lattice [46], and ours [22]), there are someindications that ¯Λ > ¯Λ (although unclear from uncertainties, due to a larger uncertainty on ¯Λ ), which thereforecould indirectly further affect the actual flavor dependence of the condensate. We shall come back in more detailbelow on this point in the phenomenological discussion in next section, after establishing our final result for the precisecondensate values. E. Four-loop O ( δ ) results We finally consider the optimization of the spectral density at four-loop order, the maximal order available atpresent. In fact, the complete standard perturbative expression of our starting expression for the condensate, i.e. thenext α S order correction to (5.1), is not fully known at present. But it obviously takes the form m h ¯ qq i − loopQCD ( m, g ) = 32 π m ( g π ) (cid:0) c L m + c L m + c L m + c L m + c (cid:1) (5.12)where L m ≡ ln( m/µ ) and we choose a convenient overall normalization with respect to lowest order terms in (5.1).Now the leading (LL), next-to-leading (NLL) and next-to-next-to-leading (NNLL) logarithms coefficients c − c areeasy to derive from RG invariance properties as being fully determined by lowest orders. The NNNLL ln m coefficient c can also be inferred by RG properties from the available anomalous dimension of the vacuum energy, calculatedby Chetyrkin and Maier [41], and related to s given in Eq. (5.5). Explicitly, we obtain c ≃ .
74 (4533 . c ≃ − . − . c ≃ . . c ≃ − . − .
5) (5.13)where the first and second numbers correspond to n f = 2 and n f = 3 respectively. (N.B.: We can obtain the genericalgebraic values of c i ( n f ) but these are rather involved and not particularly instructive, so we prefer to keep anapproximate numerical form for the relevant n f = 2 , c isactually unknown at present, and could be quite challenging to compute. But since the nonlogarithmic parts cannotcontribute to the spectral density, the latter can thus be fully determined at four loops! This gives for the exact perturbative four-loop contribution to the spectral density, after taking the logarithmic singularities according toEqs. (2.7),(3.9): − ρ MS λ ) = 32 π λ ( g π ) (cid:18) c ( n f )(2 L λ − π L λ ) + c ( n f )( 32 L λ − π c ( n f ) L λ + 12 c ( n f ) (cid:19) , (5.14)6to be added to the three-loop expression in (5.6). It allows us to calculate the spectral density and the relatedcondensate at three successive orders of the variationally modified perturbation, which gives further confidence andan important stability and convergence check of our result.We obtain at four-loop order once more a unique real common RG and OPT AF compatible solution. (The bruteoptimization results actually give several real solutions for ˜ λ, ˜ α S but there are no possible ambiguities since all solutionsare eliminated from the AF compatibility requirement, except a single one, with ˜ α S > L λ < n f = 2 , α S to more perturbative values, with respect to the three-loopresults above, as well as the corresponding decrease of ˜ L λ , meaning that ˜ µ is also larger. The stabilization/convergenceof the results is even clearer for the scale-invariant condensate h ¯ qq i RGI given in the last columns in Tables II, III,which at four-loop order has almost no variation upon RG Eq. truncations.To better appreciate the very good stability of these results, consider the basic perturbative expression of thecondensate (5.6) up to four loops in more numerical form (for n f = 2) and a more standard normalization of thecoupling: − ρ MSQCD (4 − loop )( λ ) ≃ π λ (cid:18) −
12 + 4 α S π ( L λ − .
42) + ( α S π ) (9 .
46 + (29 . − . L λ ) L λ )+( α S π ) (91 . L λ ( −
129 + L λ ( −
288 + 151 L λ )) (cid:17) . (5.15)From this one can easily appreciate that the successive perturbative terms are not small, just like in most perturbativeQCD series: At successive orders the coefficients grow rapidly (even if partly damped by the decreasing α S /π higherpowers, provided that α S remains rather moderate). In fact for the relevant above values of ˜ α S ≃ . − . L λ ≃ − ( . − . α S has perturbative coefficients that similarly grow atsuccessive orders. But the RGOPT mass and coupling optimization manage to stabilize the series, in such a waythat the discrepancies between the three- and four-loop orders in the final h ¯ qq i / results are about 2% or less. Itis important that the optimized sequence has thus clearly further stabilized from three- to four-loop order, to bemore confident in a quite precise determination of the condensate, although the variation from the lowest nontrivialtwo-loop to three-loop results were already very reasonable, by only ∼ N GN model in Table I. This is a bit surprising a priori , given that direct optimizations,not going through the spectral density, give maximal convergence for the large- N GN model [20]. In fact one canunderstand these results as follows. As explained in section IV above, the rather slow convergence for the GN spectraldensity is entirely due to the large and growing factors of π p coming from the discontinuities (3.9) at successiveorders, spoiling the simple form of the series and ‘screening’ the otherwise maximal convergence with the neat so-lution (4.6). Now although the π p coefficients from (3.9) are universal, thus the same for QCD, once combinedwith the original perturbative coefficients of (5.1) their relative contributions with respect to the other perturbativeterms remain more reasonably of the same order in the QCD case than in the GN case. This is because the originalperturbative coefficients are comparatively larger in the QCD case. More precisely inspecting the QCD spectraldensity series at three-loop order g in (5.7), the π contribution (last term in the RHS of (5.7)) is roughly twicethe other non-logarithmic contribution (first term in the RHS of (5.7)). Similarly at the next four-loop g orderfrom (5.12) the π contributions are roughly twice the other non-logarithmic contribution ( c / π p contributions remain roughly of the same order of magnitude as the original perturbative coefficients, such thatsome balance can occur from the optimization process, they should produce a moderate disturbance of the observedstability. We expect these properties to remain true at even higher orders, because the QCD original perturbationcoefficients also grow rapidly with the order.Lines 4 to 6 of Tables II, III give our direct optimization results for n f = 2 , VI. EVOLUTION TO HIGHER ENERGY AND PHENOMENOLOGICAL COMPARISON
In order to get a more precise result it is necessary to take into account the (moderate but not completely negligible)running of the condensate values, since the optimal scales obtained are somewhat different at three- and four-looporder, though reasonably perturbative, roughly of order ˜ µ > ∼ µ ≃ A. RGOPT h ¯ qq i ( µ = 2GeV) results for n f = 2 and n f = 3 The procedure to evolve perturbatively the condensate from one to another scale is straightforward since from exactRG invariance of m h ¯ qq i it is simply given by the inverse of the well-known running of the quark masses, h ¯ qq i ( µ ′ ) = h ¯ qq i ( µ ) exp[ Z g ( µ ′ ) g ( µ ) dg γ m ( g ) β ( g ) ] . (6.1)Alternatively we may take the values of the scale-invariant condensate (5.10) as obtained in the last columns ofTab. II, III and extract from those the condensate at any chosen (perturbative) scale µ ′ by using again (5.10) nowtaking g ≡ πα S ( µ ′ ). This is of course fully equivalent to performing the running from µ to µ ′ with Eq. (6.1). Sinceall relevant scales ˜ µ obtained above are in a fairly perturbative range > ∼ . We choose the highest optimized scales obtained (given by the four-loop results for both n f = 2 and n f = 3 cases) as reference low scale(s): µ ref (2) = 3 .
45 ¯Λ and µ ref (3) = 3 .
70 ¯Λ respectively. (N.B.: Given the presentvalues of ¯Λ ≃ ± µ ref (3) happens quite accidentally to be very closely below the charm quark massthreshold). For example this running gives a ≃
2% increase of the three-loop h ¯ qq i / n f =2 value given in Tab. II andquite similarly for n f = 3. Putting all this together we obtain −h ¯ qq i / n f =2 ( µ ref (2) = 3 .
45 ¯Λ ) = (0 . − . , −h ¯ qq i / n f =3 ( µ ref (3) = 3 .
70 ¯Λ ) = (0 . − . , (6.2)which are our intermediate results in terms of ¯Λ and at the respective n f = 2 , µ ref ) and direct four-loop results. We give results in the form (6.2) in view of possibly moreprecise determinations of ¯Λ , in the future. Note that for both n f = 2 , µ ref values directly from optimisation,so without possible extra uncertainties from running.Next we perform a final evolution from the low reference scales µ ref (2 ,
3) relevant for n f = 2 and n f = 3 respectivelyas given in (6.2), up to the conventional scale µ ′ = 2 GeV, where from the present world average ¯Λ above we find¯ α S (2GeV) ≃ . ± . n f = 3 we take into account the four-loop expression of the perturbative matching[45] at the crossing of the charm quark threshold. Overall this leads to an increase of the values in (6.2) for |h ¯ qq i| / ofabout ∼ .
6% for n f = 2 and 5.3% for n f = 3 (in which taking into accout the charm quark threshold with matchingrelations for α S ( µ ≃ m c ) contributes to ∼ − .
3% with respect to a more naive one-step evolution ignoring charmthreshold effects). More precisely: −h ¯ qq i / n f =2 (2GeV) = (0 . − . −h ¯ qq i / n f =3 (2GeV) = (0 . − . . (6.3)To give a more precise determination for n f = 2 one obstacle is the presently not very precisely known value of thebasic scale ¯Λ . In principle it is beyond purely perturbative reach, as it cannot be ‘perturbatively connected’ to themore precisely known ¯Λ value [44]. Our own estimate [22] of ¯Λ from the pion decay constant gave ¯Λ ≃ +42 − MeV,while some recent lattice determinations are ¯Λ ≃ ±
45 (staggered Wilson fermions [46]) and ¯Λ ≃ ±
21 (quark For n f = 2, it is implicitly understood that this evolution is performed in a simplified QCD world where the strange and heavier quarksare all infinitely massive, i.e. ‘integrated out’. Otherwise it would not make sense perturbatively to take into account the strange quarkmass threshold effects on the running. For n f = 3 we can perform a more realistic running, taking into account properly the charmquark mass threshold effects on α S ( µ ∼ m c ), see below. by ∼
15 GeV with respect to former similar determinations [48]). Since lattice uncertaintiesare mostly statistical and systematic, while ours are theoretical errors, it is not obvious how to combine all these in asensible manner. We thus prefer to keep separate estimates. For a representative illustration, combining our presentresults in (6.3) with the above quoted most precise lattice values of ¯Λ , we obtain −h ¯ qq i / n f =2 (2GeV , lattice ¯Λ ) ≃ ± ±
18 MeV; (6.4)where the first quoted error is our intrinsical theoretical error propagated from the one in (6.2), while the second largeruncertainty originates from the lattice ones on ¯Λ . Using instead solely our above quoted RGOPT determination [22]of ¯Λ gives somewhat higher values with larger uncertainties: −h ¯ qq i / n f =2 (2GeV , rgopt ¯Λ ) ≃ ± +35 − MeV . (6.5)For n f = 3 the more precisely known ¯Λ value from many different determinations allows a more precise determinationof the condensate. Taking the latest world average values [44] ¯ α S ( m Z ) = 0 . ± . wa = 340 ± −h ¯ qq i / n f =3 (2GeV , ¯Λ wa3 ) ≃ ± ± .Using alternatively solely our RGOPT determination [22] of ¯Λ = 317 +27 − MeV, gives instead slightly lower values butwith larger uncertainties: −h ¯ qq i / n f =3 (2GeV , rgopt ¯Λ ) ≃ ± +22 − MeV . (6.7)Finally rather than fixing the scale from ¯Λ, it may be more sensible to give our results for the ratio of the scale-invariant condensate with another physical scale, which is a parameter-free prediction. Taking the h ¯ qq i / RGI results inTab. II, III, and using solely our previous RGOPT results [22] for F/ ¯Λ and F / ¯Λ (where F ( F ) are the pion decayconstant for n f = 2, n f = 3 respectively in the chiral limit), we obtain: −h ¯ qq i / RGI,n f =2 F = 3 . ± . +0 . − . , (6.8)where the first error comes from the present calculation of the condensate while the second ones from taking the mostconservative range combining linearly three- and four-loop order uncertainty results for F/ ¯Λ from Eq. (4.28) of [22].A less conservative estimate may be obtained alternatively by taking the range spanned by the maximal availablefour-loop results for F/ ¯Λ correlated with the four-loop condensate results. This gives −h ¯ qq i / RGI,n f =2 F = 3 . ± . +0 . − . . (6.9)Similarly for n f = 3 we obtain: −h ¯ qq i / RGI,n f =3 F = 3 . ± . +0 . − . , (6.10)where the first theoretical error comes from the condensate with the second one from the conservative range combininglinearly three- and four-loop order uncertainty on F / ¯Λ from Eq. (4.30) of [22]. As observed above, the direct resultsfrom optimization of −h ¯ qq i ( n f ) / ¯Λ n f in Tables II, III show a moderate relative decrease, of about 2 −
3% only on h ¯ qq i / n f =3 / ¯Λ . The effect appears slightly more pronounced, about 7% relative reduction from n f = 2 to n f = 3 whencomparing the central values of (6.8) and (6.10), due to the slight 4% reduction of (central) F relative to F , althoughthis result is not clear, being affected with rather large uncertainties. We may finally combine (6.8) and (6.10) to give h ¯ qq i / RGI,n f =3 h ¯ qq i / RGI,n f =2 ≃ (0 . ± .
01) ¯Λ ¯Λ ≃ (0 . ± . ± . F F ≃ . ± . ± . ± . , (6.11)9where all our theoretical errors are combined linearly. In the results on the RHS of (6.11) the first quoted errors arethe intrinsic RGOPT errors for the present condensate calculation only, and the second larger one is propagated fromthe F/ ¯Λ and F / ¯Λ RGOPT theoretical errors. We also stress in (6.11) that our results are by construction in thestrict chiral limit m q →
0. The result given for unspecified F /F corresponds to the present sole RGOPT condensateestimate without extra input from other methods, while the last result uses the present lattice F /F estimates [2](with its own uncertainty ∼ .
05 quoted last).
B. Comparison and discussion
One may compare (6.4), (6.5) with the latest, presently most precise lattice determination, from the spectraldensity [13] for n f = 2: h ¯ qq i / n f =2 ( µ = 2GeV) = − (261 ± ± in (6.5) tending to be relatively high. Note however that the above quoted lattice valuefrom [13] was obtained by fixing the scale with the F K decay constant rather than using ¯Λ (moreover with F K beingdetermined in the quenched approximation). It is thus probably more judicious to compare our results for the RGinvariant ratio (6.8) with theirs [13]: 2 . ± . ± .
04. Overall, recent lattice determinations from various othermethods lie in the range roughly h ¯ qq i / n f =2 ∼ − (220 − n f = 2 [2], and quite similarly for n f = 3. Themost precise lattice n f = 3 determination we are aware of is h ¯ qq i / n f =3 (2GeV) = − (245 ±
16) MeV [49]. Concerningthe n f = 3 to n f = 2 condensate ratio, various lattice results have still rather large uncertainties at present [2] butsome recent results are compatible with a ratio unity [50]. Our results compare a bit better with the latest ones fromspectral sum rules [3]: h ¯ uu i / ∼ − (276 ±
7) MeV (and for the ratio [51, 52]: h ¯ ss i / h ¯ uu i = 0 . +0 . − . ). However thesum rules results [3] actually determine precisely the current quark masses, extracting then the h ¯ uu i value indirectlyfrom using the exact GMOR relation (1.1).We stress again the rather moderate n f dependence of our result. This is in some tension with the larger estimateddifference between n f = 2 and n f = 3 cases obtained by some authors [18]. Since our results are by construction validin the strict chiral limit, taken at face value they indicate that the possibly larger difference obtained by some otherdeterminations is more likely due to the explicit breaking from the large strange quark mass, rather than an intrinsic n f dependence property of the condensate in the exact chiral limit. VII. SUMMARY AND CONCLUSION
We have adapted and applied our RGOPT method on the perturbative expression of the spectral density of the Diracoperator, the latter being first obtained from the perturbative logarithmic discontinuities of the quark condensate inthe
M S -scheme. This construction allows successive sequences of optimized nontrivial results in the strict chirallimit at two-, three- and four-loop levels. These results exhibit a remarkable stability and empirical convergence.The intrinsic theoretical error of the method, taken as the difference between three- and four-loop results, is of order2%, while the final condensate value uncertainty is more affected by the present uncertainties on the basic QCDscale ¯Λ, specially with a larger uncertainty for n f = 2 flavors. The values obtained are rather compatible, withinuncertainties, with the most recent lattice and sum rules determinations for n f = 2, and our values indicate a moderateflavor dependence of the h ¯ qq i / n f / ¯Λ n f ratio. [1] M. Gell-Mann, R. J. Oakes and B. Renner, Phys. Rev. , 2195 (1968).[2] G. Colangelo et al. , Eur. Phys. J. C , 1695 (2011); S. Aoki et al. , arXiv:1310.8555 [hep-lat].[3] S. Narison, Phys. Lett. B , 346 (2014) [arXiv:1401.3689 [hep-ph]].[4] Y. Nambu and G. Jona-Lasinio, Phys. Rev. , 345 (1961).[5] S. P. Klevansky, Rev. Mod. Phys. , 649 (1992); T. Hatsuda and T. Kunihiro, Phys. Rept. , 221 (1994).[6] P. Maris, C. D. Roberts, and P. C. Tandy, Phys. Lett. B , 267 (1998); P. Maris and C. D. Roberts, Phys. Rev. C ,3369 (1997).[7] K. Langfeld, H. Markum, R. Pullirsch, C. D. Roberts, and S. M. Schmidt, Phys. Rev. C , 065206 (2003).[8] See e.g. H. G. Dosch and S. Narison, Phys. Lett. B , 173 (1998); M. Jamin, Phys. Lett. B , 71 (2002). [9] M. A. Shifman, A. I. Vainshtein, and V. I. Zakharov, Nucl. Phys. B147 , 385 (1979).[10] See e.g.
D. Becirevic and V. Lubicz, Phys. Lett. B , 83 (2004); P. Hernandez, K. Jansen, L. Lellouch, and H. Wittig,J. High Energy Phys. , 018 (2001); A. Duncan, S. Pernice, and J. Yoo, Phys. Rev. D , 094509 (2002).[11] E. Marinari, G. Parisi, and C. Rebbi, Phys. Rev. Lett. , 1795 (1981).[12] L. Del Debbio, L. Giusti, M. L¨uscher, R. Petronzio, and N. Tantalo, J. High Energy Phys. , 011 (2006); L. Giusti,and M. L¨uscher, J. High Energy Phys. 03 ( ) 013; H. Fukaya, S. Aoki, T.W. Chiu, S. Hashimoto, T. Kaneko, J. Noaki,T. Onogi, and N. Yamada, Phys. Rev. Lett. , 122002 (2010) [Erratum- ibid. , 159901 (2010)].[13] G. P. Engel, L. Giusti, S. Lottini, and R. Sommer, arXiv:1411.6386 [hep-lat].[14] T. Banks and A. Casher, Nucl. Phys. B169 , 103 (1980);[15] J. Gasser and H. Leutwyler, Annals Phys. , 142 (1984); Nucl. Phys.
B250 , 465 (1985).[16] H. Leutwyler and A. V. Smilga, Phys. Rev. D (1992) 5607.[17] A. V. Smilga and J. Stern, Phys. Lett. B , 531 (1993); K. Zyablyuk, J. High Energy Phys. , 025 (2000).[18] S. Descotes-Genon, L. Girlanda, and J. Stern, J. High Energy Phys. 01 , 274 (2000); S. Descotes-Genon, N. H. Fuchs, L. Girlanda, and J. Stern, Eur. Phys. J. C , 201 (2004);V. Bernard, S. Descotes-Genon, and G. Toucas, J. High Energy Phys. 01 ( ) 107; V. Bernard, S. Descotes-Genon,and G. Toucas, J. High Energy Phys. , 051 (2012); M. Kolesar and J. Novotny, Nucl. Part. Phys. Proc. , 90[arXiv:1409.3380 [hep-ph]].[19] J. Bijnens and T. A. Lahde, Phys. Rev. D (2005) 094502.[20] J.-L. Kneur and A. Neveu, Phys. Rev. D , 125012 (2010).[21] J.-L. Kneur and A. Neveu, Phys. Rev. D , 014005 (2012).[22] J.-L. Kneur and A. Neveu, Phys. Rev. D , 074025 (2013).[23] C. Arvanitis, F. Geniet, J. L. Kneur, and A. Neveu, Phys. Lett. B , 385 (1997).[24] J.-L. Kneur, Phys. Rev. D , 2785 (1998).[25] V.I. Yukalov, Theor. Math. Phys. , 652 (1976); W.E. Caswell, Ann. Phys. (N.Y) , 153 (1979); I.G. Halliday and P.Suranyi, Phys. Lett. B , 421 (1979); P. M. Stevenson, Phys. Rev. D , 2916 (1981); Nucl. Phys. B203 , 472 (1982); J.Killinbeck, J. Phys. A , 1005 (1981); R.P. Feynman and H. Kleinert, Phys. Rev. A , 5080 (1986); A. Okopinska, Phys.Rev. D , 1835 (1987); A. Duncan and M. Moshe, Phys. Lett. B , 352 (1988); H.F. Jones and M. Moshe, Phys. Lett.B , 492 (1990); A. Neveu, Nucl. Phys. B, Proc. Suppl. B18 , 242 (1991); V. Yukalov, J. Math. Phys (N.Y.) , 1235(1991); S. Gandhi, H.F. Jones and M. Pinto, Nucl. Phys. B359 , 429 (1991); C. M. Bender et al. , Phys. Rev. D , 1248(1992); S. Gandhi and M.B. Pinto, Phys. Rev. D , 2570 (1992); H. Yamada, Z. Phys. C , 67 (1993); K.G. Klimenko,Z. Phys. C , 677 (1993); A.N. Sissakian, I.L. Solovtsov, and O.P. Solovtsova, Phys. Lett. B , 381 (1994); H. Kleinert,Phys. Rev. D , 2264 (1998); Phys. Lett. B , 74 (1998); Phys. Rev. D , 085001 (1999); Mod. Phys. Lett. B ,1011 (2003).[26] C. Arvanitis, F. Geniet, M. Iacomi, J.-L. Kneur, and A. Neveu, Int. J. Mod. Phys. A , 3307 (1997).[27] R. Seznec and J. Zinn-Justin, J. Math. Phys. , 1398 (1979); J.C. Le Guillou and J. Zinn-Justin, Ann. Phys. , 57(1983); J. Zinn-Justin, arXiv:1001.0675.[28] R. Guida, K. Konishi, and H. Suzuki, Ann. Phys. (N.Y.) , 152 (1995); , 109 (1996).[29] J.-L. Kneur and D. Reynaud, Phys. Rev. D , 085020 (2002).[30] H. Kleinert, Mod. Phys. Lett. B , 1011 (2003); B. Kastening, Phys. Rev. A , 061601 (2003); Phys.Rev. A , 043613(2004).[31] J.-L. Kneur, A. Neveu, and M. B. Pinto, Phys. Rev. A , 053624 (2004).[32] B. Hamprecht and H. Kleinert, Phys. Rev. D , 065001 (2003).[33] J. L. Kneur, M. B. Pinto, R. O. Ramos, and E. Staudt, Phys. Rev. D , 045020 (2007).[34] J. L. Kneur, M. B. Pinto and R. O. Ramos, Phys. Rev. C , 065205 (2010)[35] J. A. M. Vermaseren, S. A. Larin, and T. van Ritbergen, Phys. Lett. B , 327 (1997). K. G. Chetyrkin, Nucl. Phys. B710 , 499 (2005); M. Czakon Nucl. Phys.
B710 , 485 (2005).[36] D. J. Gross and A. Neveu, Phys. Rev. D , 3235 (1974).[37] Mathematica, version 6, S. Wolfram Company.[38] K. G. Chetyrkin and V. P. Spiridonov, Sov. J. Nucl. Phys. , 522 (1988);[39] J. C. Collins, “Renormalization. An Introduction To Renormalization, The Renormalization Group, And The OperatorProduct Expansion” , (1984) Cambridge University Press, Cambridge, UK.[40] K. G. Chetyrkin and J. H. K¨uhn, Nucl. Phys. B432 , 337 (1994); K. G. Chetyrkin and A. Maier, J. High Energy Phys. 01( ) 092.[41] K. G. Chetyrkin and A. Maier, private communication.[42] T. Inagaki, D. Kimura and A. Kvinikhidze, Phys. Rev. D , 116004 (2008).[43] K. G. Chetyrkin, J. H. Kuhn and M. Steinhauser, Comput. Phys. Commun. , 43 (2000)[44] See the QCD review chapter by S. Bethke, G. Dissertori and G. P. Salam in K. A. Olive et al. [Particle Data GroupCollaboration], Chin. Phys. C , 090001 (2014).[45] K.G. Chetyrkin, B.A. Kniehl and M. Steinhauser, Phys. Rev. Lett. (1997) 2184 [hep-ph/9706430]; Nucl. Phys. B510 ,61 (1998); K. G. Chetyrkin, J.H K¨uhn and C. Sturm, Nucl. Phys.
B744 , 121 (2006).[46] ETM collaboration, B. Blossier et al. , Phys. Rev. D , 034510 (2010).[47] F. Karbstein, A. Peters and M. Wagner, arXiv:1407.7503 [hep-ph].[48] K. Jansen et al. , [ETM Collaboration], J. High Energy Phys. 01 ( ) 025.[49] A. Bazavov et al. , [arXiv:0910.2966]. [50] C. T. H. Davies, C. McNeile, A. Bazavov, R. J. Dowdall, K. Hornbostel, G. P. Lepage and H. Trottier, PoS ConfinementX,042 (2012).[51] C. A. Dominguez, N. F. Nasrallah, R. Rontsch and K. Schilcher, J. High Energy Phys. , 020 (2008).[52] S. Narison and R. Albuquerque, Phys. Lett. B , 217 (2010); R. M. Albuquerque, S. Narison, and M. Nielsen, Phys.Lett. B684