Program package SF-MP-HMC: stability and scalability

We have developed a new code for Schrödinger functional simulations which is based on the mass preconditioned [65] Hybrid Monte Carlo algorithm and named it SF-MP-HMC. This implementation is based on the previously mentioned DD-HMC code of Martin Lüscher and our new implementation of the mass precon-ditioned HMC that was discussed in chapter 4. The new SF-MP-HMC code cludes SSE inline-assembly optimization routines and SSE memory prefetch in-structions, mainly based on the routines from DD-HMC and MP-HMC. The ge-ometry and routines related to the Sheikholeslami-Wohlert term as well as the used Conjugate Gradient (CG) solver are taken from the new correlation function code of the ALPHA Collaboration[96].

For practical reasons, the code is parallelized only in three (spatial) directions.

We first test the scalability of the SF-MP-HMC code starting from the thermalized N_f =2 configurations of 24⁴lattices. This test is done in a setup without the back-ground field, starting from the configurations and parameter sets obtained in the GHMC code studies of the quark mass renormalization from the Schroedinger

5.6 Algorithmic challenges with SF formulation

10 100 1000 10000

2 4 8 16 32 64

Time/Traj. [s]

N_Proc

X7560(uv)

100 1000 10000 100000

2 4 8 16 32 64 128 256

Time/Traj. [s]

N_Proc

X7560 (uv) X5570 (ice2)

Figure 5.3: Scaling of the SF-MP-HMC code on a 16⁴ lattice (up) and a 24⁴ lat-tice(down) for N_f = 4 flavors of dynamical fermions. The tests are performed on the cluster of 8 core Intel machines with two types of processors: X7560 and X5570 (HLRN, ZiB Berlin). The background field given in eq. 5.12 for the computation of the coupling is switched on.

functional scheme [20]. The result of the test is shown in Figure 5.2. A very good scaling up to the maximal number of processors for this lattice size can be observed. The geometry needed for running on 96 processors has a very large surface to bulk ratio for the local lattices, therefore a more intensive communica-tion overhead is observed and makes this choice of the number of processors not suitable for longer runs. Hence, optimal choice of the number of processors for this lattice size would be 128 or even 144. In Figure 5.3 we show the scaling of the code on 16⁴ and 24⁴ lattices with N_f = 4 flavors of massless quarks. Good scalability is also observed in these runs.

The Dirichlet boundary conditions imposed on the fermions in SF induce a fi-nite gap in the spectrum of the Dirac operator and serve as a natural infrared

cutoff. Hence, the practical problems of simulating Wilson fermions on the lat-tice, such as the stability issues discussed in section 4.2.2, do not occur in the SF simulations we are performing² here. To illustrate that, we show in Figure 5.4 the history of the energy violation in a long chain of the simulation on a 24⁴ lattice. Spikes similar to those from Figure 4.4 characterizing periodic boundary conditions with Wilson fermions are not present and the energy violation does not leave the rangeΔH ∈ [−2, 2].

-2.0 0.0 2.0

0 1000 2000 3000 4000 5000

Δ H

Figure 5.4: History of the energy violationΔH in the SF-MP-HMC run on the 24⁴ lattice(right). The plotted data correspond to roughly 5100 subsequent HMC trajectories of the lengthτ =2.0.

2The smallest eigenvalue for the free Dirac operator in SF, with vanishing boundary conditions, isλ²_min= (π/2T)². For the lattice sizes needed in SF simulations of the running coupling (up to T = L = 24) we are still protected with the discussed mass gap. Once the lattice size L reaches larger values, problems similar to those for light Wilson quarks on a periodic lattice may occur.

6 Running coupling in N _f = ^{4 theory}

We have previously defined the step scaling function in continuum theory. It is also possible to define the lattice step scaling function that is dependent on the lattice resolution and equals the continuum SSF up to the cutoff effects. Here we describe the strategy for the computation of the step scaling function from lattice simulations in the theory with four dynamical fermion flavors. Due to the peculiarities in the SF formulation of theO(^a)improved Wilson theory, the way we chose different improvement coefficients has to be discussed. We review here the running coupling and the Lambda parameter results from Ref. [21, 95] and discuss the strategy for improving it. In order to achieve better precision than the authors of [21, 95], it is important to make sure that the quark mass is tuned to zero accordingly with the aimed improved accuracy. We derive the precision criteria to be fulfilled such that the only remaining dependence in the coupling is the one on the system size. After these criteria are fulfilled, we are confident that the systematic errors are negligible in comparison to the statistical ones and it is possible to proceed towards the precise computation of the running coupling for the N_f =4 approximation of QCD. We then perform a continuum extrapolation and arrive to the continuum step scaling function. Finally, we give the value of the Lambda parameter of the theory in units of the system size.

6.1 Lattice step scaling function

To employ lattice simulations for the computation of the step scaling function (SSF) defined in section 5.4, one defines the discretized step scaling function Σ(u,L/a), which has an additional dependence on the lattice resolutiona/L

Σ(u,L/a) = ^g2(2L), u =^g2(^L), (6.1) with bare coupling g₀ fixed, L/a fixed and vanishing quark massm = 0. Let us recall that the boundary conditions for the fermions given in eq. 5.7 and 5.8 intro-duce a gap into the spectrum of the Dirac operator and allow for the Schrödinger functional lattice simulations to be performed for vanishing quark masses. We are hence able to supplement the definition of the running coupling by the re-quirement m = 0. Therefore, the lattice SSF also remains independent of the quark mass. The lattice SSF gives the continuum SSF as the lattice spacing is sent

Figure 6.1: A sketch of the recursive finite-size scaling method used in the com-putation of the SF couplinggat different energy scalesμ = 1/L. The evolution of the coupling with respect to energy is computed in sev-eral steps, changing μ by a factor of 2 in each of the steps. The steps of keeping a fixed lattice resolutionaand increasingL →2L (horizon-tal direction in the figure) are alternated with the steps of keeping the lattice extension L fixed and decreasing the a (vertical direction). In this way, the non-perturbative renormalization group (expression 6.3) is implemented, such that it is not necessary to have large ratios of the scales and therefore the discretization errors are already kept small for L/a1. Illustration is taken over from [97].