Qlua: Scripting Lattice QCD

• Qlua
• Manual
• Conventions
• Running qlua
• Core language
• Complex numbers
• Dirac algebra
• Arithmetic
• The lattice
• Lattice fields
• Sequential types
• Vectors and matrices
• Input and output
• Level III solvers
• Contractions
• Lattice matrices
• Numerical helpers
• Timing and memory
• Optional back ends
• Libraries
• LHPC library
• Library names
• Index

Lattice fields

A lattice field holds one value at every site: an integer, a real, a complex number, a color vector, a color matrix, a Dirac spinor or a full propagator. There are seven of them, and they are built the same way, so what follows describes the shape once and then only what each one adds.

Three things hold for all seven.

A field belongs to exactly one lattice and never mixes with a field on another; F.lattice is the one it belongs to. A field is created zeroed unless its constructor is given something to fill it with. And every operation on a field is confined to the lattice's current subset, so the same line of code touches every site inside L:Subset("all"):where(...) and half of them inside L:Subset("even"):where(...).

The data is spread over the ranks, one piece each, but a script never sees that: a reduction returns the same number everywhere, and reading one site gives the same value everywhere. What it costs is that these are collective operations — see Parallel execution.

Every field type has these:

set and set_zero are the only field operations that modify their receiver; everything else returns a new field. Both return the receiver, so they chain: psi:set_zero():norm2() is 0, and U[1]:set(V) gives back U[1].

Shift directions

F:shift(mu, dir) moves the field one site along axis mu. The direction is a string, and there are two ways of saying each of the two shifts: one names where the data comes from, the other where it goes to.

The lattice wraps, so a site at the end of an axis takes from, or sends to, the site at the beginning.

Reading and writing one site

F[{x0, x1, ...}] reads the value at a site. The coordinate is a list with one entry per axis, counted from 0, and the result is a plain number for Int and Real and the matching sequential object otherwise. The read is collective: the rank that owns the site finds the value and broadcasts it, so every rank gets the same answer.

F[{x0, x1, ...}] = v writes it. The write happens on the owning rank, and v must be of the matching type.

Single-site access is convenient and slow — a broadcast per call. It belongs in a test or a diagnostic, not in a loop over the volume.

L:Int

L:Int()                                           --> L:Int
L:Int(k)                                          --> L:Int
L:Int(I)                                          --> L:Int

A field of integers: zero everywhere, the constant k everywhere, or a copy of another Int field.

Beyond the common members, L:Int has the bitwise and ordering operations, each taking a second Int field on the same lattice:

and and or are Lua keywords, so those two cannot be written with the method colon. Reach them through the index instead:

local both = I["and"](I, I2)
local some = I["or"](I, I2)
local diff = I:xor(I2)                            -- xor is not a keyword

I:sum() returns a number; I:sum(ms) a vector.real.

L:Real

L:Real()                                          --> L:Real
L:Real(x)                                         --> L:Real
L:Real(I)                                         --> L:Real
L:Real(R)                                         --> L:Real

A field of doubles, from nothing, a number, an Int field or another Real field.

The elementary functions are applied site by site. All of them return a new Real field unless the table says otherwise.

expi is the quickest way to turn a phase into a field. Watch the result types: ceil and floor stay in Real while round and trunc cross to Int, so R:floor() and R:trunc() agree in value on a positive field and differ in type. sign has no zero: it answers +1 at a site where the field is exactly zero.

Two more take a second Real field:

And two reduce the field to a single number:

Unlike sum and norm2, those two ignore the current subset: they look at every site whatever subset is active.

L:Complex

L:Complex()                                       --> L:Complex
L:Complex(x)                                      --> L:Complex
L:Complex(z)                                      --> L:Complex
L:Complex(I | R | C)                              --> L:Complex
L:Complex(Rre, Rim)                               --> L:Complex

A field of complex numbers. The two-field form takes the real and imaginary parts as separate Real fields.

sqrt and log take the same branch cuts as the scalar complex methods of the same name.

L:Complex:project

C:project(C2, ms)                                 --> vector.complex

Momentum projection in one call: Σ_x conj(C(x))·C2(x) accumulated into bins chosen by the multiset ms. The same answer as qcd.dot(C, C2):sum(ms), without building the intermediate field.

L:ColorVector

L:ColorVector()                                   -- zero, default colors
L:ColorVector(v)                                  -- constant fill from a color.vector
L:ColorVector(V)                                  -- copy of another field
L:ColorVector(C, {c = k})                         -- color k set to the Complex field C
L:ColorVectorN(nc, ...)                           -- the same, with nc named

A color vector at every site. The plain form uses the lattice's default number of colors; the N form takes the count as its first argument and shifts the rest along.

L:ColorMatrix

L:ColorMatrix()                                   -- zero, default colors
L:ColorMatrix(r)                                  -- r times the identity
L:ColorMatrix(z)                                  -- z times the identity
L:ColorMatrix(m)                                  -- constant fill from a color.matrix
L:ColorMatrix(M)                                  -- copy of another field
L:ColorMatrix(C, {a = i, b = j})                  -- entry (i,j) set to C, rest zero
L:ColorMatrix(V, {b = j})                         -- column j set to V, rest zero
L:ColorMatrix(V1, V2)                             -- the outer product V1[i] * conj(V2[j])
L:ColorMatrixN(nc, ...)                           -- the same, with nc named

A color matrix at every site: an Nc×Nc complex matrix acting on color vectors. The N form takes the number of colors as its first argument and shifts the rest along.

Its two indices are named rather than positional. a is the left index, the row — the one a matrix contracts on when it acts from the left. b is the right index, the column. Both count from 0, so M[{a=0, b=2}] is the entry in the first row and third column, and M * V sums over b:

(M * V)[{c = a}]  =  Σ_b  M[{a = a, b = b}] * V[{c = b}]

proj iterates towards the nearest special unitary matrix, stopping when the change falls below eps or after imax iterations. It is defined for two and three colors; for a general n it is not.

L:DiracFermion

L:DiracFermion()                                  -- zero, default colors
L:DiracFermion(d)                                 -- constant fill from a dirac.fermion
L:DiracFermion(D)                                 -- copy of another field
L:DiracFermion(C, {c = i, d = j})                 -- entry (i,j) set to C, rest zero
L:DiracFermion(V, {d = j})                        -- spin j set to the ColorVector V
L:DiracFermionN(nc, ...)                          -- the same, with nc named

A Dirac spinor at every site: Ns spin components, each a color vector, so Nc×Ns complex numbers in all. The N form takes the number of colors as its first argument and shifts the rest along.

Its two indices are named. c is the color index, running over the Nc colors; d is the Dirac index, running over the qcd.Ns spin components. Both count from 0. Giving only d picks a whole color vector; giving both picks one number.

D:gamma{mu=m} and gamma{mu=m} * D are the same thing; the method form takes the same {mu=}/{n=} table the gamma constructor does.

L:DiracPropagator

L:DiracPropagator()                               -- zero, default colors
L:DiracPropagator(p)                              -- constant fill from a dirac.propagator
L:DiracPropagator(P)                              -- copy of another field
L:DiracPropagator(M)                              -- P[.,k,.,k] = M for every spin k
L:DiracPropagator(D, {c = i, d = j})              -- column (i,j) set to D
L:DiracPropagatorN(nc)                            -- zero, nc colors
L:DiracPropagatorN(nc, D, {c = i, d = j})         -- column (i,j) set to D

A full propagator at every site: a matrix on the combined color⊗spin space, so (Nc·Ns)² complex numbers, carrying four indices. The N form takes the number of colors as its first argument.

Here c and d are both left indices — the color and the Dirac index of the row, the pair that names which DiracFermion column is being selected. The right pair, the source color and spin, is not addressable from Lua: there is no key for it. P[{c=k, d=s}] therefore gives a whole DiracFermion, the column of P that a point source at color k and spin s produced, and assigning a DiracFermion there writes that column.

for c = 0, 2 do
   for d = 0, qcd.Ns - 1 do
      P[{c=c, d=d}] = solve(src[{c=c, d=d}])      -- one column per solve
   end
end

Assigning a DiracFermion to P[{c=k,d=s}] writes that column.

Quark contractions

qcd.quarkContract12

qcd.quarkContract12(P1, P2)                       --> DiracPropagator

The ε-tensor contraction that builds baryon blocks: the color indices of P1 and P2 are contracted against two Levi-Civita symbols, and one spin index of each is contracted with the other's.

There are six, differing only in which two of the four spin indices are summed over. Number them across the pair of arguments:

The digits in the name are the two that are contracted; the remaining two become the result's left and right spin indices, in that order.

All six take the same arguments and give the same type. Both arguments must be DiracPropagators over three colors — the ε-tensor fixes that — and they must be of the same kind: two lattice fields give a lattice field, two sequential values a sequential one. Sequential arguments need 2026.10d or later.

Random number fields

L:RandomState

L:RandomState(seed, I)                            --> L:RandomState

A field of independent random-number generators, one per site. seed is a number and I is a lattice Int field giving each site its own stream — the usual choice is a unique index per site, so that the result does not depend on the layout:

local function site_index(L)
  local x = L:Int(0)
  for mu = 0, #L - 1 do x = x * L[mu] + L:pcoord(mu) end
  return x
end
local S = L:RandomState(1234, site_index(L))

The gaussian generators draw each component from a unit normal. A random state is itself a QDP field, so it can be written to and read from a QIO file (Writer:RandomState).

Shifts, scatters and gathers

F:shift(mu, dir) moves a field by one site along an axis. Arbitrary communication — permutations, restrictions, projections onto a smaller lattice — goes through qcd.scatter or qcd.gather, which are built once from an index map and then applied as often as needed.

qcd.shift

qcd.shift(F, mu, dir)                             --> field

The same as F:shift(mu, dir), as a function. Accepts the field types the method does, and in addition a table of fields, which it shifts one by one.

qcd.scatter

qcd.scatter(dstLattice, srcLattice, {I0, I1, ...})   --> scatter

Builds a scatter. The index list has one lattice Int field per axis of the source lattice, each defined on the destination lattice: at every destination site they give the source site to read from. Indices outside the source lattice are rejected when the scatter is built, not when it is used.

local s = qcd.scatter(Lsmall, Lbig, idx)
local small = s(big)

A scatter accepts Int, Real, Complex, ColorVector, ColorMatrix, DiracFermion and DiracPropagator fields on its source lattice.

qcd.gather

qcd.gather(dstLattice, srcLattice, {I0, I1, ...})    --> gather

The transpose of a scatter: the index list has one Int field per axis of the destination lattice, each defined on the source lattice, naming where each source site contributes. Several sources may land on the same destination, so a gather is applied through a reduction rather than by calling it.

The arithmetic reductions take the numeric and color field types; the bitwise ones take Int only.

A gather's reductions combine the contributions in whatever order they arrive, which depends on the layout and on message timing. For and, or, xor, min and max that does not matter — they are associative and exact. For add and mul on floating-point data it does: the same gather applied to the same field can give answers differing in the last bits from run to run, and from one rank count to another. Do not compare two such results for exact equality, and do not expect a job to reproduce its own output bit for bit across layouts.