The Stochastic Feynman-Hellmann Method

TL;DR

This paper introduces a stochastic extension of the Feynman-Hellmann method, enabling simultaneous computation of multiple operators and momenta, demonstrated through nucleon axial charge and form factors.

Contribution

It presents a novel stochastic approach to the Feynman-Hellmann method, overcoming previous limitations of single-operator and single-momentum restrictions.

Findings

Successfully reproduces nucleon axial charge

Demonstrates non-zero momentum transfer form factors

Validates the stochastic method's effectiveness

Abstract

The Feynman-Hellmann method, as implemented by Bouchard et al. [1612.06963], was recently employed successfully to determine the nucleon axial charge. A limitation of the method was the restriction to a single operator and a single momentum during the computation of each "Feynman- Hellmann" propagator. By using stochastic techniques to estimate the all-to-all propagator, we relax this constraint and demonstrate the successful implementation of this new method. We show reproduction of the axial charge on a test ensemble and non-zero momentum transfer points of the axial and vector form factors.

Equations

$$\frac{\partial E_n}{\partial \lambda }=<n|H_{\lambda }|n>$$

$$C_{\lambda }(t)=<\lambda |\mathcal{O}(t){\mathcal{O}}^{†}(0)|\lambda >$$

$${\mathcal{Z}}_{\lambda }$$

$$S_{\lambda }$$

$$\lambda j_{\mu }(x)=\bar{q}(x){\gamma }_{\mu }{\gamma }_{5}q(x)$$

$$\frac{\partial C_{\lambda }(t)}{\partial \lambda }=-\left[\frac{{\partial }_{\lambda }{\mathcal{Z}}_{\lambda }}{{\mathcal{Z}}_{\lambda }}{C}_{\lambda }(t)+\frac{1}{{\mathcal{Z}}_{\lambda }}\int \mathcal{D}\psi \mathcal{D}\bar{\psi }\mathcal{D}{A}_{\mu }{e}^{-(S+S_{\lambda })}\int d^4{x}^{\prime }j(x^{\prime })\mathcal{O}(t)\mathcal{O}(0)\right]$$

$$\frac{\partial C_{\lambda }(t)}{\partial \lambda }{|}_{\lambda =0}=C_{\lambda }(t)\int d^4{x}^{\prime }<\Omega |j(x^{\prime })|\Omega >-\int d{t}^{\prime }<\Omega |T\left[\mathcal{O}(t)J(t^{\prime })\mathcal{O}(0)\right]|\Omega >$$

$$\frac{\partial m_{\mathrm{eff}}(t,\tau )}{\partial \lambda }{|}_{\lambda =0}=\frac{1}{\tau }\left[\frac{-{\partial }_{\lambda }{C}_{\lambda }(t+\tau )}{C(t+\tau )}-\frac{-{\partial }_{\lambda }{C}_{\lambda }(t)}{C(t)}\right]\overrightarrow{t\to \infty }\frac{J_{00}}{2E_0}$$

$$S_{\Gamma }=\sum _{z}S(y,z)\Gamma (z)S(z,x)$$

$$D|\psi >$$

$$|\psi >$$

$$|\varphi >=D^{-1}\Gamma e^{iq\cdot x}|\psi >.$$

$$|{\varphi }^{\prime }>=\frac{1}{N}\sum _{i=1}^N{D}^{-1}|{\eta }_i><{\eta }_i|\Gamma e^{iq\cdot x}|\psi >$$

Full text

1 Introduction

The Feynman-Hellmann (FH) method connects matrix elements to variations in the spectrum.

[TABLE]

Above, $H$ is some perturbing Hamiltonian: $H=H_{0}+\lambda H_{\lambda}$ and $|n>$ is an eigenstate of $H_{0}$. A few variations of the FH method have been employed in Lattice QCD to compute hadronic matrix elements [1, 2, 3, 4, 5]. As pointed out in Ref. [6], these FH inspired methods are very similar to, or practically the same as those found in older publications [7, 8, 9, 10] and they can be rigorously connected to functional derivatives of the partition function through an application of the FH Theorem.

In the adaptation from [6], a two-point correlation function in the presence of an external field is considered

[TABLE]

where $j_{\mu}(x)$ is some current density.

[TABLE]

The partial derivative of $C_{\lambda}(t)$ yields an expression with two terms.

[TABLE]

The first is the vacuum matrix element, which is zero unless $j(x)$ has vacuum quantum numbers. The second term is an integral over all spacetime of the creation/annihilation operators and current. Setting $\lambda=0$ reduces $\mathcal{Z}_{\lambda}$ to the ordinary partition function in equation 5. Therefore, this method is employed in practice on an ensemble sampled with $\lambda=0$.

[TABLE]

Above, $J(t^{\prime})$ is a time-dependent current given by $J(t^{\prime})=\int d^{3}x^{\prime}\ j(t^{\prime},\overrightarrow{x}^{\prime})$. The region for which $0<t^{\prime}<t$, the second term is the hadronic matrix element of interest.

In analogy to the effective mass from a standard two-point function, the long-time limit of the partial derivative of the effective mass plateaus to the matrix element of interest.

[TABLE]

There is a beautiful symmetry between equations 1 and 7, hence the name, FH method. By studying $J(t^{\prime})$ in all regions, including where the current is inserted outside of the hadron, a functional fitting form can be derived for which the matrix element of interest is uniquely, linearly temporally-dependent. For a rigorous treatment of this derivation, we refer the reader to [6].

In practice, the derivative of the correlation function is computed by applying a current globally in spacetime to a quark propagator through a matrix product. The resulting object is then treated as a source and used as a right-hand side in solving the Dirac equation.

[TABLE]

We refer to $S_{\operatorname{\Gamma}}$ as the “Feynman-Hellman” propagator. In order to compute $\partial_{\lambda}C_{\lambda}$ the FH propagator intercepts ordinary quark lines in a hadron’s two-point function. This method, as numerically described above, was efficaciously used to calculate the axial coupling to within $1\%$ uncertainty [5].

There are several benefits to using this technique for computing correlation functions. Since the insertion time is summed over, there is only one time variable that may be varied, leading to potentially easier systematics to control. Additionally, for the cost of a single propagator, information from all timeslices is gained. Furthermore, $S_{\operatorname{\Gamma}}$ is reusable for a variety of different sink interpolators, including different choices of sink smearing and different hadrons altogether.

2 Stochastic Feynman-Hellmann

The FH method has a number of disadvantages compared to the more conventional sequential source approach to computing three-point functions [11, 12]. The current time is summed everywhere on the lattice, including over sites outside of the hadron and at the source/sink, so called contact terms. This yields extra terms that do not appear in other methods. In addition, explicit time dependence of the current is lost. This complicates the analysis for non-zero momentum transfer points, in which transition amplitudes appear due to differing initial and final state energies. Finally, each choice of $\Gamma$ and momentum injection costs an individual FH propagator, resulting in high computational costs to calculate form factors.

The aim of this work is to remedy two out of the three disadvantages. To achieve this, we employ a stochastic basis, whose vectors obey the property: $<\eta(i)\eta^{*}(j)>=\delta_{ij}$. If we consider one spin/color component of a propagator, the Dirac equation is solved as

[TABLE]

where $|\chi>$ is a spin/color component of the source and $|\psi>$ is of the propagator. A component of the FH propagator is

[TABLE]

Inserting an outer product of noise vectors factorizes the algorithm, allowing different matrix elements and momentum-transfer points to be computed without requiring additional solves.

[TABLE]

There are many possible choices for the type of noise basis, such as Gaussian, $\mathbb{Z}_{2}$, $\mathbb{Z}_{4}$, or some kind of dilution pattern. In the context of disconnected diagrams or distillation, the effectiveness of different types of noise and variance reduction techniques has been studied extensively [13, 14, 15, 16, 17, 18, 19, 20, 21, 22]. For this work we employ fully spin and color diluted hierarchical probing [19].

3 Numerical Results

The stochastic FH method is tested on a Möbius Domain Wall on gradient flowed HISQ mixed action. The lattice spacing, pion mass, and volume of the ensemble are $a=.12\ \operatorname{fm}$, $m_{\pi}=310\ \operatorname{MeV}$, and $m_{\pi}L=4.5$ respectively. Details of the action setup and ensemble can be found in [23]. Figures 1(a) and 1(b) show the calculation of $g_{A}$ and $g_{V}$ with regular and stochastic FH methods. Figure 2, 3(a), and 3(b) demonstrate the computation of other observables by reusing the same stochastic basis.

In these figures, the derivative of the effective mass is plotted by constructing the correlated difference defined in equation 7. In the long-time limit, this quantity should plateau to the matrix element. From Figures 1(a) and 1(b), it’s clear that the stochastic basis successfully reproduces the “exact” answer, albeit with larger statistical uncertainty. In the case of Figure 2, there is no exact calculation to compare to, however the signal looks clean and plateaus to a reasonable value of the tensor charge [24].

Figures 3(a) and 3(b) illustrate $g_{A}$ and $gV$ with one unit of lattice momentum. This is done by fixing the sink and current with one unit of momentum in the same direction. It should be noted however that with the stochastic FH method, arbitrary momenta may be chosen for both the current insertion and sink. An analysis that takes advantage of different momentum combinations at the current and sink is forthcoming. The statistical uncertainty for non-zero momentum is large; enhancements to the basis are also currently being explored.

4 Conclusion and Outlook

We study an algorithm that computes the FH propagator with a stochastic basis, extending the original FH method to give arbitrary observables and non-zero momentum transfer points in a feasible amount of computing time. This is achieved through a stochastic noise basis. It should be noted that this basis is required to compute disconnected diagrams, therefore amortizing the cost for the connected calculation. Additionally, the current insertion time may be exposed by introducing a delta function in time to the bilinear, thereby recovering the full dependence of the three-point function. Finally, it has been shown in [21, 22] that singular value deflation plays a synergistic role with hierarchical probing in greatly reducing variance. Incorporating deflation into this method is being explored.

Acknowledgments

This work was supported by an award of computer time by the Lawrence Livermore National Laboratory (LLNL) Multiprogrammatic and Institutional Computing program through a Tier 1 Grand Challenge award. We would like to thank Balint Joo for general support with USQCD software and the larger CalLat collaboration for useful discussions. The work of ASG was performed under the auspices of the U.S. Department of Energy by LLNL under Contract No. DE-AC52-07NA27344. This research used the Sierra computer operated by the Lawrence Livermore National Laboratory for the Office of Advanced Simulation and Computing and Institutional Research and Development, NNSA Defense Programs within the U.S. Department of Energy.