discretization.md

Gauge-invariant discretization

We do the discretization of the Usadel equation as follows. We first discretize the action S. The discrete saddle-point equation or expressions for currents are not derived analytically, but instead obtained from the computer implementation of the action via automatic differentiation (AD). This simplifies things somewhat, because proper handling of the gauge fields is simpler to do in the action formulation.

The action is

  S[Q] = \frac{i\pi\nu}{8} \mathrm{Tr}[D(\hat{\nabla}Q)^2 + \alpha F_{ij}Q\hat{\nabla}_iQ\hat{\nabla}_jQ + 4 i \Omega Q]
  \,,

where \hat{\nabla}_j Q = \partial_j Q - i[\mathcal{A}_j, Q] and \Omega = \epsilon\tau_3 + i\Delta, \Delta^\dagger = \Delta. We'd like the discretized functional to be also gauge invariant.

We subdivide the space to cells, centered on a rectangular lattice \vec{r}_j=(h j_x, h j_y, h j_z) with j_{x,y,z}\in\mathbb{Z} and h are the lattice spacings. See the figure below:

We choose Q(j) = Q(\vec{r}_j) to be values of Q at the lattice sites. We define the Wilson link matrices U_{ij} =U_{ji}^{-1} =U(\vec{r}_{i},\vec{r}_{j}) = \mathrm{Pexp}[i\int_{L(\vec{r}_i,\vec{r}_j)}d\vec{r}'\,\cdot\vec{\mathcal{A}}(\vec{r}')] where L(\vec{r}_i,\vec{r}_j) is straight line from \vec{r}_j to \vec{r}_i and \mathrm{Pexp} is the path-ordered integral. Since the gauge field \vec{\mathcal{A}} is fixed, the link matrices can be computed ahead of time.

We want a discretized theory that is invariant under the discrete gauge transformation

Q(j) \mapsto u(j) Q u(j)^{-1} \,,

and

  U_{ij}
  \mapsto
  u(i)U_{ij}u(j)^{-1}
  \,,

and in the continuum limit reduces to the original theory. Such discretization has the advantage of being insensitive to the gauge choice, and having exact discrete conservation laws associated with the symmetry.

Derivative term

We discretize the derivatives a

  D_{ij}
  :=
  \frac{1}{h}[Q(i) U_{ij} - U_{ij} Q(j)]
  \,,
  \qquad
  \Rightarrow
  Q(i) D_{ij} = - D_{ij} Q(j)
  \,,

which satisfies an anticommutation relation analogous to continuum derivative.

We then discretize

S_2
=
\mathrm{Tr} D(\hat{\nabla} Q)^2
\mapsto
-D \frac{2d}{|\mathrm{neigh}|} h^d \sum_{\mathrm{neigh}(i,j)}\mathrm{Tr} D_{ij} D_{ji}

where d is the space dimension and |neigh| the number of neighbors on the grid.

Expanding in small h we have

S_2
\simeq
\frac{h^{d}d}{|\mathrm{neigh}|}
\sum_{i}\sum_{j\in{}\mathrm{neigh}(i)}\mathrm{Tr}\Bigl(\frac{Q(i) - Q(j)}{h} - [i\mathcal{A}(\frac{\vec{r}_i+\vec{r}_j}{2}), Q(j)]\Bigr)^2
\,,

which verifies the reduction to the continuum limit.

Field strength

The link matrices have the usual continuum limit expansion

U_{ji}
=
\sum_{n=0}^\infty i^n |\vec{r}_i-\vec{r}_j|^n \int_{-1/2}^{1/2}ds_1\int_{-1/2}^{s_1}ds_2\ldots\int_{-1/2}^{s_{n-1}}ds_n\, \mathcal{A}(s_1)\cdots\mathcal{A}(s_n)

where \mathcal{A}(s)=\mathcal{A}[(\frac{1}{2} - s)\vec{r}_i + (\frac{1}{2} + s)\vec{r}_j].

As well known in lattice QCD, the field strength can be expressed in terms of the link matrices U via a plaquette loop. Consider then the plaquette (see Fig. 1), with edges j\leftarrow{}i and k\leftarrow{}l in direction \mu and l\leftarrow{}i, k\leftarrow{}j in direction \nu, and expand around \vec{r}_0=\frac{\vec{r}_i+\vec{r}_j+\vec{r}_k+\vec{r}_l}{4}:

  P_{lkji}
  =
  U_{il}U_{lk}U_{kj}U_{ji}
  \,,

  \Rightarrow
  P_{lkji}
  \simeq
  1
  - i h^2
  \bigl(
  \partial_\mu \mathcal{A}_\nu - \partial_\nu \mathcal{A}_\mu - i[\mathcal{A}_\mu,\mathcal{A}_\nu]
  \bigr)
  +
  \mathcal{O}(h^3)
  =
  1
  + i h^2 F_{\mu \nu}(\vec{r}_0)
  +
  \mathcal{O}(h^3)

which then allows expressing F_{ij} in terms of the link matrices.

Hall term

In the plaquette of Fig. 1 we have

  \mathrm{tr}
  F_{\mu\nu}Q \hat{\nabla}_\mu Q \hat{\nabla}_\nu Q
  \rvert_{\mu=ji, \nu=li}
  \simeq
  \frac{i}{h^4}
  \mathrm{tr}
  (1 - P_{lkji}) U_{ij} D_{ji} Q(i) D_{il} U_{li}

  =
  \frac{-i}{h^4}
  \mathrm{tr}
  (U_{lk} U_{kj} - U_{li} U_{ij}) D_{ji} Q(i) D_{il}

So the structure is \frac{-i}{h^4} \times (gauge factors for counterclockwise loop - return back clockwise) \times DQD for three sites in counterclockwise order. Since only Q(i), Q(j), Q(l) appear here, one can also think of this as a plaquette corner with Q(i) being the corner.

To do the sum over the \mu, \nu indices, one can average over corners of the plaquettes as follows, which gives the final discretization of the Hall term:

  S_H
  =
  \alpha
  \mathrm{Tr} F_{\mu\nu}Q\hat{\nabla}_\mu Q\hat{\nabla}_\nu Q

  S_H
  \mapsto
  \alpha
  \frac{h^d d}{4ih^4}
  \sum_{\mathrm{plaqc}(ijkl; i)}
  \mathrm{Tr}
  \Bigl(
  (U_{lk} U_{kj} - U_{li} U_{ij}) D_{ji} Q(i) D_{il}
  \Bigr)

  S_H
  =
  \frac{\alpha h^{d-4} d}{4i}
  \sum_{\mathrm{plaqc}(i; jkl)}
  \mathrm{Tr}
  \Bigl(
  (P_{lkji} - 1)(Q(i) - U_{ij} Q(j) U_{ji} Q(i) U_{il} Q(l) U_{li})
  \Bigr)
  \,,

where \mathrm{plaqc}(i;jkl) are the plaquette corners. The last line follows with direct calculation, and considering the cyclic permutations.

Matrix current

We can also calculate the matrix current j\mapsto{}i between neighboring cells:

  J_{ij}^a
  =
  \frac{\partial}{\partial\xi}
  S[
    U_{ij}
    \mapsto
    U_{ij}
    +
    T^a U_{ij}\xi
  ]\rvert_{\xi=0}

where T^a is an appropriate matrix generator, and the transformation is inserted to only one of the links.

The local gauge invariance now implies that there is an exact discrete continuity equation

  \sum_{j\in{}N_i} J_{ij}^a
  =
  \frac{\partial}{\partial\xi}
  S[Q(i)\mapsto{}e^{\xi T^a}Q(i)e^{-\xi T^a}]\rvert_{\xi=0}
  =
  \frac{\partial}{\partial\xi}
  S[\Omega(i)\mapsto{}e^{-\xi T^a}\Omega(i)e^{\xi T^a}]\rvert_{\xi=0}
  =
  R_{i}^a
  \,.

This is one main advantage of the procedure, in addition to its invariance vs. gauge fixing.

However, the definition of the matrix current is asymmetric, so that J_{ij}\ne{}-J_{ji} as they differ by a quantity of order \mathcal{O}(h^2). The asymmetry arises when T^a does not commute with U_{ij}, so it is an issue only for the SU(2) component. The problem is only in the interpretation.

The currents measured at sites i and j are different: since the parallel transport of spin between these two locations can imply rotation due to the gauge field, there's no reason why we should have J_{ij}=-J_{ji}, except in the continuum limit where the two points become close to each other.

Implementation

The above is sufficient for numerical implementation: it is not necessary to find out the saddle point equation symbolically. It can be deduced by parametrizing first so that Q^2=1 and then using automatic differentiation (AD). Only a routine evaluating the value of the action is then needed, in terms of the AD dual variables.

When looking for the saddle point in equilibrium, the problem becomes even simpler; Q has then additional restrictions, and we are looking for the minimum of the free energy S.