The confining color field in SU(3) gauge theory

We extend a previous numerical study of SU(3) Yang-Mills theory in which we measured the spatial distribution of all components of the color fields surrounding a static quark-antiquark pair for a wide range of quark-antiquark separations, and provided evidence that the simulated gauge invariant chromoelectric field can be separated into a Coulomb-like 'perturbative' field and a 'non-perturbative' field, identified as the confining part of the SU(3) flux tube field. In this paper we hypothesize that the fluctuating color fields not measured in our simulations do not contribute to the string tension. Under this assumption the string tension is determined by the color fields we measure, which form a tensor $F_{\mu \nu}$ pointing in a single direction in color space. We call this the Maxwell mechanism of confinement. We provide an additional procedure to isolate the non-perturbative (confining) field. We then extract the string tension from a stress energy-momentum tensor $T_{\mu \nu}$ having the Maxwell form, constructed from the non-perturbative part of the tensor $F_{\mu \nu}$ obtained from our simulations. To test our hypothesis we calculate the string tension from our simulations of the color fields for ten values of the quark-antiquark separation ranging from 0.37 fm to 1.2 fm. We also calculate the spatial distributions of the energy-momentum tensor $T_{\mu \nu}$ surrounding static quarks for this range of separations, and we compare these distributions with those obtained from direct simulations of the energy-momentum tensor in SU(3) Yang-Mills theory.

Received: date / Accepted: date Abstract We extend a previous numerical study of SU(3) Yang-Mills theory in which we measured the spatial distribution of all components of the color fields surrounding a static quark-antiquark pair for a wide range of quarkantiquark separations, and provided evidence that the simulated gauge invariant chromoelectric field can be separated into a Coulomb-like 'perturbative' field and a 'nonperturbative' field, identified as the confining part of the SU(3) flux tube field.
In this paper we hypothesize that the fluctuating color fields not measured in our simulations do not contribute to the string tension. Under this assumption the string tension is determined by the color fields we measure, which form a tensor F µν pointing in a single direction in color space. We call this the Maxwell mechanism of confinement.
We provide an additional procedure to isolate the nonperturbative (confining) field. We then extract the string tension from a stress energy-momentum tensor T µν having the Maxwell form, constructed from the non-perturbative part of the tensor F µν obtained from our simulations.
To test our hypothesis we calculate the string tension from our simulations of the color fields for ten values of the quark-antiquark separation ranging from 0.37 fm to 1.2 fm. We also calculate the spatial distributions of the energymomentum tensor T µν surrounding static quarks for this range of separations, and we compare these distributions with those obtained from direct simulations of the energymomentum tensor in SU(3) Yang-Mills theory.

Introduction
Quantum chromodynamics (QCD) is universally accepted as the theory of strong interactions. Nobody doubts that the well established phenomenon of confinement of quarks and gluons inside hadrons is encoded into the QCD Lagrangian. Yet, our current understanding does not go beyond that provided by a number of models of the QCD vacuum (for a review, see Refs. [1,2]). In particular, a theoretical a priori explanation of the so called area law in large size Wilson loops, which is closely related to a linear confining potential between a static quark and antiquark at large mutual distances, is still missing.
In such a challenging situation, first-principle Monte Carlo simulations of QCD on a space-time lattice represent an indispensable tool not only for checking (or ruling out) models of confinement, but also for providing new numerical "phenomenology" and possibly stimulating original insights into the mechanism of confinement.
Numerical simulations have established that there is a linear confining potential between a static quark and antiquark for distances equal to or larger than about 0.5 fm. This linear regime extends to infinite distances in SU(3) pure gauge theory, and, in the presence of dynamical quarks to distances of about 1.4 fm, where string breaking should take place [3][4][5]. The long-distance linear quark-antiquark potential is naturally associated with a tube-like structure ("flux tube") of the chromoelectric field in the longitudinal direction, i.e. along the line connecting the static quark and antiquark [6][7][8][9].
A wealth of numerical evidence of flux tubes has accumulated in SU(2) and SU(3) Yang-Mills theories .
Most of these studies concentrated on the shape of the chromoelectric field on the transverse plane at the midpoint of the line connecting the static quark and antiquark, given that the other two components of the chromoelectric field and all arXiv:1912.04739v1 [hep-lat] 8 Dec 2019 the three components of the chromomagnetic field are suppressed in that plane.
Recent times have witnessed an increasing numerical effort toward a more comprehensive numerical description of the color field around static sources, via the measurement of all components of both chromoelectric and chromomagnetic fields on all transverse planes passing through the line between the quarks [37]; of the spatial distribution of the stress energy momentum tensor [38,39]; and the flux densities for hybrid static potentials [40,41]. A more complete numerical description of the color field around the sources brings improved visualization, enabling us to grasp features otherwise less visible.
In the numerical study [37] we simulated the spatial distribution in three dimensions of all components of the chromoelectric and chromomagnetic fields generated by a static quark-antiquark pair in pure SU(3) lattice gauge theory. We found that, although the components of the simulated chromoelectric field transverse to the line connecting the pair are smaller than the simulated longitudinal chromoelectric field, these transverse components are large enough to be fit to a Coulomb-like 'perturbative' field produced by two static sources parameterized by effective charges ±Q of the sources (see Eq. (5) below).
The longitudinal component of this Coulomb-like 'perturbative' field accounts for a fraction of the simulated longitudinal chromoelectric field. We then identified the remaining longitudinal chromoelectric field as the confining 'nonperturbative' part of the simulated SU(3) flux tube field.
It is this non-perturbative part of the simulated field which contributes to the coefficient of the linear term in the heavy quark potential, the string tension.
In this paper we extend our simulations to a wider range of quark-antiquark separations. We extract the string tension from these simulations and compare our analysis with the results of recent simulations [38] of the energy-momentum tensor in SU(3) Yang-Mills theory.
We present a new procedure (the curl method) to extract a perturbative Coulomb field E C from the transverse components of the numerically simulated chromoelectric field. We avoid the use of a fitting function, directly imposing the condition that E C is irrotational (see Eq. (6) below). This provides a second method for implementing the underlying idea of our previous paper; that is, the chromoelectric field generated by a quark-antiquark pair can be separated into perturbative and non-perturbative components by a direct analysis of lattice data on the color field distributions generated by the pair.
As noted in [37], we can extract the value of the string tension from the non-perturbative field by utilizing the fact that the value of the chromoelectric field at the position of a quark is the force on the quark [42]. However, the Coulomblike field (Eq. (5)) does not give a good description of the transverse components of the chromoelectric field at distances closer than approximately two lattice steps from the sources [37], so that we must use the curl method to isolate the confining field in order to extract the string tension directly as the force.
The color fields F µν we measure, defined by the gauge invariant correlation function ρ conn W,µν (see Eqs. (1) and (2), below), point in a single direction in color space, parallel to the direction of the 'source' Wilson loop. In this paper we construct a stress energy -momentum tensor T µν having the Maxwell form from the 'measured' flux tube field tensor F µν , and extract the string tension (see Appendix A). This leads to a picture of a confining flux tube permeated with lines of force of a gauge invariant field tensor F µν carrying color charge along a single direction.
The 'Maxwell' energy-momentum tensor T µν does not account for the contribution to the quark-antiquark force from the fluctuating color fields not measured in our simulations. On the other hand, the complete Yang-Mills energy momentum tensor T YM µν simulated in Ref. [38] includes these fluctuating contributions, so that comparison of T YM µν with the 'Maxwell' energy-momentum tensor T µν constructed from the chromoelectric and chromomagnetic fields measured in our simulations provides a measure of the fluctuating contributions to the stress tensor.
We noted in our previous paper that the Coulomb-like 'perturbative' field (Eq. (5)) generated a stronger long distance Coulomb force between the heavy quarks than the Coulomb force measured in lattice simulations of the heavy quark potential [43][44][45] indicating the importance of fluctuations for the Coulomb contribution. In this paper we reexamine this issue.
The paper is organized as follows: in Section 2 we present the theoretical background and the lattice setup. In Section 3 we show some results on the spatial distribution of the color field around the two static sources and review the procedure to extract its non-perturbative part by subtraction of the Coulomb-like perturbative part identified by a fit of the transverse components of the chromoelectric field; in Section 4 we describe the new curl method to isolate the non-perturbative part; in Section 5 we show how to determine the string tension and other fundamental parameters describing the (non-perturbative) flux tube; finally, in Section 6 we discuss our results and give some ideas for future work.

Theoretical background and lattice setup
The lattice operator whose vacuum expectation value gives us access to the components of the color field generated by a static qq pair is the following connected correlator [14, 15, Here U P = U µν (x) is the plaquette in the (µ, ν) plane, connected to the Wilson loop W by a Schwinger line L, and N is the number of colors (see Fig. 1). The correlation function defined in Eq. (1) measures the field strength F µν , since in the naive continuum limit [15] ρ conn where qq denotes the average in the presence of a static qq pair, and 0 is the vacuum average. This relation is a necessary consequence of the gauge-invariance of the operator defined in Eq. (1) and of its linear dependence on the color field in the continuum limit (see Ref. [32]).
The lattice definition of the quark-antiquark fieldstrength tensor F µν is then obtained by equating the two sides of Eq. (2) for finite lattice spacing. In the particular case when the Wilson loop W lies in the plane withμ =4 andν =1 (see Fig. 1(a)) and the plaquette U P is placed in the planes41,42,43,23,31,12, we get, respectively, the color field components E x , E y , E z , B x , B y , B z , at the spatial point corresponding to the position of the center of the plaquette, up to a sign depending on the orientation of the plaquette. Because of the symmetry (Fig. 1), the color fields take on the same values at spatial points connected by rotations around the axis on which the sources are located (thê 1-or x-axis in the given example) .
As far as the color structure of the field F µν is concerned, we note that the source of F µν is the Wilson loop connected 0 0.5 −1 0    to the plaquette in Fig. 1. The role of the Schwinger lines entering Eq. (1) is to realize the color parallel transport between the source loop and the "probe" plaquette. The Wilson loop defines a direction in color space. The color field E that we measure in Eq. (2) points in that direction in the color space, i.e. in the color direction of the source. There are fluctuations of the color fields in the other color directions. We assume that these fluctuating color fields do not contribute to the string tension, so that the flux tube can be described as lines of force of the simulated field E.
The simulated flux tube field E carries color electric charge and color magnetic current along a single direction in color space. The divergence of E is equal to the color electric charge density ρ el (x) and the curl of E is equal to the color magnetic current density J M (x). The confining force is calculated from the divergence of a stress tensor T µν having the Maxwell form Eq. (A.1) .
The operator in Eq. (1) undergoes a non-trivial renormalization, which depends on x t , as discussed in a recent work [48]. The procedure outlined in that paper to properly take into account these renormalization effects is prohibitively demanding from the computational point of view for Wilson loops and Schwinger lines with linear dimension of the order of 1 fm, where the interesting physics is expected to take place. For this reason, we adopt here the traditional approach to perform smearing on the Monte Carlo ensemble configurations before taking measurements (see below for details). As shown in the Appendix A of our previous paper [37], smearing behaves as an effective renormalization, effectively pushing the system towards the continuum, where renormalization effects become negligible. The a posteriori validation of the smearing procedure is provided by the observation of continuum scaling: as carefully checked in Ref. [33], fields obtained in the same physical setup, but at different values of the coupling, are in perfect agreement in the range of parameters used in the present work. We performed all simulations in SU(3) pure gauge theory, with the standard Wilson action as the lattice discretization. A summary of the runs performed is given in Table 1. The error analysis was performed by the jackknife method over bins at different blocking levels.
We set the physical scale for the lattice spacing according to Ref. [43]: for all β values in the range 5.7 ≤ β ≤ 6.92. In this scheme the value of the square root of the string tension √ σ ≈ 0.465 GeV (see Eq. (3.5) in Ref. [43]).
The correspondence between β and the distance d shown in Table 1 was obtained from this parameterization. Note that the distance in lattice units between quark and antiquark, corresponding to the spatial size of the Wilson loop in the connected correlator of Eq. (1), was varied in the range d = 8 a to d = 16 a.
The connected correlator defined in Eq. (1) exhibits large fluctuations at the scale of the lattice spacing, which are responsible for a bad signal-to-noise ratio. To extract the physical information carried by fluctuations at the physical scale (and, therefore, at large distances in lattice units) we smoothed out configurations by a smearing procedure. Our setup consisted of (just) one step of HYP smearing [49] on the temporal links, with smearing parameters (α 1 , α 2 , α 3 ) = (1.0, 0.5, 0.5), and N APE steps of APE smearing [50] on the spatial links, with smearing parameter α APE = 0.25.

Spatial distribution of the color fields
Using Monte Carlo evaluations of the expectation value of the operator ρ W, µν over smeared ensembles, we have determined the six components of the color fields on all twodimensional planes transverse to the line joining the color sources allowed by the lattice discretization. These measurements were carried out for several values of the distance d between the static sources, in the range 0.37 fm to 1.19 fm, at values of β lying inside the continuum scaling region, as determined in Ref. [33].
We found that the chromomagnetic field is everywhere much smaller than the longitudinal chromoelectric field and is compatible with zero within statistical errors (see, e.g., Fig. 3 of Ref. [37]). As expected, the dominant component of the chromoelectric field is longitudinal, as is seen in Fig. 2, where we plot the components of the simulated chromoelectric field E at β = 6.370 as functions of their longitudinal displacement from one of the quarks, x l , and their transverse distance from the axis, x t .
While the transverse components of the chromoelectric field are also smaller than the longitudinal component, they are larger than the statistical errors in a region wide enough that we can match them to the transverse components of an effective Coulomb-like field E C produced by two static sources. For points which are not very close to the quarks, this matching can be carried out with a single fitting parameter Q, the effective charge of static quark and antiquark sources determining E C .
To the extent that we can fit the transverse components of the simulated field E to those of E C with an appropriate choice of Q, the non-perturbative difference E NP between the simulated chromoelectric field E and the effective Coulomb field E C , will be purely longitudinal. We then identify E NP as the confining field of the QCD flux tube.
To illustrate this idea, let us fix, for the sake of definiteness, β = 6.240 and put the two sources at a distance d = 16a = 1.02 fm. We then consider the plane, transverse to the longitudinal x-axis connecting the two sources, at a distance x l = 4a from one of them, and evaluate the components E x , E y and E z of the chromoelectric field in this transverse plane. The lattice determinations of E y on this plane can be fitted by the y-component of an effective Coulomb field where r Q and r −Q are the positions of the two static color sources and R 0 is the effective radius of the color source, introduced to explain, at least partially, the decrease of the field close to the sources. This fit is shown in Fig. 3 -see black dots and black solid line. Using the values of the fit parameters Q and R 0 obtained by the fit of E y , one can construct E C z and E C x and compare them to lattice data. Furthermore, the Coulomb-like content of E z fully accounts for the z-component of the chromoelectric field (see red dots and red solid line in Fig. (3)), but E C x accounts for only a fraction of the longitudinal component of the chromoelectric field (see blue dots and blue solid line in Fig. (3)). This strongly suggests that the non-perturbative component of the chromoelectric field is almost completely oriented along the longitudinal direction. It can be isolated once the parameters of the Coulomb-like component are determined by a fit to the yand/or z-components of the lattice determination of the chromoelectric field.
The procedure we have just illustrated in a specific case, can be carried out in a systematic manner. We observe that in making the fit we must take into account that the color fields are probed by a plaquette, so that the measured field value should be assigned to the center of the plaquette. This also means that the z-component of the field is probed at a distance of 1/2 lattice spacing from the x y plane, where the z-component of the Coulomb field E C z is non-zero and can be matched with the measured value E z for the same value of Q. For further details about the fitting procedure and the extraction of the fit parameters we refer to Appendix B of Ref. [37].
In Table 2, we list the values of the effective charge Q obtained from the lattice measurements of E z and E y at the values of d, the quark-antiquark separation, considered in this work. The statistical uncertainties in the quoted Q values result from the comparisons among Coulomb fits of E y and E z at the values of x l , for which we were able to get meaningful results for the fit. The values of R 0 in physical units grow with the lattice step a, while in lattice units they show more stability. This suggests that the effective size of a color charge in our case is mainly explained by lattice discretization artifacts and the smearing procedure, and is not a physical quantity (see Appendix B of Ref. [37]).
Evaluating the contribution of the field of the quark to E C in Eq. (5) at the position r −Q of the antiquark and multiplying by the charge −4πQ of the antiquark yields a Coulomb force between the quark and antiquark with coefficient −4πQ 2 . By comparison, the standard string picture of the color flux tube gives a Coulomb correction of strength −π/12 to the long distance linear potential (the universal Lüscher term arising from the long wave length transverse fluctuations of the flux tube [51]). In addition, the strength π 12 of the Luscher term is approximately equal to the strength of the Coulomb force determined from the analysis of lattice simulations of the heavy quark potential at distances down to about 0.4 fm. [43,45].
By contrast, the strength −4πQ 2 of the Coulomb force generated by E C is roughly 4 times larger than π 12 for the values of the effective charge Q listed in Table (2) and determined from our simulations of ρ conn W,µ,ν . Therefore the fluctuating color fields not measured in our simulations must be taken into account in calculating the Coulomb correction to the long distance heavy quark potential.
In Fig. 4 we plot the longitudinal component E NP x of the non-perturbative field in Eq. (4) as a function of the longitudinal and transverse displacements x l , x t at β = 6.370. As expected, E NP x is almost uniform along the flux tube at distances not too close to the static color sources. This feature is better seen in Fig. 5, where transverse cross sections of the field E NP x (x l , x t ), plotted in Fig. 4, are shown for the values of x l specified in Fig. 5. For these values of x l the shape of the non-perturbative longitudinal field is basically constant all along the axis. A similar scenario holds in the other lattice setups listed in Table 1. x l = 3.5 a x l = 4.5 a x l = 5.5 a x l = 6.5 a x l = 7.5 a x l = 8.5 a x l = 9.5 a x l = 10.5 a x l = 11.5 a x l = 12.5 a  Table 2 Values of the fit parameters Q and R 0 extracted from Coulomb fits of the transverse components of the chromoelectric field and values of the longitudinal chromoelectric fields at (d/2, 0), the midpoint between the sources and transverse distance zero, for several values of distance d. E x (d/2, 0) is the unsubtracted simulated field and E NP x (d/2, 0) is the non-perturbative chromoelectric field. For comparison, in the last column the non-perturbative chromoelectric field E NP x curl (d/2, 0) obtained using the irrotational property of the perturbative field (see section 4) is given. For the parameters of the Coulomb fit we quote, along with the statistical error, a systematic uncertainty that accounts for the variability in the values of the fit parameters extracted from all acceptable fits to E y and E z at different x l values (for more details, see Appendix B of Ref. [37]). As the distance between the sources is made smaller and smaller the quality of the Coulomb fits deteriorates and Q and R 0 cannot be reliably extracted for d ≤ 0.51 fm.  (4) In Table 2 we also compare the values of the measured longitudinal chromoelectric field E x with those of the nonperturbative field E NP x on the axis at the midpoint between the quark and antiquark, for all ten values of their separation d. Given that E NP x is almost uniform along the axis,  While the Coulomb field (5) gives a good description of the transverse components of the chromoelectric field when the distance from the sources is not too small, it does not give a good description at smaller distances, approximately two lattice spacings from the sources. This can be either the result of the non-spherical form of the effective charges, or an effect introduced by the discrete lattice.
To extract the confining part of the chromoelectric field in the data it is then preferable to have a procedure which avoids the use of an explicit fitting function, and which can work close to the quark sources. With this aim in mind we use the following two steps to separate the field into 'perturbative' and 'non perturbative' components.
1. We identify the transverse component E y of the simulated field with the transverse component E C y of the perturbative field, E C y ≡ E y . 2. We impose the condition that the perturbative field is irrotational, curl E C = 0.
Condition (1) implies that the nonperturbative field is purely longitudinal, E NP y = 0. Condition (2) will then fix the longitudinal component E C x of the perturbative field as well as the longitudinal component E NP To implement the irrotational condition (2), taking into account that the fields are measured at discrete lattice points, the sum of the measured fields along any closed lattice path is zero. For example, on a plaquette this amounts to One can easily solve this equation for E C x obtaining This of course leaves one unknown on each transverse slice of the field -the value of E C x (x, y max + 1), but if the value of y max is large enough, the perturbative field at that distance should already be small, so in our analysis we just put E C x (x, y max + 1) = 0. To check that this indeed makes a little change to our results, we have used a separate procedure in which we fixed E C x (x, y max + 1) = E x , in practice making E NP x = 0 at the largest transverse distance. This procedure gave similar results.
After the estimation of the perturbative longitudinal field E C x one can subtract it from the total field, obtaining the non-perturbative component (see Fig. 6). One can see that the non-perturbative part of the flux tube exhibits very little change along the line connecting the quark-antiquark pair; even at the smallest distances from the sources the nonperturbative field remains smooth (This is seen more clearly in Fig. 7). Table 3 The Clem parameters describing the non-perturbative field transverse section going through the midpoint between the quark and antiquark positions. The data for the fit is obtained using the curl subtraction method, taking the perturbative field at y max+1 equal to zero.  Table 4 The string tension estimated using the non-perturbative field from the curl procedure by employing different methods (from left to right: numerical integration of the field, analytical integration of the Clem function with parameters given in Table 3, estimation of fields at sources). In the last column we report also the value of the string tension obtained by numerically integrating Eq. (8) and using the non-perturbative field from the Coulomb subtraction E Coulomb  Table 5 The flux tube width estimated using the non-perturbative field from the curl procedure by employing different methods (from left to right: numerical integration of the field, analytical integration of the Clem function with parameters given in Table 3). In the last column we report also the value of the width obtained by numerically integrating Eq. (12) and using the non-perturbative field from the Coulomb subtraction E Coulomb

The string tension and the width of the flux tube
The forces between charged particles in electrodynamics are determined by a stress tensor T µν constructed from fields F µν satisfying Maxwell's equations (see Eq. (12.113) in Ref. [52]). Similarly, the force between quarks and antiquarks in Yang Mills theory is determined by the stress tensor T µν , Eq. (A.1), constructed from the field tensor F µν obtained from our simulations.
The quark-antiquark force F is then the integral of the longitudinal component T xx = (E x (x l = d/2, x t )) 2 /2 of the stress tensor over the median plane x = d/2 bisecting the line connecting the quarks, Eq. (A.8). The non-perturbative quark-antiquark force F NP = −ê x σ determining the string tension σ has the corresponding expression in terms of the non-perturbative longitudinal component of the stress tensor The square root of the string tension is then equal to We have evaluated the integral (Eq. (8) ) in two ways: 1. by direct numerical integration, using the values of E NP x determined by our simulations, and 2. analytically, by fitting the numerical data for the transverse distribution of E NP x (x t ) as in [25][26][27][28][29] to the Clem parameterization of the field surrounding a magnetic vortex in a superconductor [53]. where φ , µ and α are fitting parameters. In the dual superconducting model [54-57] λ = 1 µ is the penetration depth and is the Landau-Ginzburg parameter characterizing the type of superconductor.
Figure (8) shows an example of the fit of the data to the Clem functional form Eq. (9). for the transverse distribution of E NP x (x t ), obtained using the curl procedure.
We can obtain a second expression for the string tension by utilizing the result [42,58] that the force on a quark is equal to the value of the chromoelectric field at the position of the quark. The string tension is then equal to the corresponding value of the confining part of the chromoelectric field Eqs. (8) and (11) provide two independent ways to extract the string tension from simulations. As mentioned earlier, we must use the curl method to isolate the confining field to extract the string tension by Eq. (11). To obtain additional information about the structure of the chromoelectric flux tube we have calculated the mean square root width: .
Just as we have evaluated the integral (8) for the string tension, we have evaluated the integral (12) for the mean square root width both numerically, using the data for E NP x (x t ), and analytically, in terms of Clem parameters, fitting the longitudinal component of E NP x (X t ) in the median plane to the Clem parametrization ( (9)). The results of that fit are given in Table 3. In most cases the parametrization in Eq. (9) gives a good description of the field, shown by the values of χ 2 r , though the parameters themselves are somewhat unstable, which reflects the strong correlation between the parameter estimates.
We compared two different methods for calculating the integrals in Eqs. (8) and (12). First, we carried out the numeric integration in Eqs. (8) and (12), respectively, postulating the rotational symmetry of the field. This approach was repeated for both the non-perturbative field obtained using the "curl procedure" (resulting in √ σ int and w 2 int ) and the field obtained using Coulomb subtraction (resulting in √ σ Coulomb and w 2 Coulomb ). Next we calculated the values of the string tension and the mean square root width using the Clem parameters given in Table 3 to get the values denoted as √ σ Clem and w 2 Clem . (One remark should be made for the width -while we know that the value of E C x (x, y max + 1) in Eq. (7) is small (O(y −2 max )), in the numerator of Eq. 12 this small constant will be multiplied by y 2 (y 3 after the integration over polar angle), which will cause the error introduced to increase with y max . Indeed, the comparison with the analysis done taking E NP x (x, y max ) = 0 shows large discrepancies in this case.) Finally, we evaluated expression (11) for the string tension, σ = E NP x (x l = 0, x t = 0), using the curl method to determine the magnitude of the non-perturbative field at the sources.
Our results are gathered in Tables 4 and 5, where we use the notation σ 0 ≡ E NP y (x l = 0, x t = 0) and in Fig. (9). The data shown in Fig. (9) give a consistent value of √ σ for all values of the separation d, with scatter that increases with d as the resolution diminishes. The values of √ σ lie close to 0.465 GeV, the value used in the parameterization Ref. [43] Let us review the basis of our calculations. Our hypothesis is that the string tension is determined by the field E we measured (the 'Maxwell' mechanism). We have determined σ from both the transverse structure of the flux tube (Eq. (8)) and its longitudinal structure (Eq. 11) ) as shown in Fig. (6(c)), in which the non-perturbative field has been isolated.
We emphasize that, as discussed in Section (3), the 'Maxwell' mechanism cannot be used to obtain the Coulomb correction to the string tension. This implies that the fluctuating fields not measured in our simulations must contribute to the Coulomb force.
On the other hand, the Coulomb correction has been obtained by recent direct simulations of the stress energymomentum tensor in Yang Mills theory [38]. The Yang Mills stress tensor accounts for the contributions of fluctuating fields but cannot be directly related to measured fields, in contrast to the Maxwell stress tensor T xx = (1/2)E 2 x , determining the string tension. (See Eq. (A.9).)

Conclusions and outlook
In this paper we have determined the spatial distribution in three dimensions of all components of the color fields generated by a static quark-antiquark pair. We have found that the dominant component of the color field is the chromoelectric one in the longitudinal direction, i.e. in the direction along the axis connecting the two quark sources. This feature of the field distribution has been known for a long time. However, the accuracy of our numerical results allowed us to go far beyond this observation. First, we could confirm that, as observed in [37], all the chromomagnetic components of the color field are compatible with zero within the statistical uncertainties. Second, the chromoelectric components of the color fields in the directions transverse to the axis connecting the two sources, though strongly suppressed with respect to the longitudinal component, are sufficiently greater than the statistical uncertainties that they can be nicely reproduced by a Coulomb-like field generated by two sources with opposite charge (everywhere except in a small region around the sources).
In Ref. [37] we subtracted this Coulomb-like field from the simulated chromoelectric field to obtain a non-perturbative field E NP according to Eq. (4) and found that the dependence of the resulting longitudinal component of E NP on the distance x t from the axis is independent of the position x l along the axis, except near the sources, thus suggesting that the non-perturbative field found in this way from lattice simulations can be identified as the confining field of the QCD flux tube.
In this work we have improved the approach of Ref. [37] by presenting a new procedure to subtract the Coulomb-like field, which does not rely on any preconception about its analytic form, but is based only on the requirement that its curl is equal to zero.
Moreover, we have carefully analyzed the spatial distribution of the subtracted, non-perturbative part of the longitudinal chromoelectric field to extract from it some relevant parameters of the flux tube, such as the mean width and the string tension, both by means of a fully numerical, modelindependent procedure and by a prior interpolation with the dual version of the Clem function for the magnetic field in a superconductor.
We have also used our determinations of the color field components to construct the 'Maxwell' stress tensor. Details about its determination and a comparison with the recent literature about this topic [38] are presented in Appendix A.
In conclusion, we have shown that the separation of the chromoelectric field into perturbative and non-perturbative components can be obtained by directly analyzing lattice data on color field distributions between static quark sources, with no need of model assumptions. To the best of our knowledge, this separation between perturbative and non-perturbative components has not been carried out previously. It provides new understanding of the chromoelectric field surrounding the quarks. We have used the nonperturbative field to calculate the string tension and the spatial distribution of the energy-momentum tensor surrounding the static quarks, under the assumption that the fluctuating color fields not measured in our simulations do not contribute to the string tension. The extension of our approach to the case of QCD with dynamical fermions with physical masses and at non-zero temperature and baryon density is straightforward [32].

Appendix A: The 'Maxwell" stress tensor
In this Appendix we consider the "Maxwell" energy-momentum tensor T µν as a function of the field tensor F µν characterizing the SU(3) flux tube, which is in its turn defined in terms of the gauge invariant correlation function ρ conn W,µν of Eqs. (1) and (2) and points in a single color direction parallel to the color direction of the source (which is determined dynamically). Its six tensor components (the electric and magnetic fields E and B) correspond to the six orientations of the plaquette U P relative to that of the Wilson loop (see Fig. 1(a)).
The simulated fields E and B have the space-time symmetries of the Maxwell fields of electrodynamics, while carrying color charge in a single direction in color space. The energy-momentum tensor T µν lies in the same direction in color space of the simulated fields E and B and has the (Euclidean) Maxwell form: Its spatial components µ = i, ν = j, with i, j = 1, 2, 3 determine the Maxwell stress tensor: Taking µ = i and ν = j = i in Eq. (A.1) gives while the diagonal time component −T 44 of T µν determines the energy density, We use cylindrical coordinates (x, r, θ ), r ≡ y 2 + z 2 , tan θ ≡ z/y, and the corresponding unit vectorsê r ,ê θ : e r =ê y cos θ +ê z sin θ , (A.5) e θ = −ê y sin θ +ê z cos θ (A.6) (x is the longitudinal direction of the flux tube, i.e. the axis along which the static sources are located). The force exerted by the antiquark on the quark can be expressed, by means of the stress tensor, as a surface force F acting on the infinite plane x = d 2 bisecting the line connecting the pair: , and taking into account that the measured magnetic field B is compatible with zero The angular average over the radial vectorê r in Eq. (A.10) vanishes. Furthermore by symmetry the transverse field E r on the mid-plane x = d 2 vanishes, so that the quark-antiquark force in Eq. (A.10) becomes Replacing E x (r) by the non-perturbative field E NP x (r) in Eq. (A.11) gives the non-perturbative quark-antiquark force F NP , The remaining non-vanishing component of T µν is the energy density T 44 , Further, we note that the trace of the stress tensor T µν evaluated from Eqs. (A.13)-(A.17) vanishes independently of the simulated flux-tube fields E x (x, r) and E r (x, r): We have calculated the non-perturbative content of T rr on the symmetry plane (where T rr = T θ θ = −T 44 = −T xx ) versus r for three different values of the quark-antiquark distance: d = 0.51 fm (at β = 6.240), d = 0.69 fm (at β = 6.539) and d = 0.95 fm (at β = 6.299). Results are presented in Fig. 11, where also the full (non-perturbative plus Coulomb) content of T rr is shown.
The width of the energy density distribution T NP 44 can be obtained through Eq. (12), with E NP x replaced by T NP 44 as given in Eq. (A.18); results are presented in Table 6. Since T NP 44 is proportional to E NP x 2 the width of the T NP 44 component of the Maxwell stress tensor obtained from the nonperturbative field given in Table 6 is systematically smaller than the width of the nonperturbative part of the longitudinal chromoelectric field component E NP x given in Table 5. (The square of the field decreases more rapidly with distance than the field itself.) We now compare the above results obtained using our measured flux tube fields to evaluate the 'Maxwell' energymomentum tensor T µν with the corresponding results of recent direct simulations [38] of the expectation value of the energy momentum tensor T YM µν in the presence of a quarkantiquark pair. The latter simulations, which measure the energy and stresses in all color directions directly, were carried out in the plane midway between the quark and the antiquark, for three values of their separation.
The tensor T YM µν has the form [59] T µν = F a µα F a αν − g µν F a αβ F a αβ /4 , (A. 21) where F a µν is the Yang-Mills field tensor in the adjoint representation of SU (3), where f abc are the structure constants of the SU(3) algebra. In our definition the field is squared after color projection, whereas in Eq. (A.21) the sum over color components is taken after squaring. Moreover, the stress tensor in Ref. [38] is renormalized (this motivates the superscript R in the formulas below). In [38] the expectation value of the energy-momentum tensor in the background of a quark-antiquark pair is de-  noted by T R i j (r) QQ . We will use this notation in comparing our results (A.18) and (A. 19) with that work. In [38] the r dependence of the components of T R µν (r) QQ was plotted for the 3 values of the quarkantiquark separation for which simulations were made, and these 'noticeable' features of the results were pointed out: 1. Approximate degeneracies between temporal and longitudinal components and between radial and angular components are found for a wide range of r; T R 44 (r) QQ ≈ T R xx (r) QQ > 0 , T R rr (r) QQ ≈ T R θ θ (r) QQ > 0. (A. 23) We emphasize that the two inequalities in Eq. (A.23) are general consequences of (A.18) and (A. 19), independent of the values of the simulated field E x (x, r). In contrast, a recent study [39] of the stress tensor distribution in the + T R rr (r) QQ + T R θ θ (r) QQ < 0, is partially restored inside the flux tube. 3. Each component of the energy-momentum tensor at r=0 decreases as the separation becomes larger, while the transverse radius of the flux tube, typically about 0.2 fm, seems to increase for large separations [30,60,61], although the statistics are not sufficient to discuss the radius quantitatively. We see some indication of the increase in the width of the distributions of the diagonal components of the 'Maxwell' stress tensor in Fig. (12) and Table (6) for all ten values of the quark-antiquark separation. However, this width is greater than 0.2 fm, the transverse radius of the flux tube estimated by [38]. which can be clearly seen from Fig. 3 of Ref. [38], where the components of T R i j (r) QQ were plotted.
The 'Maxwell' stress tensor does not include the contributions to Eq. (A.21) of the fluctuating color fields not measured in our simulations. Comparison of the spatial distributions of the diagonal components of the Yang-Mills stress tensor with the corresponding distributions of the 'Maxwell' stress tensor then provides a measure of the contributions of the fluctuating color fields.