Scaling functions of the three-dimensional $Z(2)$, $O(2)$ and $O(4)$ models and their finite size dependence in an external field

We analyze scaling functions in the $3$-$d$, $Z(2)$, $O(2)$ and $O(4)$ universality classes and their finite size dependence using Monte Carlo simulations of improved $\phi^4$ models. Results for the scaling functions are fitted to the Widom-Griffiths form, using a parametrization also used in analytic calculations. We find good agreement on the level of scaling functions and the location of maxima in the universal part of susceptibilities. We also find that an earlier parametrization of the $O(4)$ scaling function, using 14 parameters, is well reproduced when using the Widom-Griffiths form with only three parameters. We furthermore show that finite size corrections to the scaling functions are distinctively different in the $Z(2)$ and $O(N)$ universality classes and determine the volume dependence of the peak locations in order parameter and mixed susceptibilities.


I. INTRODUCTION
Universal critical behavior in the 3-d, Z(2) and O(N ) universality classes plays an important role in the analysis of phase transitions in many statistical models as well as quantum field theories.In the fundamental theory of strong interactions, Quantum Chromodynamics (QCD), phase transitions that occur at finite temperature and vanishing as well as non-vanishing conserved charge chemical potentials belong to these universality classes.The spontaneous breaking of chiral symmetry in QCD is expected to exhibit universal critical behavior in the 3-d, O(4) universality class [1], and the 3-d, Z(2) universality class is expected to describe critical behavior at the so-called critical endpoint, a yet to be discovered second order phase transition that is expected to occur in QCD with non-zero quark mass values and non-vanishing baryon chemical potential.A second order phase transition in the Z(2) universality class also occurs in QCD at non-vanishing, imaginary values of chemical potentials at the so-called Roberge-Weiss endpoint [2].Also, the O(2) universality class plays a role in the studies of the phase diagram of QCD, as many calculations are performed in a discretized version of the theory, using so-called staggered fermions, in which only this smaller symmetry is realized.
Numerical studies of the phase structure of statistical models, and, in particular, of complicated theories such as QCD, are being performed on finite lattices.A good understanding of finite-size effects, thus, is generally of importance.In the limit of small external symmetry breaking fields and large volumes also these finitesize effects are universal, i.e., characteristic for a given universality class.For this purpose a powerful renormalization group framework has been developed in statistical physics which leads to a detailed finite-size scaling theory for critical behavior [3,4].This framework has been used to analyze finite-size scaling behavior of systems in the 3-d, Z(2) and O(N ) universality classes.The finite-size dependence of thermodynamic observables in 3-d, O(N ) spin models has been examined using Monte Carlo simulations [5][6][7], and finite-size scaling functions have been derived using the functional renormalization group approach [8].For the O(4) universality class, we provided an updated parametrization for the infinite volume scaling functions [9] and presented a parametrization of the finite-size scaling functions f G (z, z L ) and f χ (z, z L ) [10], which describe finite volume corrections to the singular behavior of the order parameter and its susceptibility.In this work, we will extend these studies and provide finitesize scaling functions also for the Z(2) and O(2) universality classes, by performing Monte Carlo simulations with improved Hamiltonians [11][12][13], which have been constructed to suppress contributions from correctionsto-scaling and, thus, allow for easier access to the desired scaling functions.We furthermore present a parametrization of the infinite volume scaling functions, determined from Monte Carlo simulations, using the Widom-Griffiths (WG) form [14][15][16][17] of these scaling functions.In the Z(2) [18] and O(2) [7] universality classes, the relevant parameters entering this analytic form have been determined previously using -expansion and other field theoretic methods applied directly in 3-d.This paper is organized as follows.In the next section we introduce the Z(2) and O(N ) models for which we present new Monte Carlo results and define the basic observables studied by us.In Section III we introduce basic relations for finite-size scaling functions.Section IV is devoted to the determination of the infinite volume scaling functions for the Z(2), O (2) and O(4) models, using a parametrization based on the Widom-Griffiths form.Here, we also determine the non-universal parameters for the improved Z(2) and O(2) models that are needed to introduce the scaling variables z and z L .In Section V we present our results for the finite-size scaling functions.We give our conclusions in Section VI.In Appendix A, we discuss the determination of the two non-universal scales H 0 and L 0 and in Appendix B, we give explicit expressions for the expansion coefficients d − 1 and d − 2 appearing in the scaling function f G (z, z L = 0) at asymptotically large, negative arguments.

II. LATTICE SETUP AND OBSERVABLES
We discuss here universal scaling properties for 3dimensional, N -component φ 4 models, i.e., spin models in the 3-d, Z(2) (N = 1, Ising model), O(2) (XY model) and O(4) universality classes described by the Hamiltonian, with Φ x ≡ φ x,1 for the 3-d, Z(2) spin model, (x, y) denoting nearest neighbor sites on the lattice, and Φ x ≡ (φ x,1 , ..., φ x,N ) for the 3-d, O(N ) spin models.For specific choices of λ, the above Hamiltonian is called "improved" since the quartic coupling λ appearing in the potential term of the spin models has been optimized to reduce the effect of contributions from sub-leading relevant scaling variables to universal scaling behavior of these models [11].We use the parameters λ = 1.1 [12] in the case of Z(2) and λ = 2.1 for O(2) spin models [13], respectively.In the O(4) case, we use the standard, unimproved Hamiltonian, corresponding to λ = ∞.The temperature T is defined as the inverse of the coupling J, i.e., T ≡ 1/J, and the external field coupling H controls explicit symmetry breaking in the Hamiltonian.We introduce this symmetry breaking term such that it couples only to the first component of the spin variable Φ x , defined on the sites x of a three dimensional lattice of size L 3 .Using the Hamiltonian introduced in Eq. 1 the partition functions of the 3-d, Z(2) and O(N ) models are given by, From this, one obtains the free energy density in units of temperature T , f (T, H, L) = −L −3 ln Z(T, H, L).The derivative of the free energy density with respect to the external field H defines the order parameter, M , for spontaneous symmetry breaking,  (18) 0.97917(55) 0.7686 H0 0.79522 (17) 1.36632 (28) 4.845(66) t0 0.303376(45) 0.4540 (11) 1.023 (16) Table I.Critical exponents in the 3-d, Z(2), O(2) and O(4) universality classes and non-universal parameters for the improved Hamiltonians used in our simulations.Z(2) critical exponents are taken from Zinn-Justin [18] and the O(2) values are taken from [19].Exponents used for the O(4) case are taken from [20].Other critical exponents are obtained using the hyper-scaling relations.The critical temperature (Tc = 1/Jc) of the Z(2) model with λ = 1.1 is taken from [12] and Tc for the O(2) model with λ = 2.1 is taken from [7].All other non-universal parameters have been obtained in this study.For the Z(2) model, we find results for the scales t0 and H0 that are in good agreement with previous results obtained in [21]; t0 agrees to better than 1%, and the result for H0 is smaller by about 2%.In the O(4) case, we give critical exponents and non-universal parameters used also in a previous analysis of scaling functions [9] (see the text for further references).
The (longitudinal) susceptibility χ h and the mixed susceptibility χ t are obtained as derivatives of the order parameter with respect to H and J, respectively, In the absence of explicit symmetry breaking (H = 0), the Z(2) and O(N ) spin models undergo second order phase transitions at critical temperatures T c ≡ 1/J c .The critical temperatures of the 3-d improved Z(2) [12], O(2) [7] and unimproved O(4) [20] spin models, with couplings λ, as introduced above, are well determined.We give the critical temperatures together with other universal and non-universal model parameters in Table I.
In the case of the O(4) spin model, we do not perform new MC calculations, but re-parametrize results for scaling functions already obtained in [9].We therefore give in Table I the parameters actually used in that calculation.They are consistent with analytic results [18] but differ somewhat from recent MC results [22].
For H = 0 as well as for finite lattice sizes L < ∞, pseudo-critical temperatures, T pc,o (H, L), with o = h or t, can be defined as locations of maxima in the susceptibilities χ h and χ t .
Monte Carlo simulations have been performed by us for the improved Z(2) and O(2) models.For our calculations we use a code, which has been developed and used previously in simulations of Z(2) and O(2) models1 .The statistics collected in calculations with the Z(2) and O(2) models on different size lattices is given in Tables II-IV. .These data sets were used for the determination of the scale parameters H0 and L0, discussed in Appendix A. The data sets generated on the L = 96 lattice were also used in finite-size fits discussed in Section V.
-480000 184000 120000 Table IV.Number of configurations generated per parameter set (J, H) on lattices of size L 3 in the region z < −2.These data have been primarily used for the determination of the scale parameter t0 obtained together with the determination of the infinite volume scaling functions discussed in Section IV.

III. SCALING FUNCTIONS
In order to analyze universal critical behavior in the vicinity of the second order phase transitions occurring in the 3-d, Z(2) and O(N ) spin models, the free energy is split in a singular and regular contribution, respectively, The scaling behavior of, e.g., the order parameter M and the susceptibilities χ h and χ t is derived from the renormalization group analysis of the singular part of the free energy where b is a free scale parameter and y t and y h are two relevant critical exponents2 .In Eq. 7, we introduced the reduced temperature (t), external field (h) and finite volume (l) scaling variables, They are normalized by non-universal scale parameters t 0 , H 0 and L 0 , respectively.The exponents y t and y h define the two independent critical exponents of the universality class under consideration, Here β, δ and ν are critical exponents which are related to each other through the hyper-scaling relation δ = dν/β − 1.In our current analysis we use results for the exponents β and δ as basic input.These critical exponents are well determined for the 3-d, Z(2) and O(N ) universality classes.We use here the Z(2) results obtained in [18] and the O(2) values from [19].They are given in Table I.In the O(4) case, we use critical exponents and non-universal parameters that have also been used in a previous analysis of scaling functions [9].Choosing the scale parameter b = h −1/y h , we obtain for the free energy density (10) where we have introduced the finite-size scaling function with arguments (z, z L ) defined as Using Eqs. 3 and 10 we obtain the order parameter M , and its susceptibilities with scaling functions f G (z, z L ), f G (z, z L ), and f χ (z, z L ) defined, respectively, as The finite-size scaling functions can be determined in the vicinity of the critical point (t, h, l) = (0, 0, 0), where regular contributions to the order parameter and its susceptibilities, given in Eqs.13−15, are negligible3 , The non-universal scale parameters t 0 and H 0 are fixed by the following conditions on the order parameter at infinite volume or, equivalently, in terms of the scaling function, The scale L 0 is obtained using a normalization condition for the finite-size scaling function f G (0, z L ).We define z L = 1 as the point at which the order parameter, evaluated at T c , is 30% smaller than its infinite volume value, i.e.
This differs from the choice used in [10] but has the advantage of allowing better comparison of finite-size scaling functions obtained in different universality classes.
In the following two sections, we will discuss results for the infinite and finite volume scaling functions, respectively.In order to judge which lattice sizes and external field parameters are needed to get close to the infinite volume, universal scaling regime, we first analyzed the z Ldependence of the scaling functions at some fixed values of z on different size lattices.The three scaling functions f G (z, z L ), f G (z, z L ), and f χ (z, z L ) have been calculated at a few values of z as functions of z L .Using z and z L as variables, of course, does require the determination of the non-universal scales (t 0 , H 0 , L 0 ) which we are going to discuss in the next section and in Appendix A.
Results from the calculations of the scaling functions for some fixed values of z, performed on lattices of size L 3 with L = 48 and 96, are shown in Fig. 1.As can be seen, the finite-size effects in all three scaling functions are almost negligible for z L < 0.4.This is consistent with findings obtained in calculations with the standard O(4) model [10] and can also be concluded from Fig. 8, shown in Appendix A, where we compare results for f G (z = 0, z L ) in different universality classes.
In the following, we thus use our numerical results for z L < 0.4 as approximation for infinite volume limit results.

IV. INFINITE VOLUME SCALING FUNCTIONS
In our discussion of scaling functions in the infinite volume limit, z L = 0, we suppress the second argument of the scaling functions, i.e. we introduce f G (z) ≡ f G (z, 0) and similarly for f G (z) and f χ (z).These scaling functions have been determined previously using -expansions [25] and perturbative field theoretic approaches applied directly in 3-dimensions [7,18,26], as well as in Monte Carlo (MC) simulations [5,6,21].For the Z(2) universality class it has been shown that the scaling functions, obtained in MC calculations, are in good agreement with the Widom-Griffiths form [14,15] using a resummed perturbative series for the order parameter obtained in 3d [18].However, no parametrization based on MC results has been given.The previous determination of the O(2) scaling functions, using the WG ansatz [6] has been performed using an unimproved O(2) Hamiltonian and, thus, had to take care of corrections to scaling, which, in particular, made the determination of scaling functions in the symmetry broken regime difficult.Scaling functions for the O(4) model, using Monte Carlo results obtained with the standard, unimproved Hamiltonian (corresponding to λ = ∞), have been presented in [9].In none of these cases has a parametrization of O(N ) scaling functions been presented, which uses the WG form with only three free parameters.

A. Widom-Griffiths form of scaling functions
We present here a determination of the Z(2) and O(N ) (N = 2, 4) scaling functions from Monte Carlo simulations.From our new Monte Carlo results and those obtained in [9], we determine the parameters entering a parametrization of scaling functions using the WG form  of the order parameter scaling function [14,15], where (R,θ) represents an alternate coordinate frame corresponding to the (t, h) plane [16,17].Aside from the normalization constants m 0 and h 0 , this parametrization depends on a function h(θ), which needs to be deter-  mined.For the case of Z(2), it seems that a Taylor series expansion up to O(θ 5 ) is sufficient4 [18], while in the case of O(N ) one needs to take care explicitly of the presence of Goldstone modes in the symmetry broken phase.This requires that h(θ) has a double zero at some θ 0 > 1 [7,27].We, thus, use the ansatz proposed in [7,18], The normalization constants m 0 and h 0 are determined from the conditions in Eq. 22 which gives us where θ 0 is the first positive zero of h(θ) in the Z(2) universality class and the double zero of h(θ) in the O(N ) case.Using Eqs. 12 and 19 and the above relations for the normalization constants m 0 and h 0 one can establish the relation between the WG form of the scaling function f G and the relation between the scaling variables, z = t/h 1/βδ and θ, Obviously, θ = 1 corresponds to z = 0 and θ = θ 0 corresponds to z = −∞ and these, respectively, correspond to the normalization conditions for the scaling function f G in Eq. 23.Finally, θ = 0 corresponds to z = ∞.Using Eqs. 30 and 31 we also obtain f G (z) as, The presence of Goldstone modes in the symmetry broken (z < 0) phase of O(N ) symmetric models, gives rise to a distinctively different behavior of the Z(2) and O(N ) model scaling function f G in the z → −∞ limit, , for O(N ) . ( One can arrive at the above asymptotic form for f G using Eqs.30 and 31 with the ansatz for h(θ) in Z( 2 Although earlier parametrizations of the O(4) scaling function f G (z), determined in Monte Carlo simulations [5], made use of the Widom-Griffiths form, this was done only to establish the behavior of f G (z) at large |z|.The region around z = 0 has been parametrized using polynomial ansätze.A recent parametrization of the O(4) scaling function used different fits in the small and large z regions and obtained a parametrization that uses 14 parameters [9].
In order to establish the validity of the WG form using an ansatz for the function h(θ) as suggested in [7] we reparametrized the fit results presented in [9].We used the WG form for the O(4) scaling functions as given in the previous subsection and determined optimal parameters (θ 0 , h 3 , h 5 ) in an interval around z = 0, i.e., we do not make use of the large z behavior of the scaling function given in Eq. 33.The structure of this asymptotic form is implemented already in the Widom-Griffiths ansatz and the expansion parameters d − 1 and d − 2 are determined directly from (θ 0 , h 3 , h 5 ) (see Appendix B), which can be determined from any set of z-values.We determined  4) fG(z) scaling function, obtained from our fit to results in [9], is shown in the inset in the right figure.Solid lines shown in the figures are based on fits using the Widom-Griffiths form of the scaling functions and use also data outside the parameter range shown here (see the text).For O(2), we also show an error band to the WG ansatz obtained from a bootstrap analysis.The green dashed lines show the analytic results obtained in the Z(2) [18] and O(2) [7] universality classes and MC fit results obtained in the O(4) [9] universality class, respectively.The dashed red line shows the asymptotic expansion given in Eq. 33.
these parameters in several intervals [−z max : z max ] with 1 ≤ z max ≤ 6.The resulting parameters are given in Table V, where the errors quoted there reflect the spread of results for (θ 0 , h 3 , h 5 ) obtained when varying z max .In Fig. 2, we compare the scaling functions f G (z), f G (z), and f χ (z), obtained with the WG ansatz using parameters given in Table V to that obtained in [9].As can be seen, we find excellent agreement.
Even though the scaling functions themselves are in good agreement and as such give consistent results for the positions z t and z p of the maxima of −f G (z) and f χ (z), we find different asymptotic behavior at large, negative z as shown for the case of f G (z) in the inset in Fig. 3(right).In [9], the sub-leading asymptotic correction, d − 2 has been found to vanish within errors, while we find d − 2 ∼ 0.1.This difference, however, may not be too surprising, as the earlier results for the asymptotic expansion parameters d − 1 and d − 2 have been obtained from fits in the interval z ∈ [−10, −1].We will show in the next subsection that in the O(2) case the asymptotic form is not yet valid in this z-range.
Given the good agreement between the WG parametrization of the O(4) scaling functions and the earlier results based on a 14-parameter fit to MC data we find it encouraging to analyze also the new Monte Carlo simulation results, obtained for the 3-d, Z(2) and O(2) models using a parametrization based on the WG ansatz.

O(4) Monte Carlo
Monte Carlo WG-fit to [ In order to use Eqs.30 and 31 in determination of the scaling functions f G , f G , and f χ from MC results, as given in Eqs. 19, 20 and 21, one still needs to determine the non-universal scale parameters (t 0 , H 0 , L 0 ).The non-universal scales H 0 and L 0 can be determined from the finite-size dependence of f G (z, z L ) at T c , i.e., at z = 0. We present a determination of these two scales in Appendix A. Once they have been determined from our results on different size lattices, the scale parameter t 0 can be determined from the asymptotic behavior of f G (z) in the limit z → −∞.Using Eq. 13, the second normalization condition in Eq. 22 and writing z = z 0 z b with z 0 = H 1/βδ 0 /t 0 , we obtain t 0 from As this equation relates the scale t 0 to observables calculated in the infinite volume limit, its determination can directly be incorporated into fits which we perform in the infinite volume limit for the determination of the scaling functions.We obtain t 0 and the parameters (h 3 , h 5 , θ 0 ) defining h(θ) using simultaneous fits to the scaling functions f G (z), f G (z), and f χ (z) defined in Eqs.19-21.While in the parametrization of Z(2) scaling functions θ 0 is a function of (h 3 , h 5 ), it is an additional free parameter in the O(N ) case.
Goldstone modes dominate finite-size effects at large, negative values of z, which are quite different in O(N ) universality classes from those in the Z(2) case.In the O(2) universality class finite-size effects grow rapidly with decreasing values of z.This is evident from Fig. 3, where we show results for (−z) −β f G (z) in the region z < −2.The figure shows that we had to perform MC calculations on rather large lattices to extract the scale parameter t 0 from the asymptotic behavior of the order parameter in the symmetry broken phase.In simulations of the O(2) model, lattices of size L 3 with L = 200 were needed to reach the region z ≤ −10 without suffering from finite-size effects.In the case of Z(2), lattices with L = 96 were already sufficient to perform calculations in a region down to values z −20 without observing a significant finite-size dependence in our results.
For z > −2, it was sufficient to perform calculations on lattices with L = 48 − 120.For z < −2, however, we also performed calculations with L = 160 and 200 for the O(2) model and L = 160 in the case of Z(2).The statistics collected in all parameter ranges are given in Tables II-IV.
Our Monte Carlo results obtained for the 3-d, Z(2) and O(2) models in the large volume limit, z L ≤ 0.35, are shown in Fig. 4. We performed joint fits using data in the region z L ≤ 0.35 as approximation for the infinite volume limit.All three scaling functions have then been obtained from joint fits to the WG form in the range z ∈ [−23 : 2] for Z(2) and z ∈ [−12 : 2] for the O(2) model 5 .
We summarize results for the non-universal scale parameters (t 0 , H 0 , L 0 ), determined by us, in Table I.In Table VI we give all universal fit parameters entering the definition of h(θ), and compare with results obtained in 3-d analytic calculations [7,18].In the top section of the two tables, we give the parameters (h 3 , h 5 ), entering fits performed for scaling functions in the Z(2) universality class, and (θ 0 , h 3 , h 5 ) in the O(2) case.Results for the non-universal fit parameter t 0 , obtained in the same fits, are given in Table I.The bottom part of the tables gives some universal constants derived from the Widom-Griffiths form of the scaling functions by using, on the one hand, the results from fits to our MC data and, on the other hand, the perturbative results for h(θ) and θ 0 as input.
Aside from the parameters d − 1 and d − 2 controlling the asymptotic behavior of f G (z) at large, negative z (Eq.33) we also give there the universal constants z p and z t , which are the z-values at the maxima of f χ (z, 0) and −f G (z, 0), respectively.
For the ratio of z p and z t , determining pseudo-critical temperatures in the Z(2) and O(N ) universality classes, we find In Fig. 3, we compared the MC results for f G (z) at large, negative values of z, i.e. for z < −2, with the WG form of the scaling function, given in Eq. 30, as well as with the asymptotic form given in Eq. 33.As can be seen, in the Z(2) universality class, the asymptotic expansion using the first two sub-leading corrections, gives a good approximation to the full WG form, in almost the   2)) for z < −2.5.All fits are joint fits to data for fG, f G and fχ, close to z = 0 and in the large, negative z regime.The gray data were not included in the fit.Green lines in the insets show results from analytic calculations [18] (Z(2), left) and [7] (O(2), right).entire region, z < −2.In the O(2) and O(4) universality classes, however, the first two sub-leading corrections agree with the full WG form only for z < −(8 − 10) as can be seen in Fig. 3(right) for the O(2) case and in the inset in Fig. 3(right) for the O(4) case.
Also shown in Fig. 3 are the results obtained from the 3-d, analytic calculations [18,27].While in the Z(2) case differences are insignificant, they clearly are visible in the O(2) case.However, in the asymptotic regime both parametrizations of the WG form differ by less than 1%.We observed the largest differences in the vicinity of the maximum of −f G (z), where deviations between the analytic and MC calculation amount to about 5%.This is apparent from the insets shown in Fig. 4.

V. FINITE SIZE SCALING FUNCTIONS
We now want to determine corrections to the infinite volume scaling functions in the 3-d, Z(2) and O(2) universality classes arising in a finite volume at small external field H.These corrections are universal when taking the limit (H → 0, L → ∞) at fixed z L as introduced in Eq. 12.
In the limit of small H and in the vicinity of T c , we obtain the scaling functions (f from the order parameter M and the two susceptibilities χ h and χ t using Eqs.19-21.We focus here on the region in the vicinity of T c and the pseudo-critical temperatures, T pc,h and T pc,t , determined from the maxima of the susceptibilities χ h and χ t , respectively.It is this region where correlation lengths are large and where it is of particular importance to get control over finite-size effects in the determination of pseudo-critical and critical temperatures in many models belonging to the Z(2) and O(N ) universality classes.For this reason, we determine finite-size scaling functions with parameter sets (J, H) corresponding to the interval z ∈ [−1 : 2].A similar calculation has been performed previously for finite size scaling functions in the O(4) universality class [10].
In our analysis of finite-size effects, we use a polynomial ansatz for the scaling functions which has also been used previously for calculations in the 3-d, O(4) universality class [10], For the infinite volume scaling function, f G (z, 0) ≡ f G (z) we use the parametrization determined in the previous section.Here, (n u , m l , m u ) denote the lower and upper limits of the sum over the polynomial in powers of z and z L , respectively.We take the leading order finitesize correction to be inversely proportional to the volume O(1/L 3 ), i.e., m l = 3.The upper limits n u and m u are optimized in our fits, using the Bayesian information criterion.We fix a 0m l = 0 in both universality classes; additionally we constrain the fit parameters to |a nm | < 10.
From the ansatz used for f G (z, z L ), one also obtains the parametrization of f G (z, z L ), which controls the scaling behavior of χ t , and f χ (z, z L ), which controls the scaling behavior of χ h , Using these polynomial ansätze we again perform joint fits to the MC data for the three scaling functions (f G , f G , f χ ) in the interval z ∈ [−1 : 2] and for z L ∈ [0.4 : 1.0].The data for z L < 0.4 have been excluded from these fits, as they have been used already to determine the parameters of the infinite volume scaling functions, as discussed in the previous section.
Results obtained for the finite-size scaling functions in the 3-d, Z(2) and O(2) universality classes for some fixed values of z have been shown already in Fig. 1.In Fig. 5, we show results for the scaling functions as function of z for several fixed values of z L .The fit parameters obtained with the polynomial fit ansatz, Eq. 36, are given in Table VII for the case of Z(2) and in Table VIII for the case of O(2).These fits provide a good interpolation for our data in the range z L ∈ [0.4 : 1.0].However, due to the large number of parameters involved, we cannot give significance to individual parameters entering the polynomial ansatz.We, therefore, quote our fit result without assigning errors to the fit parameters.
As can be seen, the general z L -dependence of scaling functions f G (z, z L ) and f χ (z, z L ) is similar in the Z(2) and O(2) universality classes.However, it is apparent from the upper row in Fig. 5 that finite-size effects are larger in the O(2) case than for Z(2).In the latter case, results for f G (z, z L ) are indistinguishable from the infinite volume results already for z L < 0.6, whereas in the O(2) case at z L = 0.6, deviations from the infinite volume values amount to about 3% at z = −1 and increase to 4% at z = 1 (see also the discussion of Fig. 8 in Appendix A).Furthermore, qualitative differences are evident in the z L -dependence of the scaling function f G (z, z L ).In the Z(2) case, the approach to the infinite volume limit is non-monotonic for z L < 0. A pronounced peak shows up in the symmetry broken regime (z ≤ 0) at finite z l , and the asymptotic infinite volume limit is approached from above.In the case of the O(2) universality class, f G (z, z L ) seems to approach the infinite volume limit result from below for all z.
In the case of f χ (z, z L ) the approach to the infinite volume limit is non-monotonic for z-values below the pseudo-critical scale, z < z p .As can be seen in Fig. 1, this is the case in the Z(2) as well as in the O(2) universality class.This non-monotonic behavior is not that prominently visible in Fig. 5, as it sets in only at rather large values of z L , i.e. for z L > 1.This regime is not covered in Fig. 5.
While the finite-size effects seen in the scaling functions are generally larger in the O(N ) than in the Z(2) universality class, this is not the case for the location of maxima in the scaling functions −f G (z, z L ) and f χ (z, z L ).These maxima are controlled by universal functions z t (z L ) and   z p (z L ), respectively.We determined them from the polynomial obtained from the finite-size scaling fits for Z(2) and O(2).In the case of O(4), we have used the finitesize fit given in [10].Results are shown in Fig. volume limit result are described well with ansatz with b p 4.5 and b t 5.5 as shown in Fig. 6.

VI. CONCLUSIONS
We determined the infinite volume scaling functions in the 3-d, Z(2), O(2), and O(4) universality classes using a 2 or 3 parameter parametrization based on the analytic Widom-Griffiths scaling form.We find good agreement of the O(4) parametrization with an earlier parametrization that used O(10) parameters [9].In the Z(2) case, we find excellent agreement between our parametrization based on Monte Carlo results and the analytic result obtained from a perturbative, field theoretic approach [18].The largest differences between our Monte Carlo results and analytic calculations [7] we find, in particular, for the scaling function f G (z), which controls the scaling behavior of mixed susceptibilities.
We determined the finite-size dependence of the scaling functions and showed that qualitative differences between the Z(2) and O(N ) cases show up most prominently in the scaling function f G (z, z L ) which controls pseudo-critical and critical behavior of the mixed susceptibilities.We could show that the location of the pseudocritical temperature, corresponding to z t , is less affected by finite-size effects than the pseudo-critical temperature determined by the maximum of the order parameter susceptibility (χ h ) at z p .This difference is particularly striking in the O(4) universality class.The comparison of the finite-size dependence of the scaling functions among different universality classes has been possible with our proposed normalization condition for the non-universal scale parameter L 0 .
We furthermore find, z p /z t 2, i.e. at non-zero values of the symmetry breaking parameter H deviations of the pseudo-critical temperature T pc,t from the phase transition temperature T c are about a factor of 2 smaller than that of T pc,h .All data presented in the figures of this paper can be found in Ref. [28].4) universality class are taken from [10].For this purpose, the scaling variable zL has been rescaled using L0 = 0.7686 to be consistent with the normalization condition, Eq. 24, used for the Z(2) and O(2) universality classes.3 (θ 2 0 − 1) 2 h (2) (θ 0 ) 2 6βθ 0 h (2) (θ 0 ) δ (2β − 1)θ 2 0 + 1 It should be noted that θ 0 is an independent parameter in the parametric representation of O(N ) universality class while it is a function of parameters h 3 and h 5 in the Z(2) case.

Figure 2 .
Figure2.The O(4) infinite volume scaling functions.The WG parameters have been obtained from a fit to fG, while f G and fχ have been obtained from there.Dashed red lines show results from previous calculations[9].
) and O(N ) universality classes given in Eq. 28.Explicit expressions for d − 1 and d − 2 in terms of the WG parameters h 3 , h 5 and θ 0 are given in Appendix B. B. Representation of O(4) scaling functions using the Widom-Griffiths form

Figure 3 .
Figure 3.The scaling function fG(z) in the (left) and O(2) (right) universality classes obtained from the order parameter M using Eq.19 in the region of large, negative values of z.Monte Carlo data have been obtained in simulations using the 3-d, Z(2) model with λ = 1.1 and the O(2) model with λ = 2.1.All data are from simulations at T /Tc = 0.99.The large negative z region of the O(4) fG(z) scaling function, obtained from our fit to results in[9], is shown in the inset in the right figure.Solid lines shown in the figures are based on fits using the Widom-Griffiths form of the scaling functions and use also data outside the parameter range shown here (see the text).For O(2), we also show an error band to the WG ansatz obtained from a bootstrap analysis.The green dashed lines show the analytic results obtained in the Z(2)[18] and O(2)[7] universality classes and MC fit results obtained in the O(4)[9] universality class, respectively.The dashed red line shows the asymptotic expansion given in Eq. 33.

Figure 4 .
Figure 4.The Z(2) (left column) and O(2) column) infinite volume scaling functions.Shown are Monte Carlo results obtained on lattices of size L 3 with L = 96 for |z| ≤ 2.1, L = 96, 120, 160 (Z(2), only L = 160 is shown) and L = 200 (O(2)) for z < −2.5.All fits are joint fits to data for fG, f G and fχ, close to z = 0 and in the large, negative z regime.The gray data were not included in the fit.Green lines in the insets show results from analytic calculations[18] (Z(2), left) and[7] (O(2), right).

Figure 5 .
Figure 5. Fits to data for scaling functions of in the Z(2) (left column) and O(2) (right column) universality class.All fits are joint fits to data for fG, f G and fχ, done in the intervals z ∈ [−1.0, 2.0] and zL ∈ [0.4,1.0] on lattices of size L = 48 (circle) and 96 (squares).Also shown are the infinite volume lines for zL = 0. Green crosses in the upper row mark the normalization condition fG(0, 1) = 0.7.

Figure 7 .Figure 8 .
Figure 7.The rescaled order parameter H −1/δ M versus the bare finite-size scaling variable z L,b at Tc.The left hand figure shows results for the 3-d, Z(2) model with λ = 1.1, and the right hand figure is for the O(2) model with λ = 2.1.The inset shows the region of small z L,b that is used for the determination of the non-universal scale parameter H0.The green cross marks the value of z L,b that determines the scale parameter L0 using the normalization condition given in Eq. 24.

Table VI .
Top: Fit parameters h3 and h5 for Z(2) infinite volume scaling functions in the Widom-Griffiths form appearing in the function h(θ) introduced in Eq. 28.Bottom: Fit parameters h3, h5 and θ0 for O(2) infinite volume scaling functions in the Widom-Griffiths form appearing in the function h(θ) introduced in Eq. 28.In the lower part of both tables, we give results for several universal constants computable from the WG parametrization.
Figure 6.Finite-size dependence of the location of maxima in the scaling functions fχ(z, zL) and −f G (z, zL).Shown are the universal functions zp(zL) and zt(zL) for the 3-d, Z(2) and O(N ) universality classes.Dashed lines show simple polynomial approximations (Eq.39) with parameters as given in the figure.
6.It is clearly seen that the finite-size dependence of the maxima in χ h is stronger than that of χ t in the O(N ) universality classes and vice versa in the Z(2) case.Moreover, the finite size dependence of z t and z p is stronger in the Z(2) universality class than in the O(N ) cases.Over a wide range of z L -values, the deviations from the infinite