Plaquette

A whole program: build gauge fields of varying randomness and measure them.

Contents

Now we have enough pieces to put together plaquette computation. Code pieces below make up a whole program that makes several gauge fields of varying randomness; for each of them we compute and print the plaquette value. The code pieces are interspersed with comments about the computation.

Standard libraries

require "stdlib";

We need printf() and a few other functions commonly used in QCD code. We will discuss them where they are used.

In the real QCD application, one will also say require "qcd", but in this tutorial we are dissecting most of them, so we do not require the canned routines.

Lattice definition

To be able to run this tutorial on a laptop, let us keep the lattice size small:

L = qcd.lattice{4,4,4,8};
volume = 1;
for i = 0, #L - 1 do
  volume = volume * L[i];
end

Lattice volume is a handy constant to have in the QCD code.

Random State initialization

do
  local r = L:Int(os.random());
  for i = 0, #L - 1 do
    r = r * L[i] + L:pcoord(i);
  end
  S = L:RandomState(os.random(), r);
end

Here do ... end block forms a lexical scope for local name r. Once the block is exited, r loses its binding. This is how Lua avoids name space pollution. We are not worried about reproducibility of the results, so we use true random numbers to initialize the lattice random state S.

Making an SU(n) matrix

In the following, we need to have an SU(3) matrix. One can easily construct one from an arbitrary matrix with the following function

function makeSUn(X)
  local Y = (X - X:adjoin()) / 2;
  local one = X.lattice:ColorMatrix(complex(1,0));
  return (Y - Y:trace() / X:colors() * one):exp();
end

Here we first define Y to be the antihermitian part of X, then project the result to traceless matrices (these two steps project X to the Lie algebra “su(n)”), exponentiating the result gives us an element of the group.

Notice that the function asks the matrix for everything else it needs. X:colors() is its own number of colors, rather than 3 written out: one Qlua binary can carry several values of \(N_c\) at once — this one reports colors: 2 3 N — so the number of colors is a property of the object in hand, not of the installation, and every color-carrying type has a :colors() method. X.lattice is the lattice the field lives on, which every lattice field knows, so there is no lattice argument to pass in and none to get wrong.

(There is no qcd.Nc. Earlier versions of this tutorial used one; nothing in Qlua has ever defined it.)

This is what toSUn in the gauge library does. It is written out here rather than called, and under another name, so that defining it does not shadow the library's if you require "gauge" as well.

Weak SU(n) matrix

Let us first make a random \(SU(n)\) matrix close to the group identity.

function weak_SUn(S)
  return makeSUn(0.2 * S:gaussian_ColorMatrix());
end

This code should be self-explanatory by now. If it is not, you may want to reread previous tutorials.

Strong SU(n) matrix

To make a more random SU(n) matrix we start with multiplying several random \(SU(n)\) matrices from the vicinity of the group identity and project the result back to the group for good measure. This approach allows us to avoid generating a random point with the correct Haar measure, while giving a reasonable approximation to the uniform distribution on the group manifold.

function strong_SUn(S)
  local v = weak_SUn(S);
  for i = 1, 100 do
    v = v * weak_SUn(S);
  end
  return v:proj(1e-6, 150);
end

The meaning of parameters of v:proj() is somewhat obscure. Please consult Chroma sources for their significance. Qlua tries to implement :proj() to be identical to the corresponding Chroma's routine.

Note that :proj() exists only for \(N_c\) of 2 and 3. On a color matrix of any other size it is simply not there, and the error is bad index 'proj'.

strong_SUn is here to exercise the code below it, not as a way of making a gauge field. It is a random walk on the group — a hundred and one exponentials and a projection for every link — and nothing in it measures how close to uniform it has come. More to the point, a uniformly random field is not a gauge configuration: it is what one looks like at infinite coupling, with an average plaquette near zero, as the second test at the end of this tutorial shows. Real work reads a thermalized configuration from a file, which the tutorial on QIO files covers, or generates one.

Computing plaquettes

Computing the plaquette in the \(ij\) plane is short but it deserves an explanation. First, the code

function my_plaquette(U, i, j)
  local Ui = U[j+1]:shift(i, "from_forward");
  local Uj = U[i+1]:shift(j, "from_forward");
  local Pij = (U[i+1] * Ui *
                Uj:adjoin() * U[j+1]:adjoin()):trace():real();

  return Pij:sum() / (volume * U[1]:colors());
end
  • Here we struggle with two indexing schemes: Lua starts with 1, and QDP counts from 0. We prefer the zero-based approach, hence U[j+1] and U[i+1].
  • Definitions of Ui and Uj show how Qlua does nearest neighbor communication. Each lattice field has a method called shift() which takes two arguments: the shift axis (zero based), and the direction string. The direction string could be "from_forward", "from_backward", "to_forward", or "to_backward". There is no difference between "from_forward" and "to_backward", they simply reflect different ways of looking on things.
  • The definition of Pij shows that the string of multiplications produces a lattice object, to which one can apply the trace() method.
  • Finally, Pij:sum() sums the plaquette value across the lattice. It is nice to have it normalized as shown.

Running the plaquette code on a configuration

One more routine. Here we compute the plaquette in every plane and print each one as well as the average. The inner loop runs over every ordered pair of axes, so on a four-dimensional lattice each of the six planes is visited twice — harmless, since \(\hbox{Re\,tr\,}P_{ij}=\hbox{Re\,tr\,}P_{ji}\), but it is why twelve lines are printed and not six:

function run_test(name, U)
  printf("plaquette test case: %s\n", name);
  local s = 0;
  local k = 0;
  for i in interval(0, #L-1) do
    for j in skip(i, interval(0, #L-1)) do
      local sij = my_plaquette(U, i, j);
      printf("  plaq[%d, %d] = %20.10f\n", i, j, sij);
      s = s + sij;
      k = k + 1;
    end
  end
  printf("average plaquette %20.10f\n\n", s / k);
end

Noteworthy points:

  • Once again, to make code as generic as possible, k here explicitly counts the number of 2-planes on the lattice.
  • The loop over i is somewhat unusual. Lua calls constructs like these iterators. We have defined interval(n,m) in stdlib to return an iterator producing numbers n, n+1, ..., m.
  • The loop over j uses another iterator constructor. Here skip(n,x) takes a number and an iterator and produces a new iterator that runs though all values produced by x skipping n. Using the a combination of interval() and skip() helps in making sums over subsets of directions. Function bmw_hex() in file qlib/gauge.qlua is an instructive example.
  • We compute the ij plaquette twice (once in i,j and once in j,i). This is done only to illustrate the use of the skip() iterator. In the real code one would avoid extra work by iterating j in interval(i+1,#L-1) instead.

Tests

Here are two test cases.

Weak gauge field

All plaquettes should be close to 1 in this case:

W = {};
for i = 0, #L - 1 do
  W[i+1] = weak_SUn(S);
end
run_test("weak field", W);

Strong gauge field

Assuming the number of weak matrices used to construct a strong SU(n) matrix was large enough, all plaquettes should be close to 0 now:

U = {};
for i = 0, #L - 1 do
  U[i+1] = strong_SUn(S);
end
run_test("strong field", U);