Multigrid preconditioning for the overlap operator in lattice QCD

Abstract

The overlap operator is a lattice discretization of the Dirac operator of quantum chromodynamics (QCD), the fundamental physical theory of the strong interaction between the quarks. As opposed to other discretizations, it preserves the important physical property of chiral symmetry, at the expense of requiring much more effort when solving systems posed with this operator. We present a preconditioning technique based on another lattice discretization, the Wilson-Dirac operator. The mathematical analysis precisely describes the effect of this preconditioning strategy in the case that the Wilson-Dirac operator is normal. Although this is not exactly the case in realistic settings, we show that current smearing techniques indeed drive the Wilson-Dirac operator towards normality, thus providing motivation for why our preconditioner works well in practice. Results of numerical experiments in physically relevant settings show that our preconditioning yields accelerations of more than an order of magnitude compared to unpreconditioned solvers.

Apppendix: The entries of \(D_W^HD_W-D_WD_W^H\)

To prove Proposition 2, we inspect the entries of $D_W^HD_W-D_WD_W^H$. We use the notation $\pi _\mu ^{\pm }$ for the matrices

$$\begin{aligned} \pi _\mu ^{\pm } = \tfrac{1}{2}(I_4 \pm \gamma _\mu ), \; \mu =0,\ldots ,3. \end{aligned}$$

The relations (3) between the $\gamma$-matrices show that each $\pi _\mu ^\pm$ is a projection and that, in addition,

$$\begin{aligned} \pi _\mu ^{+}\pi _\mu ^{-}= \pi _\mu ^{-}\pi _\mu ^{+} = 0, \; \mu =0,\ldots ,3. \end{aligned}$$

(28)

Considering all 12 variables at each lattice site as an entity, the graph associated with the nearest neighbor coupling represented by the matrix $D_W$ is the 4d-torus, and similarly for $D_W^H$. Table 4 repeats (6), giving the nonzero entries of a (block) row in $D_W$ and $D_W^H$ in terms of the $12 \times 12$ matrices that couple lattice site x with sites x and $x\pm \hat{\mu }$. We use m to denote $m_0+4$, with $m_0$ from (5).

Table 4 Coupling terms in $D_W$ and $D_W^H$

The product $D_W^HD_W$ represents a coupling between nearest and next-to-nearest lattice sites; the coupling $12 \times 12$ matrices are obtained as the sum over all paths of length two on the torus of the product of the respective coupling matrices in $D_W^H$ and $D_W$. A similar observation holds for $D_WD_W^H$. Table 5 reports all the entries of $D_W^HD_W$, and we now shortly discuss all the paths of length 2 that contribute to these entries of $D_W^HD_W$.

For the diagonal position (x, x), we have 9 paths of length 2, $(x,x) \rightarrow (x,x) \rightarrow (x,x)$ and $(x,x) \rightarrow (x,x\pm \hat{\mu }) \rightarrow (x,x), \mu = 0,\ldots ,3$. For a nearest neighbor $(x,x + \hat{\mu })$, we have the two paths $(x,x) \rightarrow (x,x) \rightarrow (x,x + \hat{\mu })$ and $(x,x) \rightarrow (x,x+\hat{\mu }) \rightarrow (x,x+\hat{\mu })$ and similarly in the negative directions. For a position $(x,x \pm 2\hat{\mu })$, there is only one path, $(x,x) \rightarrow (x,x \pm \hat{\mu }) \rightarrow (x,x \pm 2\hat{\mu })$, with the product of the couplings being 0 due to (28). Finally, for the other next-to-nearest neighbors, we always have two paths; for example, $(x,x) \rightarrow (x,x+\hat{\mu }) \rightarrow (x + \hat{\mu } - \hat{\nu })$ and $(x,x) \rightarrow (x,x-\hat{\nu }) \rightarrow (x + \hat{\mu } - \hat{\nu })$.

The coupling terms in $D_WD_W^H$ are identical to those for $D_W^HD_W$ except that we have to interchange all $\pi _\mu ^+$ and $\pi _\mu ^-$ as well as all $\pi _\nu ^+$ and $\pi _\nu ^-$.

Table 5 Coupling terms in $D_W^HD_W$

This shows that, in $D_W^HD_W-D_WD_W^H$, the only non-vanishing coupling terms are those at positions $(x,x+\hat{\mu } + \hat{\nu })$, $(x,x+\hat{\mu } - \hat{\nu })$ and $(x,x-\hat{\mu } - \hat{\nu })$ for $\mu \ne \nu$. They are given in Table 6, where we used the identities

$$\begin{aligned} \begin{array}{lll} \pi _{\mu }^{-}\pi _{\nu }^{-} - \pi _{\mu }^{+}\pi _{\nu }^{+} &{}=&{} \frac{1}{2}\left( - \gamma _{\mu } - \gamma _{\nu }\right) ,\\ \pi _{\mu }^{+}\pi _{\nu }^{-} - \pi _{\mu }^{-}\pi _{\nu }^{+} &{}=&{} \frac{1}{2}\left( \gamma _{\mu } - \gamma _{\nu }\right) ,\\ \pi _{\mu }^{-}\pi _{\nu }^{+} - \pi _{\mu }^{+}\pi _{\nu }^{-} &{}=&{} \frac{1}{2}\left( -\gamma _{\mu } + \gamma _{\nu }\right) ,\\ \pi _{\mu }^{+}\pi _{\nu }^{+} - \pi _{\mu }^{-}\pi _{\nu }^{-} &{}=&{} \frac{1}{2}\left( \gamma _{\mu } + \gamma _{\nu }\right) .\\ \end{array} \end{aligned}$$

By rearranging the terms, we obtain the plaquettes from (20) and (21). We exemplify this for position $(x,x+\hat{\mu } + \hat{\nu })$:

$$\begin{aligned}&\pi _\mu ^-\pi _\nu ^+ \otimes U_\mu (x)U_\nu (x+\hat{\mu }) + \pi _\nu ^-\pi _\mu ^+ \otimes U_\nu (x)U_\mu (x+\hat{\nu }) \\&\qquad - \left( \pi _\mu ^+\pi _\nu ^- \otimes U_\mu (x)U_\nu (x+\hat{\mu }) + \pi _\nu ^+\pi _\mu ^- \otimes U_\nu (x)U_\mu (x+\hat{\nu }) \right) \\&\quad = \tfrac{1}{2}( - \gamma _\mu + \gamma _\nu ) \otimes U_\mu (x)U_\nu (x+\hat{\mu }) + \tfrac{1}{2}(\gamma _\mu - \gamma _\nu ) \otimes U_\nu (x)U_\mu (x+\hat{\nu }) \\&\quad = \tfrac{1}{2}(-\gamma _\mu + \gamma _\nu ) \otimes \left( I_3 - Q_x^{\mu ,\nu } \right) U_\mu (x)U_\nu (x+\hat{\mu }) . \end{aligned}$$

Table 6 Coupling terms in $D_W^HD_W-D_WD_W^H$

Using the fact that, for the Frobenius norm, we have

$$\begin{aligned} \Vert AQ \Vert _F= & {} \Vert A\Vert _F \hbox { whenever}\ Q \ \hbox {is unitary (and}\ AQ \ \hbox {is defined)}, \\ \Vert A \otimes B \Vert _F= & {} \Vert A \Vert _F \cdot \Vert B \Vert _F \hbox { for all}\ A,B, \end{aligned}$$

we obtain the following for the squares of the Frobenius norms of all the coupling matrices from Table 6:

$$\begin{aligned} \begin{array}{ll} 2 \Vert I - Q_x^{\mu ,\nu }\Vert _F^2 &{} \text{ for } \text{ position } (x,x+\hat{\mu }+ \hat{\nu }), \\ 2 \Vert I-Q_x^{\mu ,-\nu } \Vert _F^2 &{} \text{ for } \text{ position } (x,x+\hat{\mu }-\hat{\nu }),\\ 2 \Vert I - Q_x^{-\mu ,-\nu } \Vert _F^2 &{} \text{ for } \text{ position } (x,x-\hat{\mu }-\hat{\nu }). \end{array} \end{aligned}$$

Finally, for any unitary matrix Q, we have

$$\begin{aligned} \Vert I -Q \Vert _F^2 = \mathrm {tr}\left( (I-Q^H)(I-Q)\right) = 2 \cdot \mathrm {Re}(\mathrm {tr}\left( I-Q\right) ). \end{aligned}$$

We now obtain $\Vert D_W^HD_W-D_WD_W^H\Vert _F^2$ by summing the squares of the Frobenius norms of all coupling matrices. This is a sum over 24n coupling matrices, n being the number of lattice sites. As discussed before, groups of four of these coupling matrices refer to the same plaquette $Q_x^{\mu ,\nu }$ up to conjugation in SU(3), so $\mathrm {tr}\left( I-Q\right)$ is the same for these four plaquettes Q. We can thus “normalize” to only consider all possible “first quadrant” plaquettes $Q_x^{\mu ,\nu }$ and obtain

$$\begin{aligned} \Vert D_W^HD_W - D_WD_W^H \Vert _F^2 = 4 \sum _x \sum _{\mu < \nu } 2 \cdot 2 \cdot \mathrm {Re}(\mathrm {tr}\left( I-Q_x^{\mu ,\nu }\right) ). \end{aligned}$$