Fluctuation of the phase boundary in the six-vertex model with domain wall boundary conditions: a Monte Carlo study

We consider the six-vertex model with domain wall boundary conditions in a square lattice of dimension N × N. Our main interest is the study of the fluctuations of the extremal lattice path about the arctic curves. We address the problem through Monte Carlo simulations. At Δ=0 , the fluctuations of the extremal path along any line parallel to the square diagonal were rigorously proven to follow the Tracy-Widom distribution. We provide strong numerical evidence that this is true also for other values of the anisotropy parameter Δ ( 0⩽Δ<1 ). We argue that the typical width of the fluctuations of the extremal path about the arctic curves scales as N1/3 and provide a numerical estimate for the parameters of the scaling random variable.


Introduction
In probability theory [1], the sum of a large number of random variables once divided by the number of summands might approach a deterministic value.For instance, the law of large numbers implies that the sample mean of N independent identically distributed random variables converges to their distribution average for large N .Moreover, if those random variables have finite variance, fluctuations of their sum are normally distributed with a variance of order √ N .The convergence toward smooth functions of combinations of a large number of random variables is at the root of the effectiveness of the methods of statistical mechanics and also provides a solid ground for numerical techniques such as Monte Carlo.This phenomenon is not restricted to local observables but extends to geometrical quantities such as random curves.An instructive example discussed in [2] (but see also [3]) is a one-dimensional random walk starting from the origin and conditioned in T steps to reach the point X.For X, T → ∞ with X/T = v finite, the random paths performed by the walker are with probability one within any neighborhood of the classical trajectory X = vT .Limit shapes [4] are smooth nonrandom curves that describe the thermodynamic limit of random lattice paths.Historically, one of the most studied examples is the arctic circle [5] that appears in dimer coverings of the Aztec diamond.This domino tiling problem is equivalent [6] to the six vertex model [7,8] with Domain Wall Boundary Conditions (DWBC) at its free fermionic point [9,10,11].
Let N be a positive integer, a configuration of the six vertex model on a N × N square lattice with DWBC [12,13,14,15,16,17] is specified by N non-intersecting paths traveling on the edges of the lattice.The lines enter from above and exit on the left.On each vertex, the path can have the six possible configurations given in Fig. 1; the Boltzmann weights associated with those are a for the vertices of type a 1 and a 2 , b for the vertices of type b 1 and b 2 and c for the vertices of type c 1 and c 2 .The phase diagram of the model is characterized by the ratio b/a and the parameter ∆ [9], defined as In this work, we will focus on the domain 0 ≤ ∆ < 1, which is included in the so-called disordered phase |∆| < 1 see [18] for a comprehensive study of the phase diagram.The case ∆ = 0 is the free fermionic point mentioned a few lines above.
The N random paths cannot intersect and for large enough N they can fluctuate only within a region whose boundary, after rescaling all the distances by N , is arbitrarily close to a curve, dubbed the arctic curve [5].For ∆ = 0 and a = b, such a curve is the circle (x/R) 2 + (y/R) 2 = 1, depicted in Fig. 2 and inscribed in the square lattice.Inside the arctic curve the model is in a disordered phase, outside is in a trivial ferromagnetic phase.The Cartesian coordinates are also defined in Fig. 2 and lengths, such as the circle radius R, are measured in suitable lattice units.
In this paper, we follow the convention to divide the square diagonal into N intervals, therefore x, y are semi-integers in the interval [−N/2, N/2]; for simplicity, we assume N a power of two from now on.The Cartesian equation of the arctic curve depends on the parameter ∆ and was calculated in [19,20], see also [3].A rigorous derivation at ∆ = 1/2 has been given in [21].
The limit shape, however, is a property of the model in the thermodynamic limit and it is relevant to understand how the interface between the disorder and the ferromagnetic phase fluctuates due to finite-size effects.In this case, the phases are not separated by a smooth curve; rather, the distance from a corner to the disordered phase varies according to a probability distribution.Consider, for instance, a quarter of the original lattice: the southeast (SE) quadrant on the right of Fig. 2. When N is large, there are no paths in a neighborhood of the SE corner and all the vertices are of type a 2 .For later purposes, introduce the parameter η as η := 2y/N. ( Take |η| < 1/2; as we move along a constant-y line toward the boundary (for x > 0), we cross N/4 lattice paths, see Fig. 2. Denote then by X ext (η) the x-coordinate of the intersection of a constant-y line with the furthest path from the origin.Refs.[22,23] tackled rigorously the study of the fluctuations of the extremal path in the context of dimer coverings.The same conclusions apply to the six-vertex with DWBC at ∆ = 0, [6].In particular, it can be proven that for x y The furthest path from the center of the lattice intersects a constant y line in a point with abscissa X ext (η).The picture shows X ext (η = 0).

X ext
where Z follows a Tracy-Widom (TW) distribution [24] for any |η| < 1/2.The sign in front of Eq. ( 2), is the one of X ext (η) in the coordinate system defined in Fig. 2.Moreover, the parameters Λ and Γ are O(1): ΛN is the x-coordinate of the intersection of the arctic curve (an ellipse if a ̸ = b) with the constant-y line while Γ can be derived from Ref. [25], see Sec. 3. The TW distribution characterizes fluctuations of extrema of random variables in several different contexts, among them: the spectrum of Hermitian random matrices [24], ordered sequences in random permutations [26], random growth models [22,27] and their generalizations [28,31], the KPZ [32] equation [33], quantum dynamics [34,25,35].Many of these problems can be mapped into each other [36] and one might wonder [37] how universal is the TW distribution.Ref. [29] proved for instance that the maximum reached by the extremal path at ∆ = 1/2 follows a TW distribution (although with β = 1).
In this paper, we will study numerically the probability distribution of X ext (η) for values of ∆ ̸ = 0 in the interval [0, 1).We will argue that after a rescaling analogous to Eq. ( 2), this random variable follows, for large enough N , again a TW distribution.This result is not unexpected and can be justified heuristically.Non-intersecting lattice paths can be mapped to fermionic trajectories in a Euclidean space time [25].The fermionic particles are interacting if ∆ ̸ = 0, however the particle furthest from the origin, which is also the fastest, behaves as if were free [36,38,39].This observation is also at the root of the so-called tangent method [3].
The rest of the paper is organized as follows.In Sec 2, we review briefly the numerical method.In Section 3 we discuss the comparison with the case ∆ = 0 and further analyze the lattice path fluctuations for ∆ ̸ = 0 in Sec. 4. The analysis in Sec. 4 focuses on the lattice The picture also shows an allowed state.This state has four flippable plaquettes.The plaquettes indicated in blue, are second flippable while the plaquettes indicated in red are first flippable.
diagonal (η = 0), while some comments on results at η ̸ = 0 are contained in the two conclusive appendices.Finally, in Sec. 5 we summarize our conclusions.
2 Monte Carlo study of the six-vertex model with DWBC Method.-We modify a Monte Carlo (MC) algorithm for simulating the six-vertex model with DWBC that was originally proposed in [40].The MC algorithm involves two moves called first flip and second flip; see Fig. 3.The moves modify the configuration of the paths along the edges of a plaquette of the square lattice.The state S of a plaquette is a list that contains the configurations of its four vertices: S = {v 1 , v 2 , v 3 , v 4 }, with v i given in Fig. 1.The Boltzmann weight of v i will be denoted by w(v i ) = a, b, c.The weight of the state S is then W (S) = 4 i=1 w(v i ).An MC move takes then a state S to a new state S ′ .Since in a flippable plaquette, each vertex can have two possible configurations, there are sixteen distinct plaquette states which are respectively first or second flippable.They have up to nine distinct weights.
The computer code has two lists, each containing the plaquettes flippable in one of the two different ways.We can formulate the algorithm as follows.First, randomly search among all plaquettes on the lattice, until a plaquette contained in one of the two lists is found.A flip will then be attempted, and if the attempt is successful, then both lists will be updated: at the plaquette itself and at the four plaquettes sharing a side with it.The probability P of an attempted flip that modifies the state of the plaquette S into S ′ is the Boltzmann weight of the new state, namely In the computer program, the largest of the three weights a, b, and c is 1 and therefore Eq. ( 3) is well defined.The above conditions satisfy detailed balance and ensure ergodicity.An equivalent MC algorithm with Glauber dynamics has been employed in [41] and [42] to obtain numerical estimates for the arctic curves and one-point functions of the model with several values of ∆ in the disorder phase.To ensure the system has thermalized we rely then on the results of those papers.We define a MC sweep as N 2 accepted flips.The largest value of N used in this paper is 256 whereas the longest side of the lattice used in the mentioned papers was as long as 500 units.We, therefore, consider it known at this time how long this system needs to equilibrate.From a rigorous perspective, Ref. [43] shows that the Glauber dynamics is rapidly mixing at Figure 4: Left.We show the numerical data for distributions of the random variables X ext (η = 0) and X ext (η = 1/4) along the lines y = 0 and y = N/8, see also Eq. ( 1), at ∆ = 1/2 and N = 256.Right.A realization of the random variable X ext (η) for two different values of η at ∆ = 1/2 and N = 256.∆ = 1/2, while Refs.[44,45] proved that the Markov chain exponentially slows down for large and negative ∆ 1A reliable test of the TW distribution is challenging in a statistical model, see for instance [46] and requires a much larger number of sweeps than that employed previously in Refs.[41,42].In this paper, we will construct samples with 10 6 MC sweeps.
In the disorder phase, analogous MC algorithms are nowadays commonly adopted in numerical experiments with the six vertex model [47,48,49,50,51,52].For earlier numerical approaches to ice models, see [53,54].Finally, an efficient algorithm to simulate interacting dimers has been proposed in [55] but so far it has not been implemented in the six-vertex model.The extremal lattice path.-In the coordinate system2 defined in Fig. 1 with origin in the center of the lattice, the random variable X ext (η) is the position, along a line at constant y, of the last occurrence of a vertex of different type from those in the frozen corner.The right panel of Fig. 4 illustrates a realization of X ext (η) for η = 0 and η = 1/4 on the NW corner of lattice with N = 256 and ∆ = 1/2.Nearby this corner the vertices are in the state a 1 , see Fig. 1 and Fig. 2; the black dashed line is the arctic curve [19].
The probability distribution, Prob(X ext (η) = x), will be approximated numerically at finite N by a normalized histogram which we denote by h η ; see Appendix.A for more detail about its construction.The left panel of Fig. 4 shows two examples at η = 0 and η = 1/4 for ∆ = 1/2 and N = 256.The two histograms are not symmetric (the right tail is longer than the left) and obviously have different variance and mean.Their properties will be studied in the next two Special attention must be paid to the case ∆ = 0, since it is the only point at which there are known exact results for the random process that describes the fluctuations of the coordinate X ext (η).Let us take for simplicity a = b; in [23,22], the random fluctuations were proven for large N to converge to the so-called Airy process for any |η| < 1/2.By following [25], we can derive analytically the parameters Γ and Λ in Eq. (2).It turns out that where the two branches of Eq. ( 4) depend on the sign of X ext (η).Notice that Γ(η) is positive and monotonically decreasing with Γ(1/2) = 0.When the right-hand-side of Eq. ( 5) vanishes, the intersection point X ext (1/2) = N/4 lies on the boundary of the square-lattice and its fluctuations are no longer described by the TW distribution.In this case, Ref. [56], see also [57] for a related discussion, proved that the latter are Gaussian with an horizontal width O( √ N ).We then exploit Eqs. ( 4) and ( 5) to benchmark the validity of the MC simulation.In particular, according to Eq. ( 2), the probability distribution (Prob(Z(η) = χ)) of the random variable Z can be evaluated numerically as3 Upon substituting Eqs. ( 4) and (5) into Eq.( 6), one compares the MC data at finite N with the TW distribution F ′ 2 (χ).The result of this comparison is displayed on the left of Fig. 5 for η = 1/4 and N = 256.It shows a convincing agreement already for such a value of the lattice size; it should be remarked that there are no free parameters.From a statistical viewpoint, it is possible to verify if the random variable Z is distributed according to F ′ 2 , by performing a Kolmogorov-Smirnov test [58].The latter requires to construct the empirical cumulative distribution of Z, which will be denoted by F η (χ).On the right panel of Fig. 5 is plotted F η (χ) at η = 1/4 and N = 256 together with the cumulative distribution of the TW distribution, denoted by F 2 (χ).The test statistics [59,60] is where n is the number of points in Fig. 5.For n = 14, we obtain D 14 = 1.021 . . ., which corresponds to a p-value of about 0.2 [61], well above the conventionally accepted threshold of 0.1.
4 Numerical study of the cumulants for 0 ≤ ∆ < 1 and η = 0 In this section we report numerical results for the fluctuations of the random variable X ext (η = 0), i.e. along the diagonal of the square lattice for values of the anisotropy in the domain 0 ≤ ∆ < 1.Details about the construction of the empirical probability distribution of fluctuations of the extremal lattice path are given in the Appendix A. To lighten the notation we will drop the η-dependence in the remaining part of the discussion and always assume η = 0.At large N , and for any |∆| < 1 the fluctuations of X ext shall obey (see Eq. ( 2)) where Z is a TW random variable and the parameters Λ and Γ are non-universal and O(1).In particular Λ can be estimated from the mean of X ext and Γ from its variance.It should be understood that Eq. ( 7) is an asymptotic expansion valid for large enough N .In particular, the coefficients Λ and Γ are the N → ∞ limit of non-trivial functions of the lattice size, which are not known analytically for ∆ ̸ = 0.The numerical data will be then extrapolated to N → ∞ by assuming the same functional form for the finite size-corrections of the distribution cumulants discussed in [38].A numerical study of Λ has been perfomed in detail in Ref. [41] and will be not repeated here.Skewness and Excess Kurtosis.-Sincethe random variable X ext is related to the TW random variable Z in Eq. ( 7) by a linear (affine) transformation, cumulants of order larger than two such as the skewness and excess kurtosis, are the same for the two distributions.Their study [30] provides a convenient parameter free test of the universality conjecture and will be considered first.We denote then by s N and k N the empirical skewness and excess kurtosis obtained from the Monte Carlo data for X ext at finite N .Following Ref. [38] we assume that the leading finite-size correction for large N to the cumulants is O(N −2/3 ), namely where o(N −2/3 ) represents quantities that vanish faster than N −2/3 for N → ∞.The exponent 2/3 of the first subleading correction in Eq. ( 8) was proposed in Ref. [38] in the context of the quantum XXZ spin chain when analyzing fluctuations of the density profile of the interacting fermions close to the edge.Albeit the numerical methods (Ref.[38] employed DMRG techniques) and the systems are different, the convergence patterns of cumulants toward the thermodynamic limit are resemblant.
A linear fit of the numerical data for s N and k N against N −2/3 is shown in Fig. 6 by taking into account points at N = 64, 128, 256 for ∆ = 0, 3  8 , 1 2 and 5 8 .The numerical estimates for the coefficients s and k in Eq. ( 2) correspond then to y-intercept of the regression lines; they are collected in Tab. 1 together with their statistical uncertainties.The red dashed line in Fig. 6 is the TW value.The agreement is convincing although not prefect for the skewness at ∆ ̸ = 0 and ∆ ̸ = 1/2.The data points in Fig. 6 and their error bars are obtained by averaging over several Monte Carlo histograms (5-10 histograms).Variance at the ice point.-Letus denote by σ 2 N the variance of the emprical distribution X ext at finite N and let σ 2 the variance of the TW distribution.From Eq. ( 7) and the same scaling conjecture proposed by Ref. [38] for the higher cumulants, it follows In Sec. 3, we recalled the expression for the parameter Γ 3 along the lattice diagonal at ∆ = 0, namely 2 see Eq. ( 5) for η = 0.More recently, Ref. [63] proposed a phenomenological (and not rigorous) formula for Γ 3 at arbitrary values of ∆, which involves the second derivative of the arctic curve evaluated at the intersection point with the lattice diagonal.At ∆ = 1/2, the so-called ice point, the arctic curve is an ellipse [19] and the parameter Γ 3 was conjectured to be [63]: 3 .The Ansatz elaborated in [63] relies on the KPZ [32] scaling theory applied to the fluctuations of the extremal lattice path.Fig. 7 shows a numerical estimation of the parameter Γ 3 obtained from a linear regression analogous to that performed for the skweness and excess kurtosis in the previous paragraph.The agreement with the exact [25] and conjectured [63] expressions at ∆ = 0 and ∆ = 1/2 is prefect.An equivalent result can be obtained by averaging the distribution of the extremal lattice path along lines parallel to the Figure 7: Numerical estimation of the coefficient Γ 3 in Eq. ( 7) from the variance of the empirical distribution of X ext along the lattice diagonal (η = 0).The predictions at ∆ = 0 and ∆ = 1/2 are the y-intercepts of the regression lines, see Eq. ( 9) and perfectly agree with the theoretical values represented by the dashed lines, see main text.
lattice diagonal and requiring best fit of the Monte Carlo data with the TW distribution, see Appendix B.

Conclusions
In this paper, we analyzed numerically, with Monte Carlo techniques, the fluctuations of the extremal lattice path in the six vertex model with DWBC for values of the anisotropy parameter in the domain 0 ≤ ∆ < 1. Exact results at the free fermion point (∆ = 0) show that the intersection of the furthest lattice path with the lattice diagonal follows the Tracy-Widom distribution.The same is true along any line parallel to the diagonal, as long as the contact point [56] of the arctic curve with the lattice is avoided.Our numerical simulations were tested thoroughly against these analytic results, showing convincing agreement already for lattices with side-length N = 256 and without fitting parameters.Further, we provided strong numerical evidence that the fluctuations of the extremal lattice path follow the Tracy-Widom distribution also for values of the anisotropy parameter in the domain 0 ≤ ∆ < 1.In particular, we verified the consistency of the scaling hypothesis in Eq. ( 2) by determining the skewness and excess kurtosis of the distribution of X ext along the lattice diagonal.We also verfied, see Eq. ( 9), that the width of the fluctuations about the arctic curve is of order O(N 1/3 ) with a prefactor that could be extrapolated from the variance of the Monte Carlo data and agrees perfectly with the exact value at ∆ = 0 and the recent conjecture in [63] at ∆ = 1/2.
Our finite size scaling analysis requires relatively small lattices and samples of O( 10 6 ).Therefore it is computationally less demanding than other Monte Carlo attempts to uncover the Tracy-Widom distribution [46].The results are fully consistent with those that Ref. [38] obtained in quantum spin chains through DMRG simulations.
Our conclusions should extend also to cases −1 < ∆ ≤ 0, which we did not analyze in this paper but belong to the same disordered phase [18] of the model.Finally, it would be interesting to study fluctuations of the boundary of the inner phase separation curve in the antiferromagnetic regime, see for instance or understand the statistics of interior lattice paths [65].To circumvent the problem of exponentially slow mixing of the Markov chain Monte Carlo dynamics for large and negative ∆, one might implement the algorithm proposed recently in [66].systems of Fig. 2 coincides indeed with the value of the x-coordinate in the NW corner.As mentioned in Sec. 1, the diagonal distance between two neighboring vertices is 1.In this paper, T will have the value T = 10 6 throughout.
In order to provide evidence of the validity of the scaling hypothesis in Eq. ( 7) also at η ̸ = 0 and away from the free fermionic point, we construct k max histrograms of X ext (η ̸ = 0), i.e. along lines parallel to the lattice diagonal.The vectors {v k=1 are built taking into account that the process is begun on side of the lattice and not in the corner.Therefore, there are eight simultaneous processes and not four.We thus obtain another set of vectors {v From these, we will obtain normalized, reversed histograms {h η k } if, for each k = 1, . . ., k max , we define h η k as The histogram h 1/4 (x) obtained when k = k max = N/4 is displayed for ∆ = 1/2 and N = 256 on the bottom left of Fig. 4.

B Distribution average over lines parallel to the lattice diagonal
We assume that all the histograms h η k obtained in Appendix A once shifted and rescaled to have the same mean and variance converge for large N to the same distribution.To test the scaling hypothesis in Eq. ( 2), we then average over them as described in the following.First one calculates the mean µ k and the variance σ k of the distributions h η k for all k = 0, . . ., k max .We have At this point, we will leave the actual lattice and define our distributions in a new, one-dimensional space.To perform with greater numerical accuracy the averaging procedure, it is convenient to magnify the bin width by a constant factor θ, such that 1 < θ < k max .In terms of the new variable t = θx, the probability distributions in Eq. ( 14) are then4 The value θ = 8 will be used in this paper.By following the assumption stated at the beginning of this section, we further perform for all k = 0, . . ., k max the linear change of variable which transforms the probability distributions of Eq. ( 16) into It is easy to show, by recalling Eq. ( 15), that all the probability distributions in Eq. ( 18) have zero mean and variance (θσ 0 ) 2 .According to our working hypothesis, they are then different histogram representations of the same probability distribution of the random variable In Eq. ( 19), C is a constant that ensures that Ξ has zero mean.We estimate the probability distribution of the random variable Ξ with the average in Eq. ( 20) the variable ξ belongs to the union of all the support of the histograms h η k in Eq. ( 18).Instead of calculating the cumulants as in Sec. 4, one could test the scaling hypothesis in Eq. ( 2) (and Eq. ( 7)) by directly fitting Eq. ( 20) with the TW distribution.In particular, recalling Eq. ( 19), we analyze weather h0 (ξ) with F ′ 2 , the TW distribution, by fitting at finite N the two parameters Λ N and Γ N above.As also discussed already in Sec. 4, Λ N and Γ N are non-trivial and unknown functions of the lattice size N such that (see Eq. Working with N up to 256, it becomes crucial to evaluate the leading finite size corrections to those parameters, in particular following [38] one has, see Sec.The numerical prediction for the coefficient [Γ(0)] 3 in Eq. ( 2) corresponds then to y-intercept of the regression lines.As in Sec. 4, see Fig. 7, the red dashed lines are the exact [25] and conjectured values [63] at ∆ = 0 and ∆ = 1/2 respectively.The agreement is perfect also with this approach.Finally, Eq. ( 21) is equivalent to the statement that for large enough N f0 (χ) := Γ N θ h0 Γ N θ χ + θΛ N + C (24) converges to the TW distribution for any 0 ≤ ∆ < 1.On the left of Fig. 8, we show on the same plot the numerical results for f0 (χ) at N = 256 and ∆ = 0, 3 8 , 1 2 , 5 8 .The continuous red curve underneath, which is barely visible, is the TW distribution F ′ 2 (χ).

2 Figure 1 :
Figure 1: The six possible configurations of a lattice path in the six-vertex model are denoted by thick lines.The continuous lines are lattice edges while dashed lines are edges of the dual lattice.

Figure 2 :
Figure 2: Left.Lattice path configuration in the six vertex model with DWBC for N = 64, and a = b at ∆ = 0.The possible vertex configurations have been colored differently: vertices of type a 1 are colored in red, of type a 2 in white, of type b 1 in yellow, of type b 2 in blue while of type c 1 and c 2 in green.The black circle is the arctic curve for the given choice of the parameters.Right.The SE corner of the lattice path configuration on the left for N = 64.The furthest path from the center of the lattice intersects a constant y line in a point with abscissa X ext (η).The picture shows X ext (η = 0).

Figure 3 :
Figure 3: Left.The two MC moves.Right.The six-vertex model with DWBC with N = 4.The picture also shows an allowed state.This state has four flippable plaquettes.The plaquettes indicated in blue, are second flippable while the plaquettes indicated in red are first flippable.

Figure 5 :
Figure 5: Left.Distribution of the random variable Z(η) defined in Eq. (2) calculated from the MC simulation at η = 1/4 and N = 256 for ∆ = 0.The red curve is the TW distribution F ′ 2 (χ).There are no free parameters in the comparison.Right.The empirical cumulative distribution of the random variable Z 0 (η) is plotted together with the cumulative distribution F 2 at η = 1/4 and N = 256 for ∆ = 0.

Figure 6 :
Figure 6: Left.Numerical data for the skewness of the empirical distribution X ext (η = 0) as a function of N −2/3 for N = 64, 128, 256.The intercept of the regression line is the skewness in the thermodynamic limit.The dashed red line is the TW value at β = 2 which equals 0.224 . . . .Right.Numerical data for the excess kurtosis of the empirical distribution X ext (η = 0) as a function of N −2/3 for N = 64, 128, 256.The intercept of the regression line is the excess kurtosis in the thermodynamic limit.The dashed red line is the TW value at β = 2 which equals 0.093 . . . .

Figure 8 :
Figure 8: Left.Numerical data for the distribution f0 (χ) in Eq. (24) at different values of ∆ and N = 256.The parameters Γ N and Λ β are obtained by requiring a best fit of the data with the TW distribution, which is the red curve underneath.Right.Numerical data for [Γ N ] 3 /N from Tab. 2 together with the regression lines, see Eq. (23).The red dashed lines are the analytical results for the y-intercept at ∆ = 0 and ∆ = 1/2, Sec. 4.

Table 2 :
The parameter [Γ N ] 3 as a function of N and ∆.The values of N are contained in the leftmost column and the values of ∆ are shown in the uppermost row.The number in the bracket is the estimated error in the last digit.