Improved thermodynamics of SU(2) gauge theory

In this work we present the results of our investigation about the thermodynamics of SU(2) gauge theory. We employ a Symanzik improved action to reduce strongly the discretisations effects, and we use the scaling relations to take into account the finite volume effects close to the critical temperature. We determine the beta-function for this particular theory and we use it in the determination of different thermodynamic observables. Finally we compare our results with previous works where only the standard Wilson action was considered. We confirm the relevance of using the improved action to access easily the correct continuum thermodynamics of the theory.


Introduction
Asymptotic freedom and confinement are two crucial properties of QCD. Confinement implies that the fundamental degrees of freedom of the theory, namely quarks and gluons, cannot be found as isolated particles in nature but exist only in complex bound states under normal conditions. The confinement properties of QCD-like theories are very well described in terms of a flux tube arising between quark-antiquark static charges. The vacuum quantum fluctuations of the flux tube are expected to be described by non-critical string models, and, for pure gauge theories without light quarks, string breaking does not occur and the accuracy of the predictions has been verified by many lattice Monte Carlo simulations. Despite the absence of dynamical quarks, the bound spectrum of pure gauge theories is still non-trivial due to the emergence of composite particles from the strong interactions between gluons, the so-called glueballs. There have been much effort in the past years to determine their properties at zero temperature, see Refs. [1,2,3], and also at non-zero temperature, see Refs. [4,5,6,7], although, an unambiguous experimental confirmation of their existence is still missing.
The fundamental nature of strong interactions is however quite different at very high temperature, where QCD behaves as a gas of free quarks and gluons due to asymptotic freedom. Understanding what happens to QCD for intermediate temperatures, near the deconfinement phase transition, is therefore the main reason for studying the thermodynamics of gauge theories. The chromodynamic flux tube is expected to survive below the critical temperature of quark deconfinement; various models have been developed to include thermal fluctuations to the QCD string, see for instance Refs. [8,9,10].
Quark-gluon plasma (QGP) is the phase of matter that can be probed experimentally by particle accelerators, such as RHIC and LHC, occurring at temperatures higher than ≈ 200 MeV. The properties of QGP even at quite large temperatures are compatible with those of a strongly interacting plasma that can be viewed as a perfect liquid, where color charges have long range interactions [11,12]. Because the success of the hydrodynamical description of high-energy heavy ion reactions, it is of great interest to compute the shear and bulk viscosities of the quark-gluon plasma. Because of the strongly interacting nature of the QGP, weak coupling perturbation theory is not able to capture the full thermal behavior of QCD. Lattice Monte Carlo simulations can provide a non-perturbative insight to the thermodynamics of the quark-gluon plasma, but still today it is not possible to compute the shear and bulk viscosities in full QCD and even in pure gluodynamics is an extremely complicated task, see Refs. [13,14,15].
The properties of pure gauge theories at non-zero temperature have been intensively investigated based on the idea that it is possible to describe the thermodynamics of N c = 3 QCD as a limiting case of a 1/N c expansion [16,17,18]. In particular, Feynman diagrams including quark lines give only a subleading contributions at large N c .
In the same line of research, there has been several predictions for the behavior of the quark-gluon plasma after the deconfinement phase transition in the context of the AdS/CFT correspondence. Interesting comparisons have been made in Ref. [16] with the so-called improved holographic QCD model, proposed in Ref. [19].
In Ref. [20] SU(N c ) gauge theories are instead compared with the quasi-particle ap-proach. It turns out that it gives a very good description of the interaction measure, i.e. of the trace anomaly, and of the thermodynamical quantities. All works, based on the lattice approach, which want to explore the large N c -limit, are based on simulations with N c ≥ 3. The reason is that, while the deconfinement transition for N c > 3 is a first order phase transition and for N c = 3 is a weak first order one, the case with N c = 2 is characterised by a second order phase transition. It is therefore expected that the models which describe the theories for N c ≥ 3 cannot describe the case with N c = 2 because this theory is qualitatively and quantitatively different.
Unfortunately, because simulations are missing, we do not know so much about its properties and we do not know how much SU(2) pure gauge is really different with respect to N c ≥ 3. The last systematic studies, concerning the thermodynamic properties of SU(2) gauge theory, go back to the beginning of the nineties, see e.g. Refs. [21,22] and Ref. [23] where only the energy density was considered. We believe that it is timely to perform such a study and this paper is devoted to this task.
The simulations in this work have been done using the Symanzik improved action with 6-links plaquette. This is an important aspect of our work. It is well known, see Ref. [24], that the standard Wilson action, with temporal extension N τ = 4, leads to almost 50% of corrections due to finite cut-off effects, but using our improved action this is reduced below 2%. We do not need therefore to simulate the theory with high values of N τ , making the entire work much cheaper. For example, in Ref. [25] the authors simulated SU(3) gauge theory with N τ in the range [5,8], which is by far much more expensive. From a numerical point of view, the main difficulties come from the finite volume effects close to the deconfinement transition which is due to the second order phase transition.
Some thermodynamic quantities need the evaluation of the β-function; this task is performed in Section 2. In Section 3, we measure a number of thermodynamic observables and we plot them. We compare our results with other works in Section 4 and finally we draw our conclusions in Section 5.

Scale setting and the determination of the β-function
The determination of the Callan-Symanzik β-function is relevant to set the scale, i.e. the physical temperature realised in our simulations. It is also important to determine how the physical volumes of our lattice changes when the bare coupling g changes, that is the starting point to derive the trace of the energy-momentum tensor at non-zero temperature. While there are different papers where the β-function has been determined for the case of the standard Wilson action, see e.g. Refs. [23,26,27], there are no works, to best of our knowledge, in the case of the action used in this work, therefore we must proceed to a separate calculation.
The physical measure of the lattice spacing a as a function of the inverse coupling β = 4/g 2 is performed by determining how a specified observable depends on β. In literature different observables have been considered, previous works include the plaquette, see Refs. [18,28], and the string tension, see Ref. [29]. In this work we consider three different observables: the critical temperature T c of the deconfinement phase transition, given by the critical coupling β c for different values of the lattice extent in the temporal direction N τ , the scale parameter w 0 , see Ref. [30], and the scale parameter t 0 , see Ref. [31].
The β-function has been fitted starting from the expected scaling of physical observables near the continuum limit 1 : . ( The running of the lattice spacing a as a function of g is provided by the scaling function F (g 2 ) up to corrections A(g 2 ), which takes into account the lattice artifacts [32] The scaling function F (g 2 ) is given by the product of two terms: the first one is the result of the integration of the two-loop scheme-independent weak expansion of the β-function while the term λ(g 2 ) takes into account the terms of higher order of perturbation theory. This term has been parametrised in different ways in literature; we have considered two of them. The first one, see Ref. [23]: and the second one, see Ref. [32]: We have verified that the second method gives worst results, in particular when the correction A(g 2 ) is considered. Therefore in the following we will show only the results obtained with the functional form λ(g 2 ).
The function A(g 2 ) accounts for scaling violations far away from the continuum limit driven by the running of irrelevant lattice operators at non zero lattice spacing a, see Refs. [32,33]. The form of A(g 2 ), is specified in terms of two even integer numbers ν and ν ′ , because we require that a is an even function of g. The term containing X n,ν takes into account the leading correction in a; the term Y n ′ ,ν ′ the next-to-leading one. Each term has been normalised so that X n,ν and Y n ′ ,ν ′ describe the fractional amount of scaling correction at a standard value of g = √ 2, corresponding to β = 4/g 2 = 2. The β-function β f can be expressed as a function of β = 4/g 2 : where the term ∂β/∂ log(a) can be easily determined using Eqs. (1), (2), (3): Starting from these relations we can determine three different definitions of the lattice β-function; in the following sections we present each definition and the resulting scale in detail.

Fitting the critical β
The β-function can be determined from the value of the critical bare gauge coupling g where the deconfinement phase transition occurs at different values of the lattice temporal extension N τ . To this end, we can exploit existing calculations presented in literature, see Ref. [34] and summarised in Table 1.
The points (β c , N τ ) are fitted using Eq. (1), where the lattice spacing is given by Eq. (4) and Ω ′ is the dimensionless ratio Λ Lat /T c . Note that for a given integer value of N τ , the deconfinement phase transition is located by looking for the maximum of the Polyakov loop susceptibilities as a function of β and not vice-versa. Therefore in this fitting procedure the error appears to be on the abscissa, 1.59624(13) 4 1.699(1) 6 1.8287(11) 7 1.8747(30) 8 1.920(5) Table 1: Critical value of β c = 4/g 2 for different values of N τ from Ref. [34]. Note that the value at N τ = 5 has not been considered because it was clearly too far from the interpolating function, perhaps the error cited in Ref. [34] has been underestimated.
on β, and there is no error on the ordinate, on N τ . In general, when in this work we have to combine errors in both dimensions, let us say σ x and σ y , we use the approach explained in Ref. [35], i.e. we combine the two errors as s 2 σ 2 x + σ 2 y , where s is the slope of the curve in the point we are considering. In Table 2 we report the value of χ 2 /d.o.f. for different fits of the data presented in Table 1.
According to this table, the best two fits, with the smaller χ 2 /d.o.f., are those labelled with "fit2" and "fit3". The final result is presented in Figure 1 where the error is the statistical one.
In Figure 2 we plot the ratio T /T c versus β for different values of N τ as determined from data of Figure 1. Here the value is given as average of the previous two best fits and the error takes into account the statistical error of the two fits and the systematic error given by the difference between the two curves. The values presented in this plot are the ones used in the remain part of this work every time we need a correspondence between β and T /T c .

Fitting the scale parameter w 0 and t 0
The numerical integration of the flow equations, defined from the functional derivative of the Symanzik gauge action with respect to the gauge-link variables, is performed using fourth order Runge-Kutta integrator with a discretisation of the flow time equal to δt = 0.01. We have measured the energy density, defined to be equal to the traceless antihermitian part of the clover plaquette, every ten integration steps. The scale parameter w 0 is defined as the square root of the flow time t, solution of the implicit equation We have used two different values for the reference value u (0.2 and 0.3). The observable t 0 , defined as the flow time t where the equation is fulfilled, is expected to be affected by larger discretisation effects, see Ref. [30], that will appear as scaling violations in our fitting procedure. To compute the logarithmic derivative in Eq. (15), we have performed a polynomial fit of the expectation value of the flow observable E(t) . Given that the flowed energy density is strongly correlated for flow times close to each other, we estimate the statistical error using the bootstrapping method, performing therefore a fit on each bootstrapping sample. In our tables we quote only the statistical error, without the systematic error coming from various possible fitting intervals and degrees of the interpolating polynomial. The values of the scale w 0 (β)/a and t 0 (β)/a 2 , that we have measured, can be found in Table 3. In the last column appears also the value of the "residual" non-zero temperature of the system determined by the ratio of the critical N τ plotted in Figure 1 and the size of the four-dimensional hypercubeL.
Note that, for β = 1.8, β = 1.825 and β = 1.85, we measured the scales for three different volumes, up to T /T c ≈ 0.36, and the difference never exceeds three standard deviations. Since T /T c ≈ 0.36 is also the maximum "residual" temperature corresponding to the larger value of β that we have used in our simulations, i.e. β = 2.025, we can safely assume that the finite volume effects are under control in the entire range of β. Moreover, as we will show in Sec. 3.1, within the precision of our measurements, with volumes ranging fromL = 24 toL = 56, and with a range of β corresponding to a "residual" non-zero temperature below T /T c ≈ 0.45, also the spatial plaquette is not affected by finite volume effects.
The measured data (β, w 0 (β)/a) are fitted using Eq. (2), where the lattice spacing is given by Eq. (4) and Ω = w 0 Λ Lat : Similarly, using Eq. (3), we fit the scale t 0 (β)/a 2 aŝ  (13) 3.099 (14) 7.877(68) 11.074(63) 32 0.36 Table 3: Summary of the values ofŵ 0 andt 0 used to determine the β-function. The quoted error is only statistical and does not include systematic errors arising from different choices of the fit of the flow.L is the size of the hypercube used. T /T c is the "residual" non-zero temperature of the system determined from Figure 1.
However we have to note that the value of χ 2 /d.o.f. is not of the order of one, but the best we could get is χ 2 /d.o.f. ∼ 5. Such a large deviation from the asymptotic scaling is a typical situation that occurs when fits of non-perturbative results, using lattice perturbation theory, are considered. The χ 2 /d.o.f. could be improved if more fitting parameters were included, but in our case, given the non-linear form of our fitting function, such a method provides an unstable minimum of χ 2 . An another reason for the large χ 2 /d.o.f. is that the integrated energy of the flow is estimated rather precisely, so that the statistical error is comparable to the systematic error coming from finite volume effects and the various possible fitting range and fitting polynomials of the flow observable E(t) . However, we are not concerned of a such large χ 2 /d.o.f., since the largest systematic errors are coming rather from scaling violations as lattice artefacts, i.e. when different observables are used to set the scale. Our final β-function is a combination of three different definitions, see Sec. 2.3, and the final systematic error is large enough to accommodate any possible mismatch of our fits from the pertubative scaling.
In Figure 3 we plot the scale w 0 as given in Table 3 together with the fits labelled as "fit2" and "fit6" in Table 4; overall the quality of the fit looks pretty decent. In Figure 4 we plot instead the scale t 0 of Table 3 together with the best two fits labelled as "fit9" and "fit10" in Table 5. As final result we consider the average of the two values, coming from the two fits, and as error the sum of the two statistical errors and a systematic one which comes from the difference between the two values.

Final results of the lattice β-function
In Figure 5 we plot the determination of ∂β/∂ log(a) and in Figure 6 that of the βfunction for our three observables. The perturbative dashed line present in the plots has been determined considering only Eq. (6). We have three different results which are not compatible with each other at low β, due to discretisation effects, which increase when β decreases, and that are not universal, see Ref. [32].
In this work the β-function and its error are safely defined by considering a combined final uncertainty that will arise mainly from the difference of the three possible observables used to set the scale. The final results is presented in Figure 7 and in Figure 8. In any case, it is worth to consider what is the impact of such large discrepancy induced by violation of the scaling of the physical observables at low β. From Figure 6, the stronger difference between the various β-functions appears for β 1.8, that corresponds to T /T c 1 for N τ = 5 (see Figure 2), i.e. the largest uncertainties of the β-function are in the confined phase, where thermodynamical quantities, such as pressure and energy density, are usually very small.
A method to avoid at least part of the previous uncertainties in the calculation of thermodynamical quantities must provide a direct definition of the energy-momentum tensor and of its renormalisation on the lattice, some work in this direction has been presented for instance in Ref. [37], based on the gradient flow, or in Ref. [38], based on a formulation of the thermal theory in a moving reference frame. In any case, discretisation errors are unavoidable in any lattice numerical simulations and will appear in the determination of the equation of state both in the ordinate for renormalised quantities, and in the abscissa as uncertainties in the definition of the physical temperature. As we will show in the following sections, the use of the Symanzik improved action is crucial to suppress lattice discretisation errors without requiring at the same time a demanding computational cost.

Thermodynamics
The thermodynamic quantities we are interested in are the pressure p, the energy density ǫ, the trace anomaly ∆ and the entropy density s. They are defined from the partition function according to the relations being f the Helmholtz free energy. The pressure is determined using the integral method, see Ref. [22]: where we have introduced the action density S = (T /V ) S . Note the existence of a reference β 0 which should correspond to a sufficiently small temperature where the pressure can be safely assumed to be zero (for a different approach see Ref. [39]). The trace anomaly is given by The other two quantities are determined as a function of the pressure and of the trace anomaly: The Stephan Boltzmann (SB) limit for these quantities is known to be equal to We will use these relations to normalise our final results. Note that Eqs. (29)-(32) are not taking into account the correction to the canonical partition function which is proportional to ln (N s /N τ ), see Refs. [40,41].

Simulations
Two different values of the temporal extension of the lattice, N τ = 4 and N τ = 5, have been used for the simulations at non-zero temperature, while we have fixed an aspect ratio of N s /N τ = 6, see Table 6. The remarkable result is that the second lattice, with temporal extension N τ = 5, turned out to be already close to the continuum limit, with a small correction with respect to the results coming from the lattice with N τ = 4.
For N τ = 4, we have measured the action density and the Polyakov loop in the interval 1.550 < β < 2.165 and for N τ = 5 in the interval 1.655 < β < 2.330; in both cases the measurements were done every ∆β = 0.005.
Note that finite volume effects depend on the ratioξ/N s , whereξ is the correlation length: far from the critical β this ratio goes to zero and there are small finite volume effects. On the contrary, close to the deconfinement phase transition, the correlation length will diverge for a second order phase transition, as in the case of the SU(2) Yang-Mills theory [42]. Therefore one can set a decreasing aspect ratio increasing the distance from the critical β and our value of N s /N τ = 6 is pretty arbitrary, tuned to control finite volume effects near the critical temperature.
In Table 7 we show the details of our simulations at zero temperature, generated in order to perform the subtraction of the zero temperature expectation value of the action density. We show also the "residual" non-zero temperature in each case. Its relevance and that of the finite volume effects can be seen in Figure 9 where we have compared the value of the spatial plaquette for different volumes along the entire range of β where we have used our data. We have plotted the ratio: where P (N s ) is the value of the spatial plaquette measured at the spatial volume N 4 s and ∆P is its error.Ñ s labels the value with respect to we are comparing the data and it is fixed in each single plot. From this figure it is clear that the "residual" temperature and the finite volume effects are always smaller than the statistical fluctuation of our measurements.   Table 7: Summary of the simulations used at zero temperature. The interval used to span the β interval is always ∆β = 0.005. The range in T /T c is the "residual" non-zero temperature of the system determined from Figure 1.
The error in the determination of our observables depends on the error on the action density N 4 τ [ S 0 − S T ]. It is possible to show that, for a noninteracting theory, this error is proportional to N 3.5 τ /N 1.5 s . This relation explains why the number of configurations that must be used in order to get a reasonably small statistical error increases hugely moving toward the continuum limit N τ → ∞. In our case going from N τ = 4 to N τ = 5 required already an increase by roughly a factor five of the computational cost.

Action density and finite size scaling
Since the SU(2) gauge theory is characterised by a second order phase transition, see Ref. [42], strong finite volume effects are present around the critical temperature. It is therefore necessary to simulate different volumes to extrapolate to the infinite volume limit. We have several ensembles with five different volumes both at N τ = 4, with volumes (16,20,24,28,32), and at N τ = 5 with volumes (20,25,30,35,40), see Table 6.
Close to the critical temperature, the infinite volume limit has been extrapolated using the finite-scaling approach. The action density S is a lattice operator which, under Svetitsky-Yaffe conjecture, is mapped into the energy operator of a statistical model, as shown in Ref. [43]. The scaling behaviour for S is given by: where L is the spatial extension and Q S is the scaling function for S . At t = 0 we can therefore extrapolate the action density to the infinite volume limit following the ansatz (see also Refs. [44,45]): The critical indexes for SU (2) in 4d are those of the Ising model in 3d (see section 3.2.1 of Ref. [42]): We have verified that the action density S is affected by finite size effects, non compatible with the statistical errors, in the interval T /T c ∈ [0.985, 1.005]. Since Eq. (35) is valid only at the critical point, we have tried to expand perturbatively Eq. (34) for t = 0, see Refs. [46,47,48,49], to extrapolate the results to infinite volume. Unfortunately, we have only results for five volumes that are not enough to allow a stable and reliable numerical extrapolation. Therefore, we follow the ansatz of Eq. (35) in the entire critical region T /T c ∈ [0.985, 1.005]. Anyway, to take into account the systematic error due to the sloppy infinite volume extrapolation, we verified the scaling of this quantity, plotting ( S L − S ∞ )L 1/ν−d vs tL 1/ν . Because the five curves were not compatible with each other, we increased arbitrary, in the critical region, the statistical error of S ∞ until the five curves were made compatible. At the end we tripled the statistical error at N τ = 4 and we doubled at N τ = 5. Thanks to this procedure, the value of S ∞ and its error are determined in a way which should be able to correctly estimate the presence of the systematic error.
The final value of the action density, normalised to the T = 0 value, is plotted in Figure 10. It is interesting to note that the results obtained at N τ = 4 and at N τ = 5, which correspond to a smaller lattice spacing, are compatible. Discretisation errors smaller than all the other errors are possible because we are using an improved action which brings our results already close to the continuum limit. We do not need therefore in our analysis to introduce any correction term R I (N τ ), as done for example in Ref. [16]. As a matter of fact, as discussed in Ref. [50], the expected correction for our action is, at N τ = 4, of the order of 1.35% which is smaller than our statistical error.
It is interesting to look at the scaling of the susceptibility of the Polyakov loop, which is given by a scaling function Q χ without a constant term: We show the results in Figure 11 and in Figure 12 (see for comparison Ref. [45]). Clearly the susceptibility follows the scaling relation in a very wide range of t.

Thermodynamic results
Using the relations introduced at the beginning of Sec. 3, the action density plotted in Figure 10 and the derivative of the β-function of Figure 7, we can now determine all the other thermodynamic observables. The pressure, normalised to its SB value, is plotted in Figure 13. In Figure 14, we plot the trace anomaly normalised to the SB value of the pressure (as has been done in Ref. [16]). The SB normalised energy and entropy densities can be found respectively in Figure 15 and in Figure 16.
As can be seen all our observables reach a value around 90% of the SB limit at T /T c = 5 and the results obtained at N τ = 4 and N τ = 5 are compatible, confirming that the discretisation effects are under control.

Comparison with other works
It is interesting to compare our results, in the deconfined phase, with those of Ref. [16] where results for SU(N c ), and N c ≥ 3 have been considered. Note that in that work was used only the standard Wilson action and only one volume, therefore we expect that both discretisation and finite volume effects are present. All thermodynamical observables we have considered reach the SB limit quicker than in Yang-Mills theories with N c ≥ 3. For example, in Figure 17, the value of the pressure at T /T c = 3.0 is ∼ 10% higher. The difference can be better appreciated comparing directly the trace anomaly, see Figure 18. In this case a huge difference appears around 1.5T c and the value is always above the N c ≥ 3 case.
Moreover, we can compare our results, in the confined phase, with those published in Refs. [51,52] where they simulate SU(2) pure gauge theory but, also in this case, using the standard Wilson action at fixed volume. However, they simulate different values of N τ , ranging from 5 to 10. In Figure 19 we compare the pressure with their continuum extrapolated results. Our values are always smaller and, close to the critical point, the value is about 1.4 times smaller. This difference can be explained by finite volume effects, that, as we have already seen, have a strong effect in this theory.
The results of the pressure can be affected by the choice of the reference point β 0 , see Eq. (25), therefore we compare also the trace anomaly directly with our results in Figure 20. Also in this case we can see a clear difference in the range 0.9 T /T c 1.0. In this case, the discrepancy could be given by the different choice of the β-function.
The observation that the trace anomaly falls off as 1/T 2 above T c lead to the development of many phenomenological models, see Refs. [53,54,55]. It is therefore interesting to compare SU(2) with previous studies where N c ≥ 3 was considered.
In Figure 21 we compare our results with Ref. [16]. We plot the quantity ∆/T 2 versus (T c /T ) 2 to see whether there exists a region with a linear behaviour. The figure suggests that SU (2) is compatible with the other theories only for temperature above ≈ 2 T /T c . Otherwise, for temperature from T /T c until 2T /T c the values for SU(2) are larger and not compatible with the others within the errors. The difference could be guessed already clearly from Figure 18. A real difference between SU(2) and SU(N c ) Yang-Mills theory could be claimed only after all systematic errors are carefully taken into account. Residual finite volume effects in our results, different β-functions, missing of continuum limit in the results of Ref. [16] could be at the origin of the discrepancy that we observe. However, in Figure 22 we plot the quantity ∆/(T 2 T 2 c d A ) (d A = N 2 c − 1) and our results show a better compatibility with those of Ref. [18], where both the thermodynamic and to the continuum limit have been extrapolated. We can therefore state that our results do not exclude the possibility that also for SU(2) the trace anomaly has a 1/T 2 behaviour, even if further simulations are necessary.
It is clear, by the examples considered, how much the use of an improved action is important to study the thermodynamic properties of a QCD-like theory, in particular when analytical determinations are compared to lattice results, as for the holographic model in Ref. [16], for the effective bosonic string model in Ref. [51], or for the hadronresonance-gas model in Ref. [52].

Conclusions
In this paper we have presented our results concerning the thermodynamics of SU(2) pure gauge theory. This is the first work, after almost twenty years, where a systematic study of the equation of state of this theory has been performed. The SU(2) Yang-Mills theory can still be useful to compare and test some interesting models, which go from effective string descriptions to large N c -limit results and from holografic models to quasiparticles descriptions. For our simulations we have used a Symanzik improved action so that our results, already at N τ = 5, are compatible with the continuum limit within the statistical errors.
We have performed many simulations on different volumes near the deconfinement transition, to control finite volume effects that are significant for a theory with a second order phase transition. We extrapolated our results to the thermodynamic limit following the finite size-scaling relations. We have determined non-perturbatively, employing three different methods, the β-function and later we have determined the main thermodynamic observables for the equation of state.
Finally we have compared our results, both in the confined and deconfined phase, with previous works, where clearly emerge the importance of using an improved action in this kind of measurements, in particular when we want to compare lattice results with analytic models.  Figure 1: Interpolation of N τ versus β using "fit2" and "fit3" of    Figure 3: Interpolation of w 0 versus β using "fit2" and "fit6" of Table 4.             Figure 9: In this plot we show the finite volume effects of our simulations at zero temperature. We plot the quantity [P (N s ) − P (Ñ s )]/∆P (Ñ s ), where P (N s ) is the value of the spatial plaquette measured at the spatial volume N 4 s and ∆P is its error. Comparing the results from different volumes, in the same range of β, we can see that the difference is always smaller than three standard deviations.               Figure 22: Comparison of the trace anomaly in the deconfined phase with Ref. [18].