# Unitary fermions and Lüscher’s formula on a crystal

## Abstract

We consider the low-energy particle-particle scattering properties in a periodic simple cubic crystal. In particular, we investigate the relation between the two-body scattering length and the energy shift experienced by the lowest-lying unbound state when this is placed in a periodic finite box. We introduce a continuum model for s-wave contact interactions that respects the symmetry of the Brillouin zone in its regularisation and renormalisation procedures, and corresponds to the naïve continuum limit of the Hubbard model. The energy shifts are found to be identical to those obtained in the usual spherically symmetric renormalisation scheme upon resolving an important subtlety regarding the cutoff procedure. We then particularize to the Hubbard model, and find that for large finite lattices the results are identical to those obtained in the continuum limit. The results reported here are valid in the weak, intermediate and unitary limits, and can be used for the extraction of scattering information ,via exact diagonalisation or Monte Carlo methods, of two-body systems in realistic periodic lattices.

###### keywords:

Scattering theory, Effective field theory, Lattice fermions, Finite-size effects^{1}

## 1 Introduction

The notion of low-energy universality in particle-particle collisions is a powerful concept that has been around, in one form or another, for a very long time. The general idea consists of replacing realistic, complicated interactions with much simpler ones and renormalising their bare coupling constants in favour of exact (either experimental or theoretical) scattering properties. The model interaction, irrelevant at the two-body level, can then be used to investigate the effect of interactions with higher particle numbers without the extra complications of the realistic interactions. For example, in 1957, Huang and Yang introduced their s-wave regularised pseudopotential [1], which aimed at reproducing the exact scattering length of a realistic two-body process by means of a very simple model interaction. The idea, however, was put forward much earlier by Fermi in the 1930s [2], who studied neutron-nucleon scattering, and was able to fit the scattering length in the first Born approximation.

The calculation of two-body scattering properties with rather realistic two-body interactions, if these are spin-independent, spherically symmetric, and single-channel, is not too complicated, especially with the computers we have access to nowadays. It is genuine multichannel collisions that can become quite challenging to approach directly from the corresponding Lippmann-Schwinger equation. In nuclear physics, there is a fundamental interest in obtaining nucleon-nucleon scattering properties from first principles using lattice QCD (see [3] and references therein). In condensed matter and cold atomic physics, particle-particle scattering theory in a crystal, where the incident waves are Bloch waves, and with a realistic (e.g. screened Coulomb) interaction, constitutes a formidable problem: all bands of the periodic potential are coupled by the interaction, Galilean relativity does not hold so that the collisional properties depend on the total center of mass momentum, and spherical symmetry is not present and cannot be exploited.

In situations like those described above, it is of great appeal to be able to extract scattering properties without having to directly use scattering theory. In massive quantum field theories, Lüscher showed how to extract partial-wave scattering amplitudes, i.e. phase shifts and scattering lengths, from the energy shifts due to the interactions when the system is placed in a finite volume [4], thereby generalising an older result from non-relativistic quantum mechanics [1]. However, neither Lüscher’s nor Huang-Yang’s results for the scattering lengths hold when these are unnaturally large, i.e. near the unitary limit. This problem was elegantly solved using pionless effective field theory [5] by Beane et al. in ref. [6]. In this way, an analysis of low-energy nucleon-nucleon collisions, directly from QCD, was possible soon after [7]. In the case of particle-particle scattering in a crystal, finite-size effects in the forms given in refs. [4] and [6], however, do not appear to be known. In this article we derive the relevant formulas from effective field theory on a tight-binding lattice, that is, the Hubbard model. We will focus on the three-dimensional (space) case, for the one-dimensional case is well understood [8] and Lüscher’s formula in one dimension [4] holds. We derive the pertinent expressions using Bethe-Goldstone equations, and we begin by proving these via the introduction of a more general class of methods consisting of kernel subtractions. We then introduce the Hubbard model and analyse, in first instance, its naïve continuum limit. We then return to the Hubbard model itself and show that, at low energy and large volume, the energy shifts are completely analogous to those obtained from the naïve continuum limit, which shows once more how low-energy physics, at least in the two-particle case, is universal.

## 2 Kernel subtractions and Bethe-Goldstone theory

For general interactions, one way of establishing the relation between the scattering states of a system (i.e. the positive energy states in the infinite size limit) and the energy shifts in a finite box is to use Bethe-Goldstone (BG) theory [9] in vacuum, that is, without the static Fermi sea background. In fact, even if it looks very different from it, Lüscher rederived, and used the exact same old Bethe-Goldstone’s equation in his seminal paper [4]. Here we will relate BG theory to a much more general method – kernel subtractions – for solving homogeneous integral equations ^{2}

(1) |

where is the kernel (operator), with some parameter, which in the relevant case here will correspond to an energy eigenvalue. The above equation (1) is homogeneous and therefore only has a solution for particular values of . Consider now a vector that satisfies the following integral equation

(2) |

where and are arbitrary vectors. Then, for the values of the energy for which , it holds that .

We now particularize to the physically relevant case, and derive the BG equation. Let satisfy the stationary Schrödinger equation

(3) |

We assume the Hamiltonian is of the usual form . Then, if is not in the spectrum of , the wave function satisfies the homogeneous Lippmann-Schwinger equation

(4) |

where . The kernel, Eq. (1), is therefore given by . We choose the following vectors and ,

(5) | ||||

(6) |

where is an eigenstate of , i.e. . Then, satisfies

(7) | ||||

(8) |

Eq. (7) together with the auxiliary condition (8) correspond exactly to BG theory in operator form. In order to relate it to the finite-size reaction matrix, we introduce the identity and define . After straightforward algebra, we find

(9) | ||||

(10) |

## 3 Hubbard model

We shall focus on one specific model, namely the Hubbard model on a simple cubic lattice. This is the minimal model for interacting electrons (or any other kind of spin- fermions) in a crystal, and consists of a nearest-neighbour hopping, or kinetic energy term, and an on-site interaction term. Its Hamiltonian in second quantised form is given by

(11) |

Above, is the hopping rate, denotes nearest neighbours, with , () are fermionic creation and annihilation operators, is the on-site interaction strength, is the local number operator and is the total number operator (we have set the lattice spacing ). In Hamitonian (11), the last term on the right-hand side amounts to an energy shift so that the single-particle ground state energy vanishes. In particular, the kinetic energy dispersion, of the Hubbard model (11) is given by

(12) |

where the quasi-momenta , and the Brillouin zone is given by . If the system is placed on a finite cube of side length (recall the lattice spacing ) with periodic boundary conditions, the quasi-momenta are restricted to take values (mod ), with .

## 4 Naïve continuum limit

Before studying the Hubbard model in detail, it is very instructive to consider its naïve continuum limit. The bare continuum Hamiltonian is given by

(13) |

Above, is the effective mass, is the bare interaction strength, and () is the fermionic annihilation (creation) operator in the continuum.

The above Hamiltonian, Eq. (13), is obtained as follows. Firstly, we restore the lattice spacing (which was set to in the previous section), in such a way that the action of the non-interacting Hamiltonian on a single-particle wave function, in the first quantisation, is given by

(14) |

where , with , and where we have defined the unit vectors such that (). The continuum limit of Eq. (14) is attained by setting . The continuum limit of the bare interaction strength is attained by identifying as .

The effective interaction in Eq. (13) is a zero-range ”s-wave” interaction. It is renormalisable, and accounts correctly for the s-wave scattering length with a vanishing effective range. This interaction can be renormalised in a number of ways, all equivalent to each other. For instance, minimal subtraction in the hard cutoff regularisation scheme [5], dimensional regularisation, Huang’s pseudopotential [1], or the more recent method of Tan’s distributions [10]; [11], all give the same results. In all the aforementioned methods, spherical symmetry of the contact interaction is exploited. The question arises of whether the cubic symmetry of the lattice Brillouin zone in the continuum limit changes the status of affairs. The answer to this question, which seems quite trivial, is not straightforward at all. On the one hand, the scattering properties at low energies remain unchanged. On the other hand, Lüscher’s formula heavily relies on spherical symmetry, a fact which is particularly evident in the pionless EFT considered by Beane et al. in ref. [6], and we will show that Lüscher’s formula does change significantly when the lattice symmetry is respected in the regularisation process, unless special care is taken in the renormalisation process.

### 4.1 Scattering theory

In this subsection, we shall obtain the reaction matrix (-matrix), by regularising the contact interaction in such a way that the symmetry of the Brillouin zone is manifestly preserved at all steps. The regularised contact interaction in the momentum representation for a spin-singlet fermion pair with ultraviolet (UV) cutoff is therefore given by

(15) |

where we have defined

(16) |

Above, is the Heaviside step function. The above interaction, Eq. (15), has to be contrasted with the widely used spherically symmetric regularisation [5]; [6]; [11], which reads .

After separation of center of mass and relative coordinates, and at vanishing center of mass momentum, the bare -matrix with the interaction in Eq. (15) is straightforward to calculate, and reads

(17) |

where

(18) |

with denoting Cauchy’s principal value. The UV structure of , Eq. (17), can be extracted from its value at zero energy. Numerically, we have obtained . The zero-energy on-shell -matrix can be renormalised in favour of the scattering length by means of minimal subtraction, i.e. by setting

(19) |

in which case the renormalised coupling constant takes the value . The -matrix at arbitrary energy is renormalisable if is finite for all energies , that is, if the linearly divergent part of the integrals is energy-independent. To show that this is indeed the case, write the integral as

(20) |

The second integral on the right hand side of Eq. (20) can be calculated as a series expansion that is convergent for any finite energy. We expand the integrand as

(21) |

The resulting integral in Eq. (20) is therefore given by

(22) |

where . The integral of the -th term in the series is bounded as

(23) |

Therefore, the following bound holds

(24) |

In the limit , the integral is bounded by zero from above and below, and it therefore converges to zero. Consequently, and the renormalised -matrix in this regularisation scheme, for finite energies, is constant and identical to that obtained by means of spherically symmetric regularisation.

There is a subtlety that needs to be considered in order to obtain meaningful results. The integer cutoff in lattice sums with lattice sites (which we consider odd for convenience) is , which implies that the lattice site number-dependent real cutoff is given by

(25) |

The interaction must be renormalised by using the limit of the cutoff, , since the limit is taken at the start of the calculation when doing scattering theory, that is, the renormalisation prescription reads

(26) |

where . Obviously, in the limit , we recover the result from minimal subtraction.

### 4.2 Lüscher’s formula

We place the system in a cubic box of side length with periodic boundary conditions, at vanishing pair momentum . We will consider here the lowest unbound state. From BG theory, Eq. (9) and Eq. (10), and setting , we obtain

(27) |

where the sum is assumed to be regularised and is renormalised to any given order by using the renormalisation prescription for the bare coupling constant , Eq. (26). To next-to-leading order in , the energy shift reads

(28) |

where we have defined

(29) |

In the above sums, is regularised with a cutoff (given by Eq. (25), and will be seen to be renormalised by subtracting ), calculated using the cutoff , while for the sums are regular and convergent. If we write

(30) |

we obtain, to lowest order

(31) |

where . Upon renormalisation, we obtain

(32) |

The above Equation (32) is in perfect agreement with Lüscher’s formula [4], and emphasizes the importance of the subtlety regarding the renormalisation of the effective interaction in the naive continuum limit, as discussed in the previous section.

There are several advantages of working with Eq. (32) with respect to the way one obtains Lüscher’s formula in the continuum limit of a cubic lattice [4]. Firstly, there is no need to project the interaction and the -matrix onto their -wave components, and since Eq. (32) only depends on the scattering length, this can be extracted with ease. Secondly, ”brute force” summations (i.e. direct computation of the sums, equivalent to exact diagonalisation) provide smooth convergence for all the sums in Eq. (29), and equally for the energy-dependent sums in Eq.(27). This is to be contrasted with the direct computation of the sum analogous to in the -wave channel [4]; [6]. In that case, the sums are restricted to , while in the present case these are restricted to (. In Fig. 1 we plot the value of the sums

(33) | ||||

(34) |

as functions of the cutoff , where we have subtracted , while the second sum is renormalised with the spherically-symmetric prescription of ref. [6]. Notice that the standard result for Lüscher’s sum is obtained here as . As is clearly observed in Fig. 1, direct computation of the sums using projections of the interactions onto the -wave channel gives rise to potentially high numerical errors. These would be inevitable in any computation relying on exact diagonalisation techniques, and would result in poor estimates of the scattering length.

### 4.3 Unitary limit

In the unitary limit, i.e. for , Eq. (32) is obviously incorrect. In the spherically symmetric case, Beane et al. [6] used pionless effective field theory to derive the finite-size correction to the ground state energy when the scattering length is unnaturally large. This corresponds to the limit . Here, we derive the corresponding expressions for the continuum limit of the lattice model in two different cases: (i) Beane et al.’s case () and (ii) the more natural case of in the limit .

For unnaturally large scattering lengths, Eq. (27) gives

(35) |

Clearly, we have in this case, and therefore is of the same order as the dispersion at low momenta. The energy shift in this case is given by

(36) |

Again, the leading order term is identical to the one obtained in the spherically symmetric case [6].

In the case of large scattering lengths that scale as the length of the system, , with , Eq. (27) gives

(37) |

Once more, , and by setting , Eq. (37) becomes

(38) |

The weak-coupling expansion in this case has the form

(39) |

Comparing Eqs. (36) and (39), it is easy to observe that it is not always possible to distinguish the cases (i) and (ii) referred to above, corresponding to unnaturally large and ”naturally” large scattering lengths, respectively, based on only one box size in the lowest energy shift. It can be tested, however, by choosing two different lattice sizes without changing the interaction parameters.

## 5 Lattice case

In the previous section we introduced the analog of Lüscher’s formula in the weak-coupling [4] and unitary [6] limits for the naïve continuum limit of the Hubbard model. In this section we study these expressions in the Hubbard model itself.

### 5.1 Scattering theory at low energy

The Hubbard model, Eq. (11), is a one-parameter theory, since the physics of the model only depends on the ratio . At zero energy, its -matrix gives the scattering length. It does, however, have non-zero shape parameters and, of course, higher (even) partial waves, too. On the other hand, since it only depends on all shape parameters depend exclusively on the scattering length ^{3}

The on-shell zero-energy -matrix is trivially shown to take the value

(40) |

where

(41) |

Above, the energy dispersion becomes that of the Hubbard model, Eq. (12), for the relative motion of two particles at zero total momentum, i.e.

(42) |

The zero-energy integral is known exactly and was computed by Watson [12]. It has the numerical value . Therefore, the scattering length is given by

(43) |

where is the effective mass of the particles, related to the tunnelling rate by ().

### 5.2 Energy shifts

The energy shift of the lowest unbound state is readily obtained from BG theory, Eq. (9) and Eq. (10), and reads in this case

(44) |

In terms of the scattering length, Eq. (43), Eq. (44) becomes

(45) |

In the weak-coupling limit, we therefore obtain

(46) |

The sum on the right hand side of Eq. (46) can be estimated as follows. Consider a cubic region , with being much smaller than the Brillouin zone’s extent (). The sum can be split into the two regions and . The sum over can be approximated by an integral, and therefore

(47) |

In the region , the full lattice dispersion can be approximated by

(48) |

setting , we obtain in the limit with

(49) |

where and , as calculated in the previous section. The integral in the region is easily calculated as

(50) |

Using the above results in Eq. (46), the energy shift becomes

(51) |

The above expression coincides with the corresponding formula in the naïve continuum limit, Eq. (32), as expected. However, Eq. (51) has been obtained without assuming the continuum limit at all energies, and it proves that the naïve continuum limit – with its associated subtleties – provides an adequate description of low-energy scattering in the weak-coupling limit.

The corresponding expressions in the unitary limits (both the natural and unnatural cases) are obtained in a fashion completely analogous to the one leading to Eq. (51) and are identical to the expressions obtained using the naïve continuum limit in the previous section.

## 6 Conclusions and outlook

Besides purely theoretical implications, the results reported in this work have interesting applications. For instance, two-particle scattering in a realistic three-dimensional lattice is of great relevance for condensed matter and cold atomic systems in optical lattices, as it is necessary in order to renormalise the coupling constant () in the Hubbard model. Unfortunately, the calculation of low-energy scattering properties of this system via direct use of collision theory is generally very involved, the reason being the genuine multi-channel structure of these systems. One way to estimate the scattering properties is to solve for the ground state energy of the realistic model consisting of only a few lattice wells (the number of wells being ). With the formulas derived in this article, one can fit the scattering length, and therefore the interaction coupling constant, to match the energy shift in the realistic calculation in a small box. Generalisations of the current work to study, e.g. effective interparticle interactions in graphene [13] and graphene-like optical lattices [14], are possible by using the methods presented here. In intrinsic graphene, as opposed to the cold atomic case and extrinsic graphene [15], however, the inclusion of the Coulomb interaction is necessary and the problem is more complicated, but can be done by extending the methods for continuous theories [16]. In particular, the rigorous construction of the continuum quasi-relativistic (it is not Lorentz invariant) effective field theory around the Dirac points is still an open problem ^{4}

It will be interesting to study the so-called Bertsch parameter [18], defined as the ratio between the ground-state energy of a resonantly interacting three-dimensional Fermi gas in the s-wave channel to the non-interacting ground state energy, but for the Hubbard model itself. This can be done using quantum Monte Carlo simulations [19]; [20]. Given that the energy shift, considered here, for two fermions at the lattice resonance is the same as that obtained from the pure s-wave results, it is to be expected that the lattice version of Bertsch parameter will be identical, or similar, with its continuum, pure s-wave counterpart. Its calculation would actually be much less involved since the elimination of higher-partial waves and effective range, cleverly done in ref. [19], is not necessary at all.

## Acknowledgements

We thank L. G. Phillips for illuminating discussions. The authors acknowledges support from EPSRC grant No. EP/J001392/1 and from the Danish Council for Independent Research under the Sapere Aude program.

### Footnotes

- journal: XX
- Here, integral equation is taken in the most general sense, and applies equally to discrete kernels.
- This is analogous to what happens with a hard-sphere interaction [1].
- The two-body problem in the naïve continuum limit has been studied in, e.g. ref. [17]. There, however, the interaction is purely phenomenological.

### References

- K. Huang and C. N. Yang, Phys. Rev. 105, 767 (1957).
- E. Fermi, Ricerca Sci. 7, 13 (1936).
- E. Epelbaum, H. -W. Hammer and U. -G. Meissner, Rev. Mod. Phys. 81, 1773 (2009).
- M. Lüscher, Commun. Math. Phys. 105, 153 (1986).
- S. Weinberg, Phys. Lett. B 251, 288 (1990); D. R. Phillips, S. R. Beane and T. D. Cohen, Nucl. Phys. A 631, 447 (1998); D. B. Kaplan, M. J. Savage and M. B. Wise, Phys. Lett. B 424, 390 (1998); Nucl. Phys. B 534, 329 (1998).
- S. R. Beane, P. F. Bedaque, A. Parreño and M. J. Savage, Phys. Lett. B 585, 106 (2004).
- S. R. Beane, P. F. Bedaque, K. O. Orginos and M. J. Savage, Phys. Rev. Lett. 97, 012001 (2006).
- M. Valiente, Phys. Rev. A 81, 042102 (2010); M. Valiente and D. Petrosyan, J. Phys. B: At. Mol. Opt. Phys. 42, 121001 (2009); ibid. 41, 161002 (2008); M. Valiente, M. Küster and A. Saenz, EPL 92, 10001 (2010); M. Valiente and K. Mølmer, Phys. Rev. A 84, 053628 (2011);
- A. L. Fetter and J. D. Walecka, Quantum Theory of Many-Particle Systems (Dover, N.Y., 2003).
- S. Tan, Ann. Phys. (N.Y.) 323, 2952 (2008).
- M. Valiente, Phys. Rev. A 85, 014701 (2012).
- G. N. Watson, Quart. J. Math., Oxford Ser. 2 10, 266 (1939).
- A. H. Castro Neto, F. Guinea, N. M. R. Peres, K. S. Novoselov and A. K. Geim, Rev. Mod. Phys. 81, 109 (2009).
- L. Tarruell, D. Greif, T. Uehlinger, G. Jotzu and T. Esslinger, Nature 483, 303 (2012).
- J. Hofmann, E. Barnes and S. Das Sarma, Phys. Rev. Lett. 113, 105502 (2014).
- S. R. Beane and M. J. Savage, Phys. Rev. D 90, 074511 (2014).
- C. Gaul, F. Domínguez-Adame, F. Sols and I. Zapata, Phys. Rev. B 89, 045420 (2014).
- G. A. Baker, Phys. Rev. C 60, 05311 (1999).
- M G. Endres, D. B. Kaplan, J.-W. Lee and A. Nicholson, Phys. Rev. A 87, 023615 (2013).
- J. E. Drut, Phys. Rev. A 86, 013604 (2012).