Nuclear lattice simulations: Status and perspectives

I review the present status of nuclear lattice simulations. This talk is dedicated to the memory of Gerald E. Brown.


Chiral effective field theory on a lattice
The starting point of chiral EFT for nuclei is the standard effective pion-nucleon Lagrangian supplemented by multi-nucleon field operators, that faithfully reproduces the chiral Ward identities of QCD. The latter terms are not constrained by chiral symmetry and their parameters (low-energy constants, called LECs for brevity) must be determined from a fit to few-nucleon data. Following Weinberg, the power counting is performed on the level of the effective potential between two, three and four nucleons and a properly regularized Schrödinger-type equation is solved to generate the bound and scattering states. Life becomes much simpler if one goes to Euclidean time and formulates the theory on a space-time lattice of volume (N a) 3 × (N t a t ), with a (a t ) the (temporal) lattice spacing and N (N t ) the number of grid points in the space (time) direction. In fact, the lattice spacing a serves as an UV regulator of the EFT as it implies a maximal momentum of p max ∼ π/a. For a typical value of a ≃ 2 fm, one has p max ≃ 300 MeV, i.e. a very soft interaction. A further advantage of the lattice formulation is that is automatically takes into account all possible configurations allowed by symmetries, in particular, we can have up to four nucleons on one lattice site. Such configurations naturally lead to clustering, in fact, they require some special treatment so as not to produce overbinding of α-like configurations. For that reason, the leading order (LO) interaction is smeared with a Gaussian, where the width parameter is adjusted to give the proper average S-wave NN effective range. Such smeared sources are also used frequently in lattice QCD calculations of hadron properties and dynamics. This also means that the LO lattice EFT result should not be directly compared with a corresponding continuum calculation, where such effects are generated by higher order operators in the full potential. In the lattice approach, only the LO potential is summed up whereas all higher orders are treated in perturbation theory. This is consistent with the low momentum cutoff utilized. An important new development is the so-called nuclear lattice projection Monte Carlo (MC) technique [7]. We use a larger class of initial and final states than considered in our earlier work. I give some examples: For the calculation of 4 He we use an initial state with four nucleons, each at zero momentum. For the calculation of 8 Be we use the same initial state as for 4 He, but then apply creation operators after the first time step to inject four more nucleons at zero momentum. The analogous process is done to extract four nucleons before the last step. This injection and extraction process of nucleons at zero momentum Lattice results for the 12 C spectrum at leading order. In Panel I results from three different initial states, A, B, and ∆ are shown, each approaching the ground state energy. In Panel II results starting from three other initial states, C, D, and Λ are displayed. These trace out an intermediate plateau at energy about 7 MeV above the ground state. This is nothing but the celebrated Hoyle state. helps to eliminate directional biases caused by initial and final state momenta. To probe the structure of the 12 C states, we use not only delocalized initial and final states (configurations A,B,C,D) but also cluster-type wave functions (compact triangle ∆, elongated triangle Λ). As seen in Fig. 1, while the configurations A,B and ∆ approach the ground state energy, we do not find fast convergence to the ground state for configurations C,D and Λ, thus giving the structure of the first excited 0 + state. The common thread connecting each of the initial states C, D, and Λ, is that each produces a state which has an extended or prolate geometry. The projection MC method is thus a powerful tool to unravel the structure of any considered nucleus and its excited states.
Next, we must fix the parameters (LECs) of our theory. At present, we work at next-to-next-toleading order (NNLO) in the chiral expansion, including strong and electromagnetic isospin breaking (the Coulomb force is included as described in [8]). We have 9 strong isospin-symmetric LECs. We fix these from the S-and P-waves in np scattering and the deuteron quadrupole moment, Q d . The 2 LECs related to isospin-breaking are determined from a fit to the different np, pp and nn scattering lengths, whereas the two LECs C and D of the leading three-nucleon force are fixed from the nd doublet scattering length and the triton binding energy. The first non-trivial predictions are then the momentum dependence of the pp 1 S 0 partial wave, which agrees nicely with the Nijmegen PWA result [8] and the binding energy difference of the triton and 3 He. We find E( 3 He) − E( 3 H) = 0.78(5) MeV [9], in good agreement with the empirical value of 0.76 MeV.
The ground state energies for various nuclei are summarized in Tab. 1. Consider first 4 He, 8 Be and 12 C, which show some unique trends. We observe that the smeared LO contribution is already close to the empirical numbers for these nuclei, whereas the two-nucleon forces alone underbind at NNLO. The rescue comes from the three-nucleon forces at this order, which provide just the right amount of binding. These 3N forces become more important as the atomic number A increases, as is also observed in e.g. Greens function MC [10] or no-core shell model calculations [11]. The NNLO results for these nuclei are already quite satisfactory. Of course, higher orders and improved methods to extract the energies will have to be considered in the future. Work along these lines is in progress. The situation is somewhat different for 16 O (note that these are preliminary numbers which might still change a bit). Here, the smeared leading order leads to some overbinding, which is not completely corrected from the 3N forces at NNLO. I will come back to this issue in Sec. 7.

Spectrum and structure of carbon-12
The spectrum of 12 C is one of the most challenging topics in nuclear theory. More precisely, the socalled Hoyle state, a 0 + excitation about 7.7 MeV above the ground state, has been an enigma to ab initio calculations since it was predicted by Hoyle in 1954 [12] from a calculation of the abundance of carbon and oxygen through stellar fusion and experimentally verified at Caltech a few years later [13]. Before continuing, I have to define what is precisely meant by an ab initio calculation. In a first step, one fixes the underlying forces between two and three nucleons in few-body systems, i.e. bound and scattering observables with A ≤ 4. Then, given the so determined forces, one solves the quantum A-body problem exactly, in this case using lattice Monte Carlo methods 1 . Any nuclear property for systems with A ≥ 5 is then a prediction, the precision of the calculations can be systematically improved by going to higher order in the potential and improving the statistical and systematic errors of the MC signals. As with any new method, it is only accepted by the community if one is able to solve problems that could not be resolved before. That is why as a first application of the framework of nuclear lattice simulations we focused on the spectrum of 12 C. In fact, based on standing waves properly projected onto zero total momentum as well as angular momentum and parity, the first ever ab initio calculation of the Hoyle state was reported in [16]. Since then, using the projection MC method, we have improved the calculations and gained further insight into the structure of the ground and excited state, as summarized in Tab. 2. The trends for the excited states are similar to the ground state Table 2. Low-lying spectrum of 12 C at NNLO from 2N and 3N forces [7] in comparison to experiment.
−92. 16  -about 20% of binding from the 3N forces is required to bring agreement between experiment and theory. Of particular interest besides the Hoyle state is the 2 + excitation just 2 MeV above it. This state can be interpreted as an excitation of the elongated triangle that features prominently in the structure of the Hoyle state. The predicted energy is in good agreement with recent measurements, see Tab. 2 (cf. also the group-theoretical anlysis of such cluster states in Ref. [19]). In Ref. [7], electromagnetic observables (charge radii, quadrupole moments and electromagnetic transitions among the even-parity states of 12 C) were also given at LO -all these results are in reasonable agreement with experiment. However, the calculation of the higher order corrections turned out to be more complicated than the ones for the spectrum. This requires algorithmic improvements for a precise extraction of the corresponding matrix elements which are presently under investigation.

Life on Earth -an accident?
The Hoyle state has been heralded as a prime example of the anthropic principle 2 , which states that the parameters of the fundamental interactions take values that are compatible with the existence of carbon-based life and that the Universe is old enough for such life to have formed (for a recent review, see [21]). Computer simulations allow one to address the question of how sensitive the production of carbon and oxygen in hot stars is to the fundamental parameters of QCD+QED, the gauge field theories underlying the formation of atomic nuclei. This is one of the great strengths of the framework presented here, as I will outline in the following. In fact, the energy of the Hoyle state is just 379 keV away from the 3α threshold, and the production rate of carbon in the 3α process scales with exp(−∆E/k B T ), with ∆E = E Hoyle − 3E α , T is the stellar temperature and k B Boltzmann's constant. In earlier studies, it was shown that ∆E can be modified by ±100 keV, so that one still has a sufficient production of carbon and oxygen [22]. This does not appear particularly fine-tuned, as variation of about 25% are consistent with carbon-oxygen based life. However, we have to translate this constraint on ∆E into constraints on the fundamental parameters of QCD+QCD, here the average light quark mass m q and the fine-structure constant α QED , see the recent works for more details [23,24]. Note also that the strong coupling constant α S is fixed from scale-setting through the nucleon mass, so it INPC 2013 can not be varied independently. In the following, we only consider small changes in m q and α QED , so that δE i ≃ (∂E i /∂M π )| M ph π δM π + (∂E i /∂α QED ) | α ph QED δα QED , where the superscript 'ph' refers to the physical value and we use the Gell-Mann-Oakes-Renner relation M 2 π ∼ m q to translate the quark mass dependence into the pion mass dependence. The knowledge of the pion mass dependence of the nuclear Hamiltonian is required to calculate the pertinent derivatives. In fact, consider any nuclear energy level or difference under small perturbations of the pion mass around its empirical value to leading order in the chiral expansion, whereM π refers to the explicit M π -dependence from the pion propagator in the OPE contribution, and are obtained from the quark mass dependence of the nucleon mass and the nucleon axial-vector coupling combined with the one of the pion decay constant utilizing lattice QCD data and chiral perturbation theory. The coefficients x 3 , x 4 can be mapped uniquely on the quark mass dependence of the inverse nucleon-nucleon S-wave scattering lengths, denotedĀ s ≡ ∂a −1 s /∂M π | M ph π and A t ≡ ∂a −1 t /∂M π | M ph π . These can be determined from chiral nuclear EFT, assuming some modeling based on resonance saturation [25]. The most recent analysis at NNLO givesĀ s = 0.29 +0. 25 −0.23 , A t = −0.18 ± 0.10 [26]. This introduces some amount of uncertainty, as discussed below. The various partial derivatives in Eq. (1) can be straightforwardly and precisely computed using the auxiliary field quantum Monte Carlo techniques (for details, see Ref. [24]). Independent of the precise values of A s andĀ t , we find strong evidence that the physics of the 3α-process is driven by α-clustering, as speculated before (see e.g. [27]). The condition that |∆E| ≤ 100 keV can be turned into an inequality that connects the quark mass change δm q /m q withĀ s andĀ t , i.e. (c 1Ās + c 2Āt + c 3 ) · (δm q /m q ) ≤ |∆E| (where the c i are calculated using our MC techniques), as depicted in Fig. 2, the so-called end-of-theworld plot. From that figure one sees that shifts in the light quark mass at the ≃ 2 − 3% level are unlikely to be detrimental to the development of life. Tolerance against much larger changes cannot be ruled out at present, given the relatively limited knowledge of the quark mass dependence of the two-nucleon S-wave scattering parametersĀ s,t . Lattice QCD is expected to provide refined estimates of the scattering parameters in the future. Further, variations of the fine-structure constant α QED up to ±2.5% are consistent with the requirement of sufficient carbon and oxygen production in stars.  End-of-the-world plot. For relative quark mass changes of 0.5% (wide band), 1% (middle band) and 5% (inner most band), one can read off the allowed values forĀ s andĀ t (both varied within natural ranges from −1 to +1), that are consistent with the requirements [22] for carbon-oxygen based life. The recent determination ofĀ s andĀ t from Ref. [26] is depicted by the black cross.
Since the formalism outlined above does particularly well for α-cluster states (note, however, that it is not restricted to such type of nuclei), it is natural to ask what happens if one extends these simulations to the next few α-cluster nuclei, that are 16 O, 20 Ne, 24 Mg, and 28 Si. Here, I report some preliminary numbers for the ground state energies and the low-lying spectra of these nuclei. First, it is interesting to note that for this mass range the CPU time scales approximately as A 2 (whereas naively one would expect an A 3 -scaling from the determinant calculation), so that the required computing time for 28 Si is only 5.4 times the one for the 12 C ground state. Second, we find that these heavier nuclei are overbound at NNLO, with the overbinding increasing with mass number, cf. Fig. 3. Let us consider the 16 O nucleus. Using the same action as described before leads to E g.s. ( 16 O) = (−141 ± 1) MeV, which is about 14 MeV overbinding (for an earlier variational MC calculation based on non-chiral forces, see [28]). This effect can be traced back to the new α-cluster geometries (the square and the tetrahedral structure, to be precise) that are possible with four or more alpha particles. These geometries generate too much attraction, similar to the four-nucleon configurations on one lattice site that required Gaussian smearing as discussed above. The simplest method to overcome this effect is to introduce a smeared four-nucleon operator, which is formally of higher order in the chiral expansion (for an earlier use of such type of operator, see e.g. Refs. [8,9]). The strength of this operator can be fixed from the ground state energy of 16 O. As can be seen from Fig. 3, this also leads to a good description of the ground state energies of all α-cluster nuclei up to 28 Si. Therefore, we believe that we have captured correctly the physics behind the binding energies of these nuclei. Clearly, these studies need to be backed by higher order calculations. In fact, we are presently extending the chiral potentials to next-to-next-to-next-to-leading order (N 3 LO). These exist in the continuum (see e.g. the N 3 LO three-nucleon forces in Refs. [29][30][31]) and just need to be adopted to the lattice notation.  Including a smeared four-nucleon contact interaction with its strength fixed from the 16 24 Mg, and 28 Si, respectively, and the value of the experimental binding energy is given in the square brackets.
It is also interesting to study the low-lying spectrum of these nuclei. Consider again 16 O. The tetrahedral structure is the dominant component in its ground state wave function with J π = 0 + and by its geometry it naturally leads to a 3 − excitation. The square configuration is dominant in the first 0 + excitation, that is located just below the 3 + -state and has a 2 + configuration as its lowest excitation. These configurations are exactly matching the experimentally found low-lying states in 16 O (see also [32]). In fact, we are presently performing a detailed analysis of this spectrum and the influence of the breaking of rotational symmetry of such configurations on the coarse lattices with a ≃ 2 fm used.

Summary and outlook
To summarize, nuclear lattice simulations are a novel scheme to tackle the nuclear quantum manybody problem, firmly rooted in the so successful chiral effective field theory of the nuclear forces. To make further progress, methodological investigations in three directions are required. First, one has to address the topic of corrections beyond NNLO. Second, systematic studies of the lattice spacing dependence and the related issue of rotational symmetry breaking are required. Third, the hybrid MC algorithm has to be improved to get better signals for electromagnetic and weak operators. Work in all these directions is under way. Finally, I mention a few physics topics that are presently also under investigation: • Reaction theory on the lattice. A general strategy for investigating capture reactions on the lattice as been outlined in Ref. [33]. The method consists of two major steps. First, projection MC is used to determine a multi-channel adiabatic lattice Hamiltonian for the participating nuclei. The second part is the calculation of the effective two-body capture reaction at finite volume using the adiabatic Hamiltonian. As a first application, the leading M1 contribution for p(n, γ)d is worked out in [33]. See also Ref. [34] for an ab initio approach for fusion reactions using chiral forces.
• Topological volume corrections. Scattering processes on a lattice can be analyzed using Lüscher's method [35]. However, if the bound states (like nuclei) move in a periodic volume, additional topological volume corrections arise that are sensitive to the number and mass of the constituents [36] (and its extension to QFT in [37]). The importance of these corrections has been explicitely demonstrated in Ref. [38] for elastic scattering between a fermion and a bound dimer in the shallow binding limit (and is presently investigated in the study of nd scattering [39]).
• Equation of state of neutron matter. The equation of state of neutron matter is currently a hot topic, see Steiner's talk at this conference. We are presently performing a NNLO analysis of the energy of neutron matter as a function of density, extending the work in Ref. [40], investigating also the pairing gap and possible P-wave pairing. For related works employing chiral EFT, see e.g.