Effects of a strong magnetic field on the QCD flux tube

In this work we investigate the effect of an external magnetic field B on the shape of flux tubes in QCD by means of lattice simulations, performed with N_f=2+1 flavors of stout improved dynamical staggered quarks with physical masses. After having discussed some difficulties in the practical definition of the flux tube at B=0, we show that these ambiguities do not affect the determination of the flux tube modifications induced by the magnetic field. Different results are obtained depending on the relative orientations of the flux tube and of the magnetic field: they confirm that the magnetic field acts as transverse confinement catalyser and longitudinal confinement inhibitor; moreover, the flux tube itself loses its axial symmetry when it is not directed along the magnetic background.


I. INTRODUCTION
Despite the fact that a formal proof of color confinement in Quantum Chromodynamics (QCD) is still lacking, Lattice QCD simulations provided an overwhelming amount of numerical evidence that color confinement is encoded in the QCD Lagrangian. Particularly important in this respect was the observation [1][2][3][4][5][6][7] that the field generated by two opposite static color sources is not a dipole field like in quantum electrodynamics: most of the field energy density is concentrated in a linear structure that connects the two static sources. The features of this linear structure closely resemble the ones of the flux tubes experimentally observed in type II superconductors, and for this reason the linear structure was called color flux tube. Analogous linear networks are observed when three or more static color sources are used [8][9][10][11].
The existence of the color flux tube provides an intuitive explanation for the linearly rising potential between two opposite static color sources: the slope of the potential (the string tension σ) is nothing but the energy density per unit length of the flux tube. As a consequence the study of the color flux tube imposed itself as a tool to investigate the origin of the confining potential in QCD in a way that is independent of the details of the confining mechanism, even though the very idea of flux tube emerges very naturally within the dual superconductor scenario for color confinement [12,13].
The purpose of this paper is to provide a first Lattice QCD investigation of flux tubes in the presence of a magnetic background field. Various lattice studies have shown that such an external magnetic field has a strong influence on the confining properties of QCD [14][15][16][17], with the string tension in the direction parallel to the magnetic field that is strongly reduced and which, for large enough magnetic fields, could even disappear. Magnetic field induced anisotropies could play a relevant role, especially at the level of heavy quark phenomenology [18][19][20][21][22][23][24][25][26][27][28][29]. A possible interpretation of such results can be found in various model computations [30][31][32][33][34][35][36][37][38][39][40][41][42][43], and looking at the flux tube provides a way to achieve a better comprehension of the specific way in which the magnetic field influences the confining properties of QCD.
A priori various phenomena could indeed take place: the magnetic field could change the strength of the color field within the flux tube, but it could also modify the shape of the flux tube itself, which could even become anisotropic and lose its axial symmetry when the quarkantiquark separation is not collinear with the magnetic field. Even if the tube profile remains cylindrical, one can look at the characteristic lengths that characterize the flux tube and inquire how they are modified by the presence of an external magnetic field.
After the introduction in Sec. II of the numerical setup adopted and of the observables used, a preliminary part of our study is dedicated to the investigation of the flux tube at zero magnetic field (Sec. III). Indeed, most studies in the literature have been performed in pure gauge theories and only recently results for full QCD appeared [44]. Therefore a study of the flux tube properties in N f = 2 + 1 QCD at the physical point, using a discretization different from the one adopted in Ref. [44], is interesting by itself.
This preliminary part will give us the opportunity of discussing some ambiguities related to the flux tube defi-nition, that are associated with the smoothing procedure adopted to improve the signal to noise ratio. Remarkably these ambiguities are much less severe (and indeed practically absent) if one is interested only in modifications of the flux tube induced by the magnetic field. This will be shown in Sec. IV, where the main results of this paper will be presented 1 . Finally, in Sec. V we report our conclusions.

II. NUMERICAL SETUP
A. Lattice discretization of N f = 2 + 1 QCD with a magnetic background In this work we simulate 2 + 1 flavour QCD making use of the stout improved rooted staggered fermion discretization and the Symanzik tree-level improved gauge action. More explicitly, the partition function is written as where DU stands for the product of the SU (3) Haar measure of all the links of the lattice. The gauge action S YM is given by [46,47] where the symbols P 1×1 i; µν and P 1×2 i; µν denote the real part of the trace of 1×1 and 1×2 Wilson loops. The staggered Dirac matrix is where the η i; ν s are the usual staggered phases, U (2) i; µ stands for the two times stout-smeared link in position i and direction µ [48] (with ρ = 0.15 as isotropic smearing parameter) and u f i; µ is the abelian field parallel transporter.
The transporters corresponding to a uniform magnetic field B z directed alongẑ can be written as where q f is the quark charge and all the other abelian link variables are set to 1 (N k is the lattice extent in the 1 Preliminary results have been presented at the 35 th International Symposium on Lattice Field Theory (Lattice 2017) [45].
directionk and 1 ≤ i k ≤ N k ). However in this expression the value of B z cannot be arbitrary: for Eq. (4) to describe a uniform magnetic field on a lattice with periodic boundary conditions, the value B z has to satisfy the quantization condition [49][50][51] where b is an integer number. Let us stress that the magnetic field is external: abelian transporters u f i; µ are not updated, so quarks interact with the external magnetic field but they do not back-react on it. In this way we are neglecting the direct quark-quark electromagnetic interactions, while we are properly taking into account the effect of the external field on the quark loops.
Bare parameters have been chosen in such a way that simulations stay on a line of constant physics with physical values of the quark masses. Since the lattice spacing is independent of the magnetic field [52], we could use for this purpose the values reported in Refs. [53][54][55]. Gauge configurations have been sampled by using the Rational Hybrid Monte-Carlo (RHMC) [56][57][58] algorithm; simulation parameters and details are reported in Table I.

B. Observables
To study the color flux tube in a lattice simulation two basic ingredients are needed: two static color sources of opposite charge, that are usually introduced by means of a Wilson loop (or two Polyakov loops if the temperature is non-vanishing), and a probe to investigate the field structure, that is usually a plaquette, i.e. a Wilson loop of size 1 × 1. The field profile around the static charges is extracted from the correlation of the plaquette with the 1.20 10 TABLE I: Simulation details: the bare coupling β, the bare light quark mass m l (ms/m l was always fixed to 28.15), the lattice size, the value of the lattice spacing (with a systematic uncertainty of 2-3%, see [53][54][55]), the magnetic quantum b, the magnetic field intensity eB and the number of independent configurations N confs analyzed, separated by 10 ÷ 25 molecular dynamics trajectories.
in the second case [59,60] the connected correlator is used instead where N c is the number of color and L is a parallel transporter that connects the Wilson loop and the plaquette (which is often called "Schwinger line"), see Fig. 1 for a pictorial representation. Note that in both cases, given that in the present case the electromagnetic field is not dynamical, a possible U (1) charge of the static pair is irrelevant, since it would just add a phase to the Wilson loop which cancels in the ratios, while the plaquette operator is by definition limited to the SU (3) part, since we are interested in the modifications of the non-Abelian fields in the flux tube induced by the external magnetic background. These two lattice implementations are not completely equivalent, as they probe different quantities: in the naive continuum limit the disconnected correlator ρ disc reduces to where µ and ν identify the plaquette orientation (no sum on µ, ν is intended in the r.h.s), F µν is the continuum euclidean field strength, a is the lattice spacing and g 0 is the bare coupling constant. The subscripts QQ and 0 are used in the previous expression to denote the cases in which two opposite static color charges are present in the background or not. When using the connected correlator ρ conn , the plaquette operator appears in the same trace of the Wilson loop operator, so in the naive continuum limit the term that is linear in the nonabelian field strength at the position of the probe gives the leading contribution an expression that in the literature is often denoted, for the sake of simplicity but with a clear abuse of notation, by a 2 g 0 F µν QQ . While ρ conn and ρ disc display similar behaviors as a function of the transverse displacement x t , the connected correlator has a significantly larger signal to noise ratio. This is related to the fact that using ρ conn we access the field strength and not its square; as a consequence ρ conn is much less sensitive to the fluctuations, which also means that it is expected to be more stable under smoothing of the gauge fields. In some cases one is interested precisely in fluctuation effects (e.g. in the study of the fluctuation-induced broadening of the flux tube) and it is then mandatory to use ρ disc , which in pure gauge theories can be precisely estimated using standard noise reduction techniques [61][62][63][64][65]; if this is not the case and large statistics are not available ρ conn is a more convenient choice [44,[66][67][68][69]. For these reasons and to directly compare with the recent results [44] we used the connected correlator ρ conn defined in Eq. (7).
The geometry of the connected correlator ρ conn is depicted in Fig. 1: the Schwinger line L is attached to the square Wilson loop W in the midpoint of its temporal extent, it reaches half the distance between the static color sources and then it moves x t lattice spacings in a direction orthogonal to the Wilson loop plane. Since ρ conn is a purely gluonic observable, the magnetic field can affect its value only by loop (sea) effects. In particular, to study different values of B, different sets of configurations have to be generated. When evaluating ρ conn it is important to realize that the external magnetic field breaks the lattice octahedral symmetry and Wilson loops oriented along different directions in general will not be equivalent; analogously, the two directions that are orthogonal to the plane of the Wilson loop will not be equivalent for generic magnetic field orientations (see later discussion).
In order to reduce the UV noise and improve the signal to noise ratio we adopted, as usual, smearing. With the aim of simplifying the comparison with previous results in the literature, we chose the same smearing procedure adopted in Ref. [44]: a single HYP smearing step [72] has been applied to all the temporal links, with parameters (α 1 , α 2 , α 3 ) = (1.0, 0.5, 0.5), then several APE smearing steps [73] have been applied to the spatial links, according to the definition where α AP E was fixed to 1/6 and S µ (x) is the sum of all the spatial staples associated to the spatial link U µ (x).
In order to investigate the profile of the flux tube we measured the longitudinal (i.e. directed along the flux tube) component of the chromoelectric field E l , since all previous studies showed this component of the field strength to be the dominant one. If we denote byμ the axis of the relative separation between the two color charges, the longitudinal chromoelectric field is given by where d is the distance between the charges, x t is the transverse distance at which the flux profile is probed and (t, µ) is the plaquette orientation.
Since smearing is used to reduce the UV noise, it is important to study the dependence of the results on the amount of smoothing adopted, in order to asses the reliability of the results. As an example, in Fig. 2 we report the values obtained for E l (d, x t ) using a 48 3 × 96 lattice, with a ≃ 0.0989 fm and d ≃ 0.7 fm: results are shown, as function of the number of APE smearing steps N AP E , for three different values of the transverse separation x t . It is clear that a non-trivial dependence on N APE is present and a prescription is needed to fix the value of ρ conn .
In Refs. [44,[67][68][69] the prescription adopted was to take the value at the plateau (or at the maximum) which is reached after some smearing steps; this is analogous to what has been done in the literature for similar quantities, like the gauge-invariant field strength correlators [15,70,71]. This prescription implies that the field value has to be taken after different numbers of smearing steps for different values of x t , since the plateau (or  (12)). maximum) is reached for larger values of N AP E when increasing x t , as can be seen in Fig. 2.
In this study we explore also a different prescription, in which all field strength values at different x t (and also at different values of the lattice spacing a) are measured keeping constant the smearing radius R s in physical units. The continuum limit is then taken at fixed R s and results obtained using different smoothing radii can be compared among them. Such a prescription is similar to the one adopted to compute renormalized observables using the gradient flow as a regulator (see e.g. [74]), and the similarity is even more striking in view of the practical equivalence between the smoothing techniques [75][76][77][78].
We fixed the value of the smearing radius R s in physical units according to the following relation (obtained in Ref. [77]) This expression slightly differs from the one reported in [77], since in the present work we adopt a normalization of the parameter entering the APE smearing that is different from the one used in the original derivation: A further difference is that in this study we use only the spatial staples in the APE smearing and we update only the spatial links, while in Ref. [77] a four dimensional smearing was studied. This would results in a further multiplicative factor in Eq. (12), independent of N AP E and a. Since our aim is to perform the continuum extrapolation at fixed R s , this multiplicative factor is practically irrelevant and in the following we will just use Eq. (12). In Fig. 3 we show some results for the flux tube profile extracted from simulations performed at lattice spacing a = 0.1249 fm (on a 40 4 lattice) using a Wilson loop of physical size d ≃ 0.7 fm: different symbols correspond to results obtained using different values of the smoothing radius R s and values extracted from the plateaux are also shown for comparison. From Fig. 3 it can be seen that not only the absolute scale of the flux tube depends on R s but also its shape. Similar behaviors are observed for all the values of a and d explored in this work.
Two different strategies can be adopted in order to keep the size of the Wilson loop constant in physical units while changing the lattice spacing. An approach consists in fixing a priori the extent of the Wilson loop in lattice units and the value of lattice spacing, imposing the constraint of constant physical size. A different possibility is to perform measures, for each value of the lattice spacing, using several Wilson loop extents, in order to be able to interpolate the results on a wide range of sizes. This second possibility requires some more care during the analysis, to check for possible systematics induced by the interpolation procedure; on the other hand it is much more flexible, since one can a posteriori decide the optimal size to be used in order to have small systematics and good signal to noise ratio. For this reason we adopted the second possibility. We verified that the interpolation in the smearing radius is stable (i.e. independent of the interpolating function adopted) to a high level of accuracy, as could have been expected given the fact that the dependence on R s is quite smooth. More care is needed for the interpolation in the distance d between the static color charges, since we have less data points available: for the coarsest lattice we measured Wilson loop from 4 × 4 up to 7 × 7, reaching 11 × 11 for the finest lattice. However, starting from d 0.6 fm, results are independent of the order of the spline interpolation adopted (linear, quadratic and cubic splines were tested). Using these interpolations we can compare the results obtained at different lattice spacings and, as an example, in Fig. 4 we show the outcome of this analysis for d ≃ 0.7 fm and R s ≃ 0.5 fm: finite cut-off corrections seem to be small and more visible for small values of x t ; similar considerations apply to results obtained for different values of R s and for the plateau method.
In order to perform the continuum limit of our results in a model independent way, data at all lattice spacings should be available for each value of x t . As one can see from Fig. 4, because of the discrete nature of the lattice distances this is not the case. In principle, we could think of adopting the same strategy used above and interpolate, for each lattice spacing, results obtained at different values of x t . However, in this case we find it more convenient to make use of a well known ansatz for the flux tube profile.
In particular, a parametrization that was shown in previous studies [44,[67][68][69] to well describe the data for the longitudinal component of the electric field is the Clem form where φ, α and µ are fit parameters and K 0 , K 1 are modified Bessel functions of the second kind. This parametrization of the longitudinal (chromo-)electric field is inspired by a similar parametrization of the longitudinal magnetic field around a vortex in type II superconductors, that was obtained in Ref. [79] by variationally improving an ansatz solution of the Ginzburg-Landau equations. We checked that Eq. (14) is consistent with the observed flux tube profiles for all the values of the lattice spacing, of the smoothing radius and of the quark-antiquark distance explored in the present work. For the purpose of performing the continuum limit, Eq. (14) is just a reasonably simple functional form that well describes data using three parameters. From a broader perspective, however, the fact that Eq. (14) well describes lattice data supports the dual superconductivity picture of confinement. In this picture condensation of the chromomagnetic degrees of freedom is expected to happen in the vacuum, confinement is the chromoelectic analogue of the standard Meissner effect and, in complete analogy with ordinary superconductors, the vacuum is characterized by two length scales: the penetration length λ and the coherence length ξ. These scales are related to the parameters entering Eq. (14) by the relations (see Ref. [79])  where κ = λ/ξ is the Ginzburg-Landau parameter, whose value discriminates between type I superconductors, corresponding to κ < 1/ √ 2, and type II superconductors, for which κ > 1/ √ 2 (see, e.g., Ref. [80]). To investigate the lattice spacing dependence of the results and to extract their continuum limit we used Eq. (14) with a-dependent parameters. Since in our lattice discretization the leading lattice artefacts are O(a 2 ), we performed a global fit to all the data corresponding to fixed values of d (quark-antiquark separation) and R s (smoothing radius), using the functional form in Eq. (14) with the substitutions Quantities denoted by the "0" subscript are the continuum values of φ, µ and α, while quantities with subscript "1" parametrize lattice artifacts.
The global fit works reasonably well for all explored values of R s as well as for the plateau method. For instance, for the data in Fig. 4 we obtain a value of the χ 2 /d.o.f. around 1.3; moreover the fit parameters are stable, within errors, when one eliminates data at the coarsest lattice spacing from the fit. The continuum values of the parameters φ, µ and α that are obtained in this way are shown in Figs. 5, 6 and 7 for three different values of the distance between the static color charges (d) and several values of the smoothing radius (R s ). In all cases a sizable dependence of the results on d and R s can be seen. Results computed at fixed smoothing radius converge, for large R s , to the values extracted using the plateau method, which displays only a weak dependence on the value of d used. Also the values of derived quantities like the Ginzburg-Landau parameter κ are dependent on the specific values of d and R s used, as shown  in Fig. 8. In principle, one could try an extrapolation of continuum results to zero smearing radiuds, R s = 0, however our present accuracy does not permit to perform that reliably.
For all the combinations of d and R s values studied we did not found results incompatible with κ 1/ √ 2: only for d ≃ 0.66 fm, using the plateau method or large smoothing radii, results for κ larger than this critical value are obtained, which are however compatible with it at one standard deviation. This behavior favors the interpretation of the QCD vacuum as a type I superconductor, in accordance with previous results in the literature [44], however the strong dependence of κ on the values of d and R s makes it very difficult to draw firm conclusions on the actual value of κ in QCD.  II: Equivalence classes of flux tubes, identified by the relative orientations of the magnetic field (always assumed to be directed alongẑ), the Wilson loop (in the plane (t, µ)) and the transverse directionρ, see Fig. 1.

IV. RESULTS IN A MAGNETIC BACKGROUND FIELD
The numerical results that will be presented in this section refer, unless otherwise explicitly stated, to the 48 3 ×96 lattice with lattice spacing a ≃ 0.0989 fm, that is the lattice size on which the largest number of simulations with B = 0 have been performed (see Table I).
The presence of the external magnetic field explicitly breaks some of the symmetries of the QCD Lagrangian, both in the continuum and on the lattice. Particularly important for our purpose is the breaking of the rotation symmetry: the details of the flux tube profile will depend on the relative orientations of the magnetic field, the Wilson loop and the transverse direction chosen to probe the chromoelectric field.
In this work the magnetic field will always be directed along one of the lattice axes, which we can assume to be theẑ direction. The plane of the Wilson loop is identified by the couple of indices (t, µ), with µ ∈ {x, y, z}, and we denote byρ the spatial direction, orthogonal to the Wilson loop, along which the chromoelectric field is evaluated (see Fig. 1). If the color sources are separated  Table II. along the magnetic field, i.e.μ =ẑ, the theory is invariant under rotations around theẑ axis and the two possible choicesρ =x andρ =ŷ of transverse direction are equivalent. If insteadμ =x, the orthogonal directionŝ ρ =ŷ andρ =ẑ are not equivalent, and an analogous situation happens forμ =ŷ. It is however simple to verify that the two combinationsμ =x,ρ =ŷ andμ =ŷ, ρ =x can be mapped into each other by using a rotation along theẑ axes and a reflection with respect to a plane containing theẑ axis, which are symmetry transformations also when a non-vanishing magnetic field is present. In a similar way it can be shown that the choiceμ =x, ρ =ẑ is equivalent toμ =ŷ,ρ =ẑ. We thus have, under the residual symmetry that is present when B = 0, three equivalence classes of flux tubes, that will be denoted by the shorthand L, TT, TL (T stands for transverse and L for longitudinal with respect to the magnetic field) and are reported in Table II for later reference. The flux tube shapes obtained in the three inequivalent geometries of Table II are shown in Fig. 9 for eB ≃ 3.12 GeV 2 (using R s ≃ 0.5 fm). The intensity of the longitudinal electric field when the flux tube is directed along the magnetic field is strongly reduced with respect to the case in which it is transverse to it. Moreover, in the transverse case, the flux tube loses its axial symmetry, as seen from the fact that the results for the cases TT and TL are different from each other. Notice that, in all cases, the symmetry under the transformation x t → −x t is preserved: in principle, one could expect asymmetries when the magnetic field is orthogonal to the quark-antiquark separation, however that would imply that the magnetic field induces a component of the chromoelectric field parallel to it, and this is protected by the CP symmetry, which is not broken by the magnetic background.
Also in the presence of a non-vanishing magnetic field the details of the flux tube strongly depend on the value of the smoothing radius R s . However the ratio of the chromoelectric fields with and without the magnetic field (i.e. E l (d, x t , R s , B)/E l (d, x t , R s , B = 0)) is remarkably insensitive to the value of R s . This is true for all the inequivalent classes of Table II and for all the values of the transverse distance x t studied, with some examples shown in Fig. 10. This means that we can study the effect of the magnetic field on the flux tube in an unambiguous way and for this reason the dependence on R s will be dropped in the following. Fig. 11 shows the changes in the flux tube induced by the magnetic field. The most striking effect that can be seen is the strong decrease of the chromoelectric field when the tube is collinear with the magnetic field (case L), while it slightly increases in the TT and TL cases. Two less prominent but still significant effects that are due to the magnetic field are the following: • in the longitudinal case the flux tube gets squeezed, since E(d, x t , B)/E(d, x t , B = 0) is a decreasing function of x t , • in the transverse case the flux tube loses its cylindrical symmetry, since the results for the cases TT and TL are not equal to each other.
This behavior is consistent with the general picture that emerges from previous studies of static potential and screening masses [14,16,17]: the magnetic field acts as transverse confinement catalyser and longitudinal confinement inhibitor. However previous studies (with the possible exception of Ref. [15]) investigated "integrated" quantities (like e.g. the static potential) and thus could not resolve the difference between the TT and TL cases.
In the remaining part of this section we will concentrate on the properties of the flux tube in the L case, i.e. the case in which the separation between the two color sources is collinear with the magnetic field. In this setup the cylindrical symmetry of the flux tube is preserved and it is reasonable to expect the Clem parametrization in Eq. (14) to well describe the numerical data, which indeed turned out to be the case.
The previously noted fact that the ratios E l (d, x t , R s , B)/E l (d, x t , R s , B = 0) are independent of the smoothing radius R s does not a priori implies the ratios of the fit parameters entering the Clem expression Eq. (14), like e.g. φ(d, R s , B)/φ(d, R s , B = 0), to be also independent of R s . This independence is however numerically observed to hold true with reasonable accuracy: deviations from a constant do not exceed the 5% for the range of parameters explored, an example being shown in Fig. 12. As a consequence, also the ratios of fit parameters computed with and without the magnetic field can be used to extract reliable informations on the effect of B on flux tubes.
To characterize the properties of the flux tube it is convenient to use, instead of the fit parameters φ, α, µ, some numerical parameters of more direct physical and geometrical interpretation. Two such parameters are the average square width of the flux tube w 2 and its energy density per unit length ǫ, defined by the expressions: If we assume for the longitudinal chromoelectric field the expression in Eq. (14), using known integrals of the modified Bessel functions (see, e.g., Eq. 5.52.1 and Eq. 5.54.2 in Ref. [81]) it is not difficult to prove the relations [44] that can be used to estimate w 2 and ǫ without having to numerically perform the integrals on the transverse directions. This is highly desirable due to the specific form of the integrands of Eq. (17): they are very small everywhere but for a sharp peak at intermediate values of x t and this makes the numerical integration unstable with the available numerical precision. The average square width of the flux tube w 2 is not strongly dependent on the magnetic field and slightly decreases by increasing B, being reduced by about 10% for the largest value of the magnetic field explored, eB ≃ 3.12 GeV 2 , see Fig. 13. This is consistent with the previously noted fact that (in the longitudinal case L) the flux tube gets squeezed by the magnetic field, see Fig. 11.
Since both the peak value and the square mean width of the flux tube get smaller in the presence of an external magnetic field, the energy density per unit length of the tube ǫ(B) will be a decreasing function of B. In a simple classical picture, the energy density per unit length of the flux tube is nothing but the string tension and it is thus interesting to compare the behavior of ǫ(B) with that of the string tension σ(B), that has been previously investigated in Ref. [16].
A direct comparison of ǫ(B) and σ(B) cannot however be performed, due to the dependence of ǫ(B) on the smoothing radius R s . This dependence disappears in the ratio ǫ(B)/ǫ(B = 0) and for this reason we show in Fig. 14  and σ(B)/σ(B = 0) (from [16]). Taking into account the fact that these two sets of data have completely different systematics and that in the determination of ǫ(B) there is also a theoretical bias (the form Eq. (14) of the flux tube was explicitly used in Eq. (18)), the agreement is reasonable. In future studies it will be highly desirable to have more precise data available, in order to estimate in an unbiased way the energy per unit length of the flux tube.

V. CONCLUSIONS
In this work we studied the dependence on a background magnetic field B of the chromoelectric flux tubes between two static color sources in QCD, using simulations performed with N f = 2 + 1 dynamical quarks of physical masses.
As a preliminary step we investigated the B = 0 case, pointing out that the smoothing procedure used to improve the signal to noise ratio can lead to significant systematics. We proposed, in alternative to the "plateau" method that is often used in the literature, the "fixed smoothing scale" method, in which the smoothing radius R s in physical units is kept fixed as the continuum is approached. Continuum results, however, still display a significant dependence on R s and a further extrapolation R s → 0 is desirable to obtain physically sensible results, but our data are not precise enough for this further extrapolation to be performed reliably. This makes it very difficult to provide firm results for interesting quantities like the value of the Ginzburg-Landau parameter κ in QCD.
Luckily enough, these problems are absent if one is just interested in studying the modifications of the flux tube induced by the external magnetic field B, as we showed in Sec. IV. Since for B = 0 rotational invariance is explicitly broken, three different cases have been studied, corresponding to the inequivalent relative orientations of the magnetic field, of the Wilson loop and of the transverse direction.
When the distance between the two static color charges is collinear with the direction of B, the two transverse directions are equivalent, the flux tube stays cylindrical and both the intensity of the chromoelectric field and the average square width of the tube decrease with the magnetic field. When the static charge separation and the magnetic field are perpendicular, the cylindrical symmetry of the flux tube is instead broken and the chromoelectric field inside the flux tube is a growing function of

B.
When the flux tube keeps its cylindrical symmetry the tube profile is still well described by the functional form Eq. (14), like in the B = 0 case. Using this fact we could estimate the energy density per unit length ǫ of the flux tube and we verified that, at least at a semi-quantitative level, the dependence of ǫ on B is consistent with the dependence of the string tension σ on B, as determined previously in Ref. [16].
The general picture that emerges from several studies [14,16,17], in which the magnetic field disfavours confinement in the longitudinal direction and enhances confinement in the transverse directions, is thus fully consistent with the results of the present work. Since the properties of the confining potential boil down to properties of the color flux tube, our results can in fact be seen as the microscopical origin of the macroscopic effects observed in previous works.
A natural extension of this work would be the study of other field components of the flux tube: while at vanishing magnetic field the longitudinal (chromo-)electric field is by far the dominant one, it is conceivable that at eB = 0 also other field components could be activated, making the field structure within the flux tube more complicated.