$\Lambda_{\bar{\textrm{MS}}}^{(n_f=2)}$ from a momentum space analysis of the quark-antiquark static potential

We determine $\Lambda_{\bar{\textrm{MS}}}^{(n_f=2)}$ by fitting perturbative expressions for the quark-antiquark static potential to lattice results for QCD with $n_f=2$ dynamical quark flavors. To this end we use the perturbative static potential at the presently best known accuracy, i.e. up to ${\cal O}(\alpha_s^4)$, in momentum space. The lattice potential is computed on a fine lattice with $a \approx 0.042 \, \textrm{fm}$ in position space. To allow for a comparison and matching of both results, the lattice potential is transformed into momentum space by means of a discrete Fourier transform. The value of $\Lambda_{\bar{\textrm{MS}}}^{(n_f=2)}$ is extracted in momentum space. All sources of statistical and systematic errors are discussed. The uncertainty in the value of $\Lambda_{\bar{\textrm{MS}}}^{(n_f=2)}$ is found to be smaller than that obtained in a recent position space analysis of the static potential based on the same lattice data.


Introduction
In this paper we aim at determining Λ MS by comparing lattice and perturbative results for the quark-antiquark (QQ) static potential 1 in momentum space. More precisely, we restrict ourselves to quantum chromodynamics (QCD) with n f = 2 dynamical quark flavors, i.e. exclusively focus on Λ The QQ static potential V amounts to the interaction energy of the color-singlet state made up of a static quark Q and its antiquarkQ separated by a distance r = | r| (in position space), or characterized by a momentum transfer p = | p| (in momentum space), respectively. Accordingly, we denote the potential in position space by V (r) and in momentum space by V (p).
While many studies aiming at the extraction of Λ MS from the QQ static potential have been performed so far, e.g. [1][2][3][4][5][6][7][8], to the best of our knowledge we are the first to transform the lattice data to momentum space and compare and match perturbative and lattice results in momentum space. As will be discussed in detail below, our present study is mainly motivated by the significantly worse convergence behavior of the perturbative potential in position space as compared to momentum space [9][10][11][12].
Higher-order perturbative calculations in QCD are most conveniently performed in momentum space. This is particularly true for the perturbative QQ static potential. Hence, the perturbative expression for the static potential at the highest current accuracy, which encompasses all contributions up to O(α 4 s ), is directly accessible in momentum space. Due to the fact that QCD is asymptotically free, perturbative calculations in QCD are viable only at large momentum transfers p ≫ Λ QCD , with Λ QCD denoting the QCD (momentum) scale, which can be seen as separating the regimes of perturbative and non-perturbative physics. However, note that in standard perturbation theory loop diagrams come along with integrations of the loop four-momentum over the full momentum regime, implying that such loops naturally also receive contributions from momenta Λ QCD for which perturbation theory is no longer trustworthy. The leading uncontrolled contribution δV (p) contained in the perturbative potential in momentum space arising from this kind of diagrams is quadratic in Λ QCD and scales as δV (p) ∼ − 4παs [12].
Taking into account that V (p) ∼ − 4παs p 2 1 + O(α s ) and p ≫ Λ QCD , i.e. Λ QCD p ≪ 1, this only amounts to a tiny correction. For completeness, note that the strong coupling α s depends on an a priory arbitrarily chosen renormalization momentum scale µ, i.e. α s = α s (µ). The µ dependence is such that α s (µ) ≪ 1 only for µ ≫ Λ QCD . A particularly obvious choice of the 1 In agreement with the prevalent notation, particularly in the field of lattice QCD, we use the terms static potential and static energy synonymously. Note, however, that in the literature sometimes a distinction is made and these terms refer to different quantities. 2 Exclusively focusing on QCD with just two dynamical quark flavors, in the framework of the present study we favor the specification of Λ MS rather than αs(MZ). A reasonable extraction of αs(MZ) from an n f = 2 study would at least require the discussion of flavor thresholds and extrapolations from n f → n f + 1, and in a sense obscure the main aim of our study. renormalization scale for the couplings in V (p) is µ = p, corresponding to an identification of µ with the typical momentum scale of the quantity under consideration.
Conversely, within lattice QCD the QQ static potential is naturally computed in position space. It can be extracted straightforwardly by studying the exponential decay of the rectangular Wilson loop as a function of its temporal extension [30]. Of course, lattice simulations at a given lattice spacing a cannot resolve arbitrarily small separations. Moreover, the minimum attainable lattice spacing is limited by the available computing power. The behavior of the static potential at small QQ separations is intimately related to its behavior at large momenta. Consequently, after a Fourier transform to momentum space, lattice results are expected to allow for reliable insights only below a certain threshold momentum, which is the maximum momentum p = π/a that can be resolved on a lattice of spacing a (cf. Sec. 2 below).
Taking into account the above constraints, a comparison and matching of perturbative and lattice results in momentum space is limited to momenta p fulfilling Λ QCD ≪ p ≪ p.
Analogous considerations can be invoked to delimit the fitting interval in position space (cf. e.g. [5]), employing that the typical momentum scale that can be attributed to a relative distance r scales as p ∼ 1/r. In position space the manifestly perturbative regime is thus characterized by 1/r ≫ Λ QCD . However, an important difference is that the leading uncontrolled contribution δV (r) to the static potential in position space as determined by a standard Fourier transform from momentum to position space is more pronounced than in momentum space: it is linear in Λ QCD , scales as δV (r) ∼ αs r (rΛ QCD ) = α s Λ QCD , and arises exclusively from the low momentum part of the Fourier integral over momenta Λ QCD , for which perturbation theory is no longer trustworthy [12] 3 . For the position space potential a seemingly obvious choice of the renormalization scale is µ = 1/r, which amounts to the typical momentum scale to be associated with the static quarks separated by a distance r.
However, it has been recognized long ago that in particular for the identification µ = 1/r the convergence of the perturbative potential in position space as defined by a standard Fourier transform from momentum space is spoiled [31], while it can be significantly improved by subtracting just the uncontrolled contribution δV (r) linear in Λ QCD [12]. Unfortunately this necessitates the specification of an additional subtraction scale and thereby increases the number of free parameters. Without any subtraction procedure a meaningful fit of perturbative expressions for the QQ static potential to lattice results in position space is not possible: the perturbative expressions only show a controlled convergence behavior for separations, which are much smaller than the minimum accessible separations on state-of-the-art lattice simulations (cf. e.g. [32][33][34][35]). Correspondingly, aiming at an accurate determination of Λ MS by comparing the results from lattice QCD simulations with perturbative calculations of the QQ static potential in position space, considerable efforts are needed to cope with this issue. Various strategies to retain or restore the significantly better convergence of V (p) also for V (r) have been devised in the literature [3,12,[32][33][34][35][36]].
An extraction of Λ MS directly in momentum space of course does not involve a Fourier transform of the perturbative potential, such that one may hope to circumvent most of these limitations. On the other hand, one now has to Fourier transform the lattice potential, which might even-tually lead to similar problems. As we will argue and demonstrate in detail in this paper, most favorable for us the latter concerns are not substantiated.
Our paper is organized as follows. Section 2 is devoted to the QQ static potential on the lattice. After briefly reviewing its computation in position space, we point out how we transform it to momentum space. In Sec. 3 we summarize the present knowledge of the perturbative QQ static potential and detail on the role of Λ MS . Special emphasis is put on the convergence behavior of both the perturbative QQ static potential for n f = 2 in momentum space and the QCD β-function. Section 4 constitutes the main section of our paper. Here we describe in detail our momentum space analysis to extract Λ 2 The QQ static potential in momentum space from lattice QCD

Gauge link configurations
We use the same n f = 2 gauge link configurations as for a recent determination of Λ MS , where, in contrast to this work, the lattice results and perturbative expressions for the static potential were compared and matched in position space [5]. These gauge link configurations were generated by the European Twisted Mass Collaboration (ETMC) [37][38][39] with the tree-level Symanzik improved gauge action [40], with b 0 = 1 − 8b 1 and b 1 = −1/12 and the Wilson twisted mass quark action [41][42][43][44], with Here a denotes the lattice spacing, ∇ µ and ∇ * µ are the gauge covariant forward and backward derivatives, m 0 and µ q are the bare untwisted and twisted quark masses, τ 3 is the third Pauli matrix acting in flavor space, and χ = (χ (u) , χ (d) ) represents the quark fields in the so-called twisted basis.
The twist angle ω is given by ω = arctan(µ R /m R ), where µ R and m R denote the renormalized twisted and untwisted quark masses. For the ensembles of gauge link configurations considered in the present study (cf. Table 1) ω has been tuned to π/2 by adjusting m 0 appropriately. This ensures automatic O(a) improvement for many observables including the static potential (cf. [38] for details).
The considered gauge link configurations cover several different values of the lattice spacing, the pion mass m PS and the spacetime volume L 3 × T ; cf. Table 1, which also provides the number of gauge link configurations, used for the computation of the static potential, for each ensemble. The lattice spacing in physical units has been set via the pion mass and the pion decay constant, using chiral perturbation theory. The resulting value for the hadronic scale 4 r 0 is r 0 = 0.420 (14) fm (cf. Sec. 5 of [38] and Tab. 8 of [39]). For further details on the generation of these gauge field configurations as well as on the computation and the analysis of standard quantities (e.g. lattice spacing and pion mass) we refer the reader to [38,39].

Computation of the QQ static potential in position space
First we determine the QQ static potential V ( r) in position space. In a second step, V ( r) is transformed to momentum space by means of a discrete Fourier transform. 4 The hadronic scale r0 is defined via r 2 0 F (r0) = 1.65, with F (r) = dV (r)/dr [45].  To be able to perform this Fourier transform numerically we need the static potential for all r = na inside a finite periodic spatial volume L ′3 of side length L ′ . This is achieved as follows: First, we compute V ( na) with n i ∈ {0, 1, . . . , N ′ /2}, where N ′ is even and defined as N ′ ≡ L ′ /a. Second, we realize the periodicity by defining V (n x a, n y a, n z a) ≡ V (|n x |a, |n y |a, |n z |a), where now To keep finite volume effects on a negligible level, the spatial volume L ′3 needs to be sufficiently large (cf. also point (2) below). While a lattice computation of the static potential V ( na) for all n x , n y , n z = 0, 1, . . . , N ′ /2 is possible in principle, it is extremely computer time consuming in practice: One needs to generate gauge link configurations for such large volumes and has to compute both on-and off-axis Wilson loops for all possible quark-antiquark separations r = na. On the other hand the shape of the static potential at large separations is known to be accurately described by V (r) = A 0 + σr + A 1 /r, where A 0 denotes a constant offset parameter, σ is the string tension and A 1 −π/12 [46][47][48].
Several clarifying comments are in order here. Accounting for a nonvanishing number n f of light dynamical quark flavors, there exists a certain threshold distance r c 1 fm such that for separations r > r c of the static quark Q and its antiquarkQ the QQ state is energetically disfavored in comparison to a pair of static-light mesons, BB [49]. This effect is known as string breaking. The ground state of the system is made up of dynamical quarks, gluons and static quarks Q andQ separated by a distance r. It scales string-like ∼ σr in the nonperturbative regime below r c , but becomes completely independent of r for r > r c , and saturates at about twice the B meson mass. However, the QQ state scaling ∼ σr can still be traced for r > r c also, where it corresponds to an excited state: Apart from resulting in a narrow mini-gap feature at r ≈ r c , i.e. where the energy of the QQ state equals that of the BB state, mixing effects between the Q and B sectors are tiny. Wilson loops as also employed here to extract the static potential on the lattice [cf. Eq. (5) below] are particularly insensitive to the B sector. On the other hand, standard perturbation theory for the QQ static potential manifestly focuses on the Q sector of the theory: While it accounts for dynamical light quarks in loop diagrams, in this framework the virtual light quarks can never become real. For these reasons, even though there occurs string breaking for n f = 0, when focusing on the QQ potential V (r) we manifestly limit ourselves to the Q sector, i.e. focus on the r dependent, linearly rising component of the potential also for r > r c .
Moreover, the large distance behavior of V (r) is expected to have a rather weak effect on the Fourier transformed potential for p ≫ Λ QCD , i.e. the momentum regime used for the matching to perturbation theory and the Λ MS determination in this work. 5 Therefore, we stick to the following strategy: For quark-antiquark separations | n|a ≤ r max we extract V ( r) from the exponential decay of Wilson loop averages W ( r, t) with respect to their temporal extent t, while keeping their spatial extent r fixed. To this end we first compute Somewhat arbitrarily, we choose r max ≈ 0.42 fm corresponding to r max = 10a for our smallest lattice spacing (cf. Table 1). In a second step the t-independent quantity V ( r) is obtained by performing an uncorrelated χ 2 minimizing fit to V (effective) ( r, t) in a suitable t range. This range is chosen such that excited states are strongly suppressed, while statistical errors are still small.
We use the ensembles listed in Table 1 and consider on-and off-axis Wilson loops formed by APE smeared spatial links (N APE = 60, α APE = 0.5 for all our gauge link ensembles) and ordinary, i.e . unsmeared temporal links. For a detailed explanation regarding the construction of off-axis Wilson loops cf. [5]. For a definition of APE smearing we refer to [50].
For quark-antiquark separations | n|a > r max while |n i |a ≤ L ′ /2 we model the lattice potential by For our finest lattice spacing a ≈ 0.0420 fm we use L ′ /a = 256, obviously fulfilling L ′ ≫ L (cf. Table 1). To ensure that the extracted value of Λ MS is independent of the choice for L ′ , we also performed computations with L ′ /a = 128 and L ′ /a = 512 and found essentially identical results for Λ MS (the deviations are below 1 MeV).
In Sec. 4 different values of M ∈ {1, 2, 3, 4} are used to quantify systematic errors associated with this modeling of the long range part of the lattice potential. The string tension is fixed to σ = 1550 MeV/fm (corresponding to r 0 = 0.420 fm determined on the same gauge link configurations we are using throughout this work [39]). While A 1 = −π/12 in the bosonic string picture [46,47], lattice simulations with n f = 2 quark flavors yield a larger value [48]. We determine A m with m ∈ {0, 1, . . . , M } by a χ 2 minimizing fit of Eq. (5) to the lattice results determined in step (1) in the region r min ≤ r ≤ r max . In order to have χ 2 1, for our smallest lattice spacing a = 0.0420 fm we choose r min = 8a for M = 1, r min = 6a for M = 2, and r min = 4a for M = {3, 4}.
The resulting function V 3 (r) for M = 3 is shown in Figure 1 together with the lattice results for V ( r).  Figure 1: Position space potential V 3 (r) (green curve) obtained by a χ 2 minimizing fit to the lattice results V ( r) (red dots) in the region 4a ≤ | r| ≤ 10a with a ≈ 0.0420 fm. Note that the statistical errors (actually also depicted here) are so tiny that the curve for V 3 (r) and its error band fall on top of each other and are indiscernible by eye.

Computation of the QQ static potential in momentum space
We define the QQ static potential in momentum space, V ( p) with p = 2π k/L ′ , by the discrete Fourier transform of V ( r), where the sum is also over all possible values of In the continuum and in infinite volume V ( p) is rotationally symmetric, i.e. V ( p) = V (p). Since both a cubic lattice discretization and a cubic periodic volume break rotational symmetry, this is no longer true for the lattice potential (6). The deviations from the rotationally invariant continuum and infinite volume potential are expected to be particularly small, when restricting the dimensionless lattice momenta k to values inside a cylinder of unit radius around the lattice diagonal, i.e. demanding This so-called cylinder cut is frequently used in lattice momentum space computations of propagators. For details cf. [51].
The maximal momentum value along each of the three axes is p = π/a ≈ 15 GeV for our smallest lattice spacing a ≈ 0.0420 fm. Similar to position space, where the minimum on-axis separation is a and the static potential is essentially free of discretization errors for r 3a, one might expect rather small discretization errors for p p/3 ≈ 5 GeV. These expectations will be confirmed in Sec  : V (p) for a ≈ 0.0420 fm obtained from V 3 (r) after applying the cylinder cut. 6 Throughout this paper all computations are performed in units of the lattice spacing a, e.g. we work with quantities aV (r), V (p)/a 2 and aΛ MS , which are independent of any potential errors or uncertainties regarding scale setting; for a recent review cf. [52]. Nevertheless, the axes in the presented plots as well as the numbers quoted in the main text often given in units of MeV or fm. This is intended to make the physical scales more obvious. To this end we use a = 0.0420 fm for our smallest lattice spacing corresponding to β = 4.35 (cf. also Table 1).
3 Perturbative insights into the QQ static potential 3.1 The QQ static potential in perturbation theory To allow for easier reference and to keep this paper self-contained, we briefly summarize the present status of the QQ static potential in perturbation theory. Quantities which depend on the particular renormalization scheme used are given in the MS scheme [53,54].
In the perturbative momentum regime, i.e. for momenta p = | p| ≫ Λ QCD , the QQ static potential is conveniently represented as with C F the eigenvalue of the quadratic Casimir operator for the fundamental representation of the gauge group; C F = 4/3 for SU (3). Pulling out the overall factor ∼ 1/p 2 the entire non-trivial structure of V (p) can be encoded in the dimensionless quantity α V , which in turn is a function of both the coupling α s (µ) and L ≡ L(µ, p) = ln µ 2 p 2 . α s (µ) is evaluated at an a priori arbitrarily chosen renormalization scale µ in the perturbative regime, i.e. µ ≫ Λ QCD , guaranteeing that α s (µ) ≪ 1.
The running of the coupling α s (µ) as a function of the renormalization scale µ is governed by the QCD β-function defined as whose series expansion in powers of α s is presently known with the following accuracy, While the expansion coefficients β 0 and β 1 are independent of the renormalization scheme, β 2 and β 3 are scheme-dependent. They have been determined for arbitrary compact semi-simple Lie groups in the MS scheme [55]. For SU(3) with n f = 2 massless dynamical quark flavors they read The same quantities for arbitrary values of n f can be found e.g. in [5].
As the static potential is a physical observable, it should of course be independent of the explicit value of the renormalization scale µ and form a renormalization group (RG) invariant, i.e. fulfill One might wonder, how this can come about with α V in Eq. (8) being a function of the two µ-dependent quantities α s (µ) and L(µ, p). However, knowing V (p) -and thus α V -at a certain accuracy in perturbation theory, e.g. up to O(αk s ), Eq. (12) only has to hold to this order, i.e.
Presently, all terms are known explicitly fork = 4, and α V is of the following form, While α V has a strict power-series expansion in α s up to O(α 3 s ), beyond this order one also encounters logarithmic contributions in α s [56], the first such term being ∼ α 4 s ln α s [57]. Equation (13) constrains the P k (L) accounting for the entire L dependence of α V in Eq. (14) to be polynomials in L of degree k, i.e.
with dimensionless expansion coefficients ρ km , and implies that, apart from the explicit values of a k ≡ ρ k0 = P k (0) with a 0 = 1, the ρ km are fully determined by the coefficients of the β-function [35,58]. Their explicit expressions for k ≤ 3 are given in our notations in [35].
Therewith, all coefficients in the perturbative expansion of the static potential in momentum space up to order α 4 s have been assembled. We emphasize again that the resulting expression is independent of the explicit choice for µ, in the sense that different choices for µ only lead to deviations at O(α 5 s ) [cf. Eq. (13)], which is beyond the accuracy of the contributions taken into account in Eq. (14).
In order to prevent the logarithms L in Eq. (14) from becoming large and thereby spoil the perturbative expansion, it is desirable to ensure that µ does not deviate much from p. Hence, a particularly convenient choice for the renormalization scale is µ ≡ p, implying L = 0. We will exclusively stick to this choice throughout the remainder of this paper. Adopting this choice, α V simplifies significantly. It becomes a function of α s (p) only and reads For completeness, note that -as briefly mentioned in the introduction -an analogous choice in position space, i.e. setting µ ≡ 1/r completely spoils the convergence of the perturbative potential V (r) [31]. As will also become obvious below, the choice µ ≡ p does not lead to any problems, the reason being the much less pronounced IR sensitivity of the perturbative potential V (p) in momentum space [12].
Let us here also mention the papers by [68,69] who argue that a particularly convenient choice of the RG scale µ is given by µ = p e −5/6 , rendering Eq. (14) independent of n f up to O(α 2 s ). While this choice is especially convenient when comparing the results for α V as an expansion in powers of α s (µ) for different numbers of dynamical quarks as the n f dependency is relegated to higher order expansion coefficients, for our analysis this choice has no advantages and does not provide a handle to improve the results: Even though the expansion coefficients of powers of α 1+n s (µ) show a slightly less pronounced increase with n for this choice (cf. also our detailed analysis for µ = p presented below), for a given value of p the explicit numerical value of α s (p e −5/6 ) is substantially increased as compared to α s (p). In the momentum regime where both lattice simulations and perturbative calculations for the QQ static potential are viable, the combined effect of these two opposite tendencies does not favor µ = p e −5/6 in comparison to µ = p.

An alternative representation of Eq. (19) is
α s (p) = 1 + x a 1 + (x a 1 ) 2 a 2 a 2 1 + (x a 1 ) 3 a 3 + a 3ln ln(4π/a 1 ) a 3 1 1 + a 3ln a 3 + a 3ln ln(4π/a 1 ) ln(xa 1 ) + . . . , (20) where we employed the shortcut notation x ≡ αs(p) 4π . In order to allow for more insights into the structure of α V [α s (p)] for n f = 2, we insert the explicit numerical values for the coefficients a 1 , a 2 , a 3 and a 3ln into Eq. (20), and rewrite it as with a 1 ≈ 8.11. It is easy to check that for xa 1 1 5 ↔ α s (p) 4π 5a 1 ≈ 0.31 the contributions ∼ (xa 1 ) n in Eq. (21) are ordered in the sense that they become increasingly less important, when increasing n from 0 to 3. For larger values of α s (p) 0.31 this ordering is spoiled. Correspondingly, the perturbative expansion of α V [α s (p)] with n f = 2 can in particular be considered as well-behaved and controllable for α s (p) ≪ 0.31. As the highest order contribution ∼ (xa 1 ) 3 in Eq. (21) is still significantly smaller than the next-to-highest one ∼ (xa 1 ) 2 , in the present paper we will consider the perturbative expansion as trustworthy even up to α s (p) 0.31. Taking into account the value of Λ MS as determined in [5], Λ

The scale Λ MS
So far we have not discussed, how an explicit numerical value can be attributed to the strong coupling α s (µ) evaluated at a given momentum µ. We emphasize that such an identification is renormalization scheme dependent, the most widely used scheme being the MS scheme [53,54], which we also adopt here.
A straightforward integration of Eq. (9) yields where the conventional definition of Λ in the MS scheme, i.e. Λ → Λ MS , corresponds to the choice of C = β 1 Employing some elementary manipulations and rearrangements, Eq. (22) adopted to the MS scheme can be written as involving a definite integral over α s (cf. e.g. [71]). The additional terms apart from the factor of 1/β(α s ) have been included in the integrand to make the finiteness of the integral over the interval from 0 to α s (µ) manifest.
In general Eq. (23) cannot be solved explicitly to provide the coupling α s (µ) at a given momentum scale µ as a function of the ratio µ/Λ MS . An exact closed form solution is only possible at leading order -i.e. taking into account only the leading contribution of the β-function (10), β(α s ) ≈ − αs 2π β 0 -and reads α s (µ) = β 0 2π ln(µ/Λ) −1 . Note, however, that approximate results for α s (µ) as a function of the ratio µ/Λ MS are available also for higher order contributions: the derivation of such expressions involves expansions in terms of ln(µ/Λ MS ); see e.g. [70,72].
Aiming at the determination of Λ MS by fitting the perturbative expression for the static potential (8) to numerical data from lattice simulations we do not see any reason to resort to these further approximations. Our strategy is rather to solve the implicit equation (23) for α s (µ) numerically. Hence, we will always use the relation (23) between α s (µ) and µ/Λ MS at the presently best available accuracy, irrespectively of the order of the expansion of the perturbative potential in α s (p).
However, note that there is still some freedom left in deciding, how to proceed with the evaluation of Eq. (23). We can (I) either plug in the perturbative expression of the β-function (10) at the presently best known accuracy and then do the integration over α s numerically, (II) or adopt a Taylor expansion of the integrand in Eq. (23) and do the integral analytically, keeping only those terms, whose coefficients are known explicitly, Of course, limiting ourselves to the truly perturbative regime, i.e. the regime, where the perturbative expansion of the β-function (10) is such that higher order corrections become increasingly less important, both choices should be equally justified.
As pointed out in Sec. 3, when fitting the perturbative static potential (8) in momentum space to lattice data we will always stick to the identification µ ≡ p. Correspondingly, higher order corrections in Eq. (10) should remain small throughout the integration interval from 0 to α(p) in Eq. (23) and Eq. (24), respectively.
To get a feeling, how the perturbative expansion of the β-function behaves, we employ the same strategy as adopted in the context of Eq. (21) and study its behavior at the upper integration limit in Eq. (23), i.e. at α s ≡ α s (p). First we rewrite Eq. (10) as with x ≡ αs(p) 4π . Second we insert the explicit numerical coefficients of the β-function for n f = 2 into this equation, resulting in with β 1 β 0 ≈ 7.93. Using the same reasoning as in the context of Eq. (21) above, we find that for the contributions ∼ (xβ 1 /β 0 ) n to become increasingly less important with n, we have to demand xβ 1 /β 0 0.4, which corresponds to α s (p) 0.63.
Hence, in particular for α s (p) 0.31 -which was the value found to crudely delimit the range of validity of a perturbative expansion of the static potential in momentum space for n f = 2 [cf. in the context of Eq. (21) above] -higher order terms in the perturbative expansion of the βfunction (10) for n f = 2 are expected to become much less important. In turn both possibilities (I) and (II) discussed above to numerically solve Eq. (23) for α s (p) as a function of p/Λ MS should yield very similar results.

Determination of Λ MS
In this section, we determine Λ MS in units of the lattice spacing, i.e. aΛ MS , by fitting perturbative expressions for the static potential in momentum space (cf. Sec. 3) to corresponding lattice results (cf. Sec. 2). Using the values of the lattice spacing listed in Table 1, these results can easily be converted to physical units, i.e. MeV. Unless explicitly stated otherwise, the errors provided for Λ MS do not account for the uncertainties associated with the lattice spacing a. These errors will however be accounted for when quoting our final result for Λ MS [cf. Eq. (31), below]. For completeness, we will also quote r 0 Λ MS , which is dimensionless and, hence, unaffected by any potential uncertainty in a (uncertainties in r 0 /a, which are collected in Table 1, are, of course, included).

Fitting procedures
The perturbative QQ potential to be fitted to lattice results is given by where we included an overall constant offset V 0 of the potential. We use different orders in the expansion of V (p) in powers of α s (p) to test and judge the convergence behavior of our results; the abbreviation N n LO stands for (next-to-) n leading-order. As only energy differences are measurable V 0 does not have any observable consequences. However, it is necessary to allow for a meaningful matching of the perturbative potential to lattice results. We employ an uncorrelated χ 2 minimizing fit with two degrees of freedom, V 0 and Λ MS , to fit the perturbative QQ potential (27) to the lattice potential in momentum space. Our fitting interval is delimited by p min and p max . Note that the lattice potential in momentum space originates from position space results whose large distance behavior is modeled by Eq.

Individual variation of input parameters
We investigate the stability of Λ MS with respect to variations of the input parameters p min , p max and M , using our finest lattice spacing a ≈ 0.042 fm and (L ′ /a) 3 = 256 3 lattice sites. To this end the input parameters delimiting the fitting range are varied in the following intervals: • p min = 1500 MeV . . . 2250 MeV: As discussed in detail in Sec. 3, for p min 1500 MeV the validity of perturbative expressions for the static potential is rather questionable.
• p max = 2250 MeV . . . 3000 MeV: The maximum momentum on our finest lattice along an axis is p = π/a ≈ 15 GeV. For p p/3, it seems reasonable to expect rather small discretization errors. These expectations have also been confirmed by numerical investigations. Fitting the perturbative potential to lattice results we obtain an essentially stable value of Λ MS up to 3000 MeV . . . Below, we demonstrate that the fit results for Λ MS are rather stable with respect to these parameter variations, i.e. we confirm that a meaningful and rather precise matching of perturbative expressions for the QQ static potential and lattice results is possible.
Exemplary fits of the perturbative static potential (27) at LO, NLO, NNLO and NNNLO to lattice results are depicted in Figure 3. Noteworthily, all four orders are suited to describe the lattice potential within statistical errors, i.e. fulfill χ 2 /dof 1.0. Similar χ 2 /dof are obtained when varying input parameters (cf. below).
To understand, how the value of Λ MS depends on the input parameters p min and p max , we vary them individually. Our results are summarized graphically in Figure 4  To allow for a meaningful determination of Λ MS , it is important to ensure that the fit results for Λ MS only exhibit a weak -preferably negligible -dependence on p min and p max . In plots of Λ MS as a function of p min and p max this should manifest itself in the formation of plateaus. This behavior is clearly visible in Figure 4 for all six plots and at all orders in the expansion (27). Moreover, the curves shown in Figure 4 clearly reveal a convergence behavior when increasing the order in Eq. (27)  2). Another source of uncertainty in our analysis is the modeling (5) used to extrapolate the lattice potential computed from Wilson loops to larger r. In order to scrutinize this uncertainty, we study the dependence of Λ MS on M ∈ {1, 2, 3, 4}, resorting to fitting procedure (I) and keeping p min = 1875 MeV and p max = 2625 MeV fixed. The corresponding results, which are collected in Table 2, indicate that Λ MS is also rather stable with respect to variations of M .   The variation of Λ MS as a function of the input parameters p min , p max and M (cf. also Figure 4 and Table 2) when fitting (27), both at NNLO (red) and at NNNLO (green), to lattice results is visualized in Figure 5. As systematic uncertainty one could e.g. quote the whole range of values covered by the NNNLO variations, In particular note the drastic improvement in stability as compared to the previous position space analysis [5], where a systematic uncertainty larger by a factor of ≈ 4 has been obtained, Λ The sources of the systematic error discussed above might, however, be correlated. A method to determine the overall systematic error accounting for potential correlation is discussed in the following subsection.

Consideration of correlations between different systematic error sources
To account for possible correlations, we perform a large number of fits, with the input parameters chosen randomly and uniformly in the intervals specified in Sec. 4.2.1. As systematic error we then take the variance of the fit results. This procedure is analogous to that used for a Λ MS determination from the QQ static potential in position space based upon the same lattice data [5] to which we will compare our results in the following.
In these fits we randomly vary p min and p max in the intervals p min = 1500 MeV . . .
Again it is most instructive to compare the error with that obtained for a position space analysis of the QQ static potential based upon the same lattice data [5]: • The final result of Ref. [5] is Λ (n f =2) MS = 315 (26) MeV. Hence, the momentum space analysis pursued in the present paper yields a result which is more precise by a factor of ≈ 2.
• While fitting procedure (A) used in [5] is rather similar to (I) and (II) of this paper, fitting procedure (B) used in [5] is somewhat different. 10 Hence, in order to allow for a consistent and fair comparison of our present momentum space analysis and the previous position space analysis, only fitting procedure (A) should be used. We carried out such an analysis, resulting in Λ (20) MeV. While there is now perfect agreement with respect to the average, the error associated with a momentum space analysis is still smaller by a factor of ≈ 1.5 than the error obtained in a position space analysis.  (13) MeV. In this case the error of Λ MS as obtained from a momentum space analysis is even smaller by a factor of ≈ 4 as compared to an analogous position space analysis. This rather drastic difference can be attributed to the almost perfect stability of Λ MS at NNNLO with respect to variations of the fitting range in momentum space (cf. Figure 4). In position space, a similar stability has not been observed [5].
Let us emphasize again that all these results for Λ Further systematic uncertainties, exclusively originating from the lattice static potential, are discussed in the next section.

Lattice discretization errors
To estimate the order of magnitude of lattice discretization errors, we utilize the lattice ensembles listed in Table 1, featuring pions of roughly of the same mass ≈ 284 MeV . . . 352 MeV. Keeping the input parameters p min = 1300 MeV, p max = 2050 MeV and M = 3 fixed, we extract values for both r 0 Λ MS and Λ MS by performing the corresponding fits. For these investigations we exclusively adopt fitting procedure (I) and use the perturbative QQ static potential (27) at the best known accuracy, i.e. NNNLO. 11 In Figure 6 we depict our results for r 0 Λ MS and Λ MS as a function of a 2 . The results determined at the four available lattice spacings perfectly agree within errors, i.e. there is no indication of any sizable lattice discretization errors. For completeness note that the errors associated with the extracted values for Λ MS in Figure 6 are rather large. The reason for this is that we have here included the errors associated with the lattice spacings and also with r 0 /a (cf. Table 1). This is absolutely essential when comparing results obtained at different values of the lattice spacing.

Finite volume effects and non-vanishing light quark mass corrections
Finite volume effects were investigated in detail in Ref. [5]. In this analysis such effects were found to be negligible compared to the other errors discussed above.
Similarly, potential corrections on Λ MS due to non-vanishing light quark masses were examined in detail in [5] by studying the variation of Λ MS for different pion masses in the range m PS ≈ 325 MeV . . . 517 MeV at fixed lattice spacing and spacetime volume. In the quark mass region investigated, Λ MS was found to be constant within tiny statistical errors of ≈ ±1 MeV. Therefore, we do not expect the non-vanishing light quark masses on the lattice to induce any significant deviations from the zero quark mass limit, for which the perturbative expressions in Sec. 3 were derived. In other words, we consider the systematic error introduced by comparing a lattice computation featuring massive light quarks with a perturbative calculation in the zero quark mass limit negligible compared to the other uncertainties discussed above.
For the above reasons we do not take any potential errors arising from finite volume effects and due to non-vanishing light quark mass corrections into account when quoting our final result for Λ MS in this paper.

Final results for Λ MS
In the following, we present our final results for Λ Let us emphasize again that the errors quoted here and in Eq. (31) also account for contributions which are not directly related to the extraction of Λ MS , as e.g. the uncertainties of a and of r 0 /a. Hence, with regard to a comparison of the results for Λ MS obtained from an analysis of the QQ static potential in momentum space as opposed to a similar analysis in position space, we consider the comparison in Sec. 4.2.2 -not accounting for these additional uncertainties -as most relevant and meaningful.
For easy reference, we collect all error contributions for our momentum space analysis in Table 3, where we also confront them to the analogous error contributions for a position space analysis [5].  Table 3: Error contributions for a momentum space analysis (this work) confronted to those for a position space analysis (Ref. [5]).

Conclusions
We have determined Λ (n f =2) MS by fitting perturbative expressions for the QQ static potential at the presently best know accuracy, i.e. up to O(α 4 s ) in momentum space, to lattice results obtained at a rather fine lattice spacing a ≈ 0.042 fm.
In contrast to previous works in this direction (cf. e.g. [1][2][3][4][5][6][7]) we have employed a discrete Fourier transform to transform the lattice results for the QQ static potential to momentum space. The extraction of Λ MS by fitting perturbative expressions of the static potential to lattice results has exclusively been performed in momentum space.
Resorting to a previous position space extraction of Λ (n f =2) MS based on the same lattice data [5], we could show that the momentum space analysis allows for a more precise determination of Λ MS . We could reduce the associated errors by a factor of ≈ 1.5 . . . 4, depending on the details of the fitting procedures and the estimate of the error associated with the truncation of the perturbative expressions used (cf. in particular Sec. 4.2.2). This improvement can mainly be attributed to the nearly perfect stability of the results for Λ MS with respect to variations of the momentum fitting range (cf. Figure 4). Such behavior has not been observed in position space, where comparably rather strong variations have been observed (cf. also Figures 3 and 4 of [5]).
In units of the hadronic scale r 0 our final result for Λ MS reads The quoted uncertainties include all sources of systematic error: the neglect of higher orders in the perturbative expansion, the dependence of the fit results on the fitting range p min . . . p max and the parameter M , lattice discretization errors, finite volume effects and uncertainties due to nonvanishing quark masses on the lattice. The uncertainties are dominated by the variations of Λ MS , when switching from NNLO to NNNLO, and in the case of Eq. (33) by the uncertainties associated with the lattice spacing.
In the future it would be interesting to adopt a similar momentum space analysis to lattice results for the QQ static potential with n f = 0, n f = 2 + 1 and n f = 2 + 1 + 1 dynamical quark flavors. For all these cases we expect a momentum space analysis to benefit from the better convergence behavior of the perturbative static potential in momentum space as compared to position space (cf. also the detailed discussion in the introduction, Sec. 1), which -as demonstrated in the present work -can reduce the error on the extracted values of Λ MS as compared to a position space analysis of the QQ static potential.