Phys. Rev. B 47, 11897 (1993)

Home | Retrieve | Search | Browse | Feedback | Help


Previous article | Next article | Up to abstract


---------------



This ASCII version was used to create the search index. No attempt has been made to render any mathematics or special characters correctly and errors in the conversion are quite common. If you need a more accurate rendering of the article, please use the page images or the PDF.

View Page ImagesPDF (693 kB), or ASCII text (41 kB)


Two-hole d-wave binding in the physical region of the t-J model: A Green's-function Monte Carlo study

Massimo Boninsegni and Efstratios Manousakis
Department of Physics, Center for Materials Research and Technology, Supercomputer Computations Research Institute, Florida State University, Tallahassee, Florida 32306
(Received 10 December 1992)

We investigate numerically the ground state of two holes with dx2-y2 symmetry in the two-dimensional t-J model, on lattices of significantly larger sizes than the ones studied so far. The Green's-function Monte Carlo method, in conjunction with suitable initial states and guidance functions, is used to estimate the ground-state energy of one and two holes, as well as to compute the two-hole distribution function. Our results show a significant decrease of the two-hole binding energy as the size of the lattice is increased from 4 x 4. A critical value Jc~0.27t is found such that d-wave hole binding no longer occurs for J<Jc.

©1993 The American Physical Society.

PACS: 74.20.Mn 02.70.Lq 75.10.Jm

I. INTRODUCTION

The problem of two mobile holes in the two dimensional (2D) t \(hy J model^{1}has elicited significant attention for its possible relevance to superconductivity in the cuprous oxides. The 2D t \(hy J model attempts to describe mobile holes in a 2D quantum antiferromagnet by adding a simple nearest neighbor hole hopping term to the spin cf12 antiferromagnetic Heisenberg model (AFHM):

H hat = - t ^ sum from { langle ij rangle s } ^ ( a hat_{is^{dag}a}hat_{js}+ H.c.) + J ^ sum from { langle ij rangle } ^ ( s vec_{i}cdot s vec_{j}- case 1 over 4 n hat_{i}n hat_{j}) ,   L (1)

where a hat_{is^{dag}=}c hat_{{i}s}^{dag}( 1 - n hat_{i-s})_, c hat_{is^{dag_}being}the creation operator for an electron with spin projection s at lattice site i, and where n hat_{i}= tsum_{s}c hat_{is^{dag}c}hat_{is_}is a number operator; a hat_{is^{dag_}creates}an electron only on an empty site, thus avoiding double occupancy. s vec_{i_}is the spin operator associated with site i and is defined as s vec_{i}= cf12 tsum_{{}alpha beta } c hat_{{}i alpha }^{dag}sigma vec_{{}alpha beta } c hat_{{}i beta }_, where sigma vec is a vector of Pauli matrices. We consider a square lattice of N = L times L sites with periodic boundary conditions. Several analytical approaches have been employed to study the Hamiltonian (1), leading to important theoretical predictions mainly regarding the physics of a single hole.^{2}Hole pairing has been investigated mostly numerically, both in the t \(hy J model^{3-4}and in the related strong coupling limit of the Hubbard model,^{5}by means of exact diagonalization on lattices of small size, such as 4 times 4. These studies have yielded some indications that the formation of a bound state of two holes may be energetically favorable within certain ranges of the parameters of the models. However, there are reasons to suspect that the relatively small size of the lattices studied may not afford a reliable prediction of the physics of the infinite system. Therefore the numerical investigation must be extended to lattices of significantly larger sizes than the ones accessible to exact diagonalization, which is limited by computer memory constraints. Monte Carlo techniques, with their relatively modest memory requirements, are a valid computational alternative; in particular, the Green's-function Monte Carlo (GFMC) method^{6}has been successfully used to compute ground state expectation values for the no hole case of the Hamiltonian (1), i.e., the spin cf12 AFHM.^{7}This method consists of projecting out the ground state of a given Hamiltonian from a starting trial state by acting iteratively with a suitable projection operator. Its application to the t \(hy J model in the case of mobile holes is complicated by the presence of the well known ``minus'' sign problem, which affects Monte Carlo simulations of systems with fermionic character.^{6}Because of increasingly larger statistical fluctuations, it is usually impossible to iterate the algorithm until convergence to the ground state is observed.

However, as we have recently shown for the case of a single hole,^{8}such complications can be overcome if a sufficiently accurate initial trial state and a suitable guidance function are used. In that case, convergence to the ground state can be actually observed before the statistical fluctuations become too large and accurate ground state estimates can be obtained by performing a transient estimation.

In this paper we extend our implementation of the GFMC method to carry out a numerical study of the ground state of two mobile holes in the Hamiltonian (1) on lattices of significantly larger sizes than the ones studied so far. We investigate the binding of two holes as a function of the parameter J/t of the model by computing the ground state two hole binding energy DELTA = delta E_{2_} - 2 delta E_{1_,}where delta E_{M}= E_{M}- E_{0_,}E_{M_}being the ground state energy of the Hamiltonian (1) with M holes. We gain further insight into the occurrence of binding by studying the ground state two hole distribution function.

Let us begin by reviewing what is known about the physics of the 2D t \(hy J Hamiltonian in the presence of two holes. In the large-J/t limit, two holes in the ground state will reside on nearest neighboring sites, as this configuration features the least number of broken antiferromagnetic bonds; this was observed in the static (t=0) limit by exact diagonalization on a 4 times 4 lattice.^{9}As a result of the same mechanism, if more than two holes are present, they will aggregate in a single cluster (phase separation). Monte Carlo simulations for two and four holes in the static limit^{10}and for relatively large values of J/t (J/t > 5) (Ref. 11) have confirmed these intuitive predictions. It should be noted that in this limit the ``sign'' problem is not too important (it actually disappears at t=0). In the opposite limit (J/t < 1), the higher hole mobility renders the situation much less clear; this limit is considered physically more realistic. It is argued by some authors^{12}that phase separation will take place at \fIany\fR value of J/t; on the other hand, numerical results obtained by means of high temperature expansion^{13}and exact diagonalization^{14}show evidence of hole binding and no phase separation at J/t <= 1.

The results that we present in this paper for the ground state of two holes in the t \(hy J model indicate that the magnitude of the two hole binding energy DELTA on an infinite lattice is considerably reduced from its value on a 4 times 4 lattice; for instance, at J/t = 0.4, which is a value in the physical regime of the t \(hy J model, we estimate DELTA to be app 0.1 t in magnitude on an infinite lattice, as opposed to the value of 0.349t found on a 4 times 4 lattice. We find that the finite size effects on DELTA are mainly due to the single hole energy delta E_{1_,}which shows a more pronounced dependence on the size of the lattice than the two hole energy delta E_{2_.}Results for the two hole distribution C(r) give further indications of a weakened binding of the holes at J/t < 1 as the lattice size is increased. Our results also indicate that a critical value J_{c_}of the parameter J exists such that hole binding only takes place for J >= J_{c_;}we find J_{c}app 0.27_ (with t=1).

In the next section, we will briefly sketch the implementation of the GFMC method used in this paper (for a more detailed description, see Ref. 8); in Sec. III we will introduce the two hole variational state which we use as the initial state in our GFMC simulation, whereas Sec. IV will be devoted to the discussion of our results.

II. THE METHOD

The GFMC method permits the calculation of the ground state expectation value scrO of an observable scrO hat as the limit of the sequence

scrO^{(n)}= { langle PSI_{T}| G hat^{n}scrO hat G hat^{n}| PSI_{T}rangle } over { langle PSI_{T}| G hat^{2n}| PSI_{T}rangle } ,   (2)

as n -> inf. Here | PSI_{T}rangle_ is an initial state not orthogonal to the ground state | PSI_{0}rangle_ and G hat is a suitable projection operator, whose role is enhancing the component of | PSI_{T}rangle_ along | PSI_{0}rangle_. If the spectrum of the Hamiltonian H hat of the system is bounded, as is the case of (1) on a finite lattice, a convenient choice for G hat is G hat = E - H hat with E >= E_{max_,}the largest eigenvalue of H hat .^{8}For the case of two holes in the t \(hy J model, it is E_{max}= 8t_, which corresponds to the maximum kinetic energy of two free holes in a ferromagnetic background.

A Monte Carlo implementation of this scheme consists of evaluating (2) as an average over the random walks of scrN independent ``walkers'' through configuration space. The walkers are initially assigned scrN starting points, corresponding to the scrN configurations | c_{1}rangle , c_{2}rangle , .^.^.^ , | c_{scrN}rangle_ of the system. These configurations are stochastically drawn from the probability distribution | PSI_{T}(c) | PSI_{G}(c)_, where PSI_{T}(c)_ is the wave function of the trial state | PSI_{T}rangle_ and PSI_{G}(c)_ is a positive definite ``guidance'' function, based upon some physical insight on the system and whose role is to guide the random walk toward the most important configurations.

The random walk of each walker consists of a succession of transitions to new configurations; a transition from a given configuration | x rangle to a new configuration | x prime rangle is stochastically sampled from a probability distribution w ( x -> x prime ) proportional to | langle x prime | G hat | x rangle | PSI_{G}( x prime ) /_\p PSI_{G}(x)_. In order to obtain a statistical estimate of the quantity scrO^{(n),}each walker must be allowed to perform 2n transitions, since every transition corresponds to one operation of G hat.

The quantity scrO^{(n)}must be computed for increasingly larger n until convergence is observed, within statistical error bars. In this paper we calculate the ground state energy of a quantum antiferromagnet with mobile holes by following the same procedure explained in Ref. 8 for the single hole calculation; we use a correlated two hole variational state as the initial state | PSI_{T}rangle_. We also compute the ground state two hole distribution function by means of a scheme known as ``forward walking'' (see Ref. 6). The difficulty, in the case of the t \(hy J model and other fermionic problems, arises from the matrix elements langle x prime | G hat | x rangle not being always positive. This results in what is known as the ``minus'' sign problem, i.e., in large statistical fluctuations of the quantity scrO^{(n)}as n increases, which clearly make it very difficult to perform a large number of iterations because of the increasing statistical uncertainty affecting the estimates. Therefore it is very important for the success of the method to have a sufficiently accurate initial state so that convergence to the ground state can be achieved within a relatively small number n of iterations, before the statistical fluctuations grow exceedingly. The algorithm can be also significantly enhanced by the use of a suitable guidance function, which by directing the random walk toward the most important configurations can prevent the statistical variance of the estimates from growing too rapidly, thereby improving the change of observing convergence.

III. TRIAL STATE FOR TWO HOLES

Our GFMC calculation for two holes is based upon an initial state which is a generalization of the ``string''-based variational state used in Ref. 8 in the GFMC calculation for a single hole. Let us therefore recall the single hole variational ansatz first (see also Ref. 17):

| PSI_{T}( k vec ) rangle = ^ sum from { R vec c } ^ e^{{}- i k vec cdot R vec } F hat ( k vec ) exp ( - 1 over 2 ^ sum from { i < j } ^ u_{ij}s hat_{i^{z}s}hat_{j^{z}right}) a hat_{{}R vec s } | c rangle ,  L

where k vec is the hole momentum and where the sum runs over lattice sites R vec and over all lattice spin configurations lcurl | c rangle rcurl, specified by assigning the value of the spin projection s^{z}at all lattice sites; the sum in (3) is restricted to configurations | c rangle with a value of the z component of the total spin S^{z}= 0. The operator exp ( - cf12 tsum_{{}i < j } u_{ij}s hat_{i^{z}s}hat_{j^{z})_}is a background antiferromagnetic spin spin correlation operator, and the function u_{ij_}depends on the distance between the two sites i and j. Finally, F hat ( k vec ) = 1 + ^ sum from a vec ^ f_{a}vec ( k vec ) scrP hat_{a}vec + ^ sum from { a vec a vec prime } ^ f_{{a}vec a vec prime} ( k vec ) scrP hat_{{}a vec prime } scrP_{a}vec ,   (3)

with a vec = +- x bhat , +- y bhat connecting two nearest neighboring sites and scrP hat_{a}vec = tsum_{{}R vec s } a hat_{{}R vec s }^{dag}a hat_{{}R vec + a vec s }_. F hat ( k vec ) is a spin hole correlation operator; f_{a}vec ( k vec )_ and f_{{}a vec a vec prime } ( k vec )_ are variational parameters. This ansatz describe ``strings'' of spins displaced by one site along the path. In order to generalize the state (3) to the case of two holes, consider first a general, translationally invariant state of two holes with opposite spin projections in an antiferromagnetic spin background:

| PSI_{0}( Q vec ) rangle = ^ sum from c ^ sum from { R vec , r vec } ^ mark e^{{}- i Q vec cdot ( R vec + r vec / 2 ) } exp ( - cf12 ^ sum from { i < j } ^ u_{ij}s hat_{i^{z}s}hat_{j^{z}right})  

lineup times g ( r vec ) a_{{}R vec uarrow } a_{{}R vec + r vec darrow } | c rangle ,   (4)

where Q vec is the total momentum of the state and r vec is the relative distance of the two holes. Analogously to the single hole case, we introduce the ``string'' correlation operator

F hat ( Q vec , r vec ) = 1 + ^ sum from a vec ^ f_{a}vec ( Q vec , r vec ) scrP hat_{a}vec + ^ sum from { a vec a vec prime } ^ f_{{}a vec a vec prime } ( Q vec , r vec ) scrP hat_{{}a vec prime } scrP hat_{{}a vec } ,  L (5)

where f_{a}vec ( Q vec , r vec )_ and f_{{}a vec a vec prime } ( Q vec , r vec )_ are variational parameters. By acting on the state (4) with F hat ( Q vec , r vec ), we obtain our variational state for two holes:

| PSI_{T}( Q vec ) rangle = ^ sum from c ^ sum from { R vec , r vec } ^ mark e^{{}- i Q vec cdot ( R vec + r vec / 2 ) } F hat ( Q vec , r vec )  

lineup times exp ( - 1 over 2 sum from { i < j } ^ u_{ij}s hat_{i^{z}s}hat_{j^{z}right})  

lineup times g ( r vec ) a_{{}R vec uarrow } a_{{}R vec + r vec darrow } | c rangle .   (6)

The expectation value of the energy is minimized by taking g ( r vec ) to be nonzero for nearest neighbor distances only and with d-wave spatial symmetry, i.e., g ( +- x bhat )\p = - g ( +- y bhat ), corresponding to a singlet state of the two holes. With this choice of g, the variational parameters f_{a}vec ( Q vec , r vec )_ and f_{{}a vec a vec prime } ( Q vec , r vec )_ can be computed approximately analytically by minimizing the energy expectation value with u=0. In our calculation we allowed for strings of length 1 only; i.e., we set f_{{}a vec a vec prime } ( Q vec , r vec ) = 0_. We found that allowing for strings of length 2 can improve the initial variational energy estimate, but has little effect on the convergence to the ground state; namely, the ground state estimates are the same, within statistical error bars, with the ones obtained by including strings of length 1 only. We restricted our calculation to a state with Q vec = (0,0); in this case, it is f_{a}vec ( Q vec , r vec ) = f_{0_,}with f_{0_}real. This choice of variational state is consistent with exact diagonalization results according to which the lowest energy state of two holes is a d-wave singlet with momentum Q vec = (0,0), although degeneracies exist which are attributed to peculiar geometrical properties of the lattice studied.^{3}

As we did in the single hole calculation, we used the function

PSI_{G}= exp ( - 1 over 2 ^ sum from { i < j } ^ u_{ij}s_{i^{z}s}sub j^{z}right )   (7)

as a guidance function. We compute the two hole energy delta E_{2}= E_{2}= E_{0_}by performing a transient estimation of the ground state energy of (1) with two holes, E_{2_,}from which we subtract the value of the ground state energy E_{0_}of the no hole case (corresponding to the spin cf12 AFHM), also calculated by the GFMC method; note that in the no hole case no ``minus'' sign problem arises, so that E_{0_}can be estimated rather straightforwardly and with remarkable accuracy (see Ref. 7).

IV. RESULTS AND DISCUSSION

Let E_{M^{(n)_}be}the estimate of the ground state energy of (1) with M holes obtained at the nth GFMC iteration, and let E_{M_}be the extrapolated value in the limit n -> inf. In Fig. 1 the estimates delta E_{2^{(n)}=}E_{2^{(n)}-}E_{0_}at different n, on a 4 times 4 lattice (open circles), are compared to the exact diagonalization results (dashed line) for J/t = 1. At this relatively large value of J/t, it is possible to obtain sufficient evidence of convergence; upon averaging the data from the last three iterations shown, we estimate delta E_{2_}on a 4 times 4 lattice to be delta E_{2}/ t = 0 +- 0.01_, in agreement with the exact value of 0.42, which corresponds to the dashed line in Fig. 1. In Fig. 2 we show analogous results for delta E_{2_}on an 8 times 8 lattice; convergence to the ground state can be observed, and if we average the values of the last two iterations we obtain delta E_{2}/ t = 0 +- 0.02_. The dotted lines in Figs. 1 and 2 correspond to a fit of the data with the expression

E_{2^{(n)}approx}E_{2^{inf}+}beta ^ exp ( - kappa n ) ,   (8)

using E_{2^{inf_,}beta,}and kappa as fitting parameters. We find that E_{2^{inf_,}i.e.,}the extrapolated estimate of E_{2_,}is the best determined of the three parameters. The extrapolated values for delta E_{2_}for the two lattices are the same as the ones given above, within statistical error bars.

As we shall see below, our results indicate that the two hole energy estimates at J/t = 0.4 do not change significantly beyond 8 times 8; we expect this to be true at J/t = 1 as well, since the holes are less mobile in this case and therefore finite size effects will be even smaller than at J/t = 0.4. Thus we did not devote our computational resources to study lattices larger than 8 times 8 at J/t = 1.

Because we are ultimately interested in a bound state of two holes, we consider the binding energy DELTA = delta E_{2}- 2 delta E_{1_,}where delta E_{1}= E_{1}- E_{0_}is the single hole ground state energy. We compute delta E_{1_}in the same way as delta E_{2_;}we use the single hole trial state (3) with momentum k vec = ( pi / 2 , pi / 2 ) (where the single hole energy band attains its minimum) as the initial state; we allow for strings of length 1 and find values of delta E_{1_}in agreement with the ones reported in Ref. 8. By using the estimates obtained above for delta E_{2_,}we obtain DELTA / t = - 0 +- 0.03 for the 4 times 4 lattice, in agreement with the exact value, and DELTA / t = - 0 +- 0.03 for the 8 times 8 lattice. Since we know that the single hole energy estimates remain unchanged^{8}beyond 8 times 8 in the range 0.2 <= J/t <= 5, we conclude that the value of DELTA / t found on an 8 times 8 lattice should be sufficiently close to the value in the infinite lattice limit and that binding occurs in this limit at J/t = 1, with a binding energy of about - 0.6 t. A further indication of binding as well as of the limited finite size effects affecting the calculation at this value of J/t comes from the very similar values found for the rms hole separation on the two lattices, namely, R_{roman}rms = 1 +- 0.02_ (in agreement with the exact result^{4}) on a 4 times 4 lattice and R_{roman}rms = 1 +- 0.03_ on an 8 times 8 lattice, the lattice constant being set equal to 1.

We now turn to a more interesting and physically more relevant case, i.e., J / t = 0.4. In Fig. 3 we show transient estimation results for delta E_{2_}on a 4 times 4 lattice. The tendency of the succession of energy estimates to approach the exact value is clear, and if we average the values from the last two iterations, we obtain delta E_{2}/ t = - 2 +- 0.03_, in good agreement with the exact value of - 2.993. However, the rapid increasing of the error bar renders it much more difficult than in the single hole case to determine the convergence value precisely. The problem is present on larger lattices as well, as shown in Fig. 4 for delta E_{2_}(open circles) on an 8 times 8 lattice. The worsening of the sign problem passing from one to two holes is expected, particularly at J/t < 1, as a larger hole concentration results in an increased rate of sign changing transitions in the random walk; such transitions are responsible for the statistical fluctuations. Since the amount of information produced decreases dramatically with the number of iterations n, it is important to extract the maximum amount of information from the data generated at ``early time,'' that is, at small n, when the size of the error bars is not too large.

A simple scheme^{15}which effectively permits the acceleration of the convergence to the ground state with no additional computational cost consists of considering the expectation values

scrE^{(n)}= { langle chi_{n}| H hat | chi_{n}rangle } over { langle chi | chi_{n}rangle } ,   (9)

where

| chi_{n}rangle = ( 1 + lambda_{n}G hat^{n}) | PSI_{T}rangle ,   (10)

lambda_{n_}being a variational parameter determined by minimizing the value of scrE^{(n).}scrE^{(n)}yields a better ground state estimate than

E_{2^{(n)}=}langle PSI_{T}| G hat^{n}H hat G hat^{n}| PSI_{T}rangle / langle PSI_{T}| G hat^{2n}| PSI_{T}rangle ,  

owing to the greater variational freedom of the state | chi_{n}rangle_ compared to G hat^{n}| PSI_{T}rangle_. This Lanczos type procedure can be applied to expectation values other than the energy, although the parameter lambda_{n_}has to be determined by minimizing the energy expectation value.

The improvement is particularly important at small n, when the difference between successive values of E_{2^{(n)_}is}larger; this is shown in Figs. 3 and 4, where we compare the values of delta scrE^{(n)}= scrE^{(n)}-E_{0_}and delta E_{2^{(n)}=}E_{2^{(n)}-}E_{0_}(circles). We can obtain extrapolated estimates for delta E_{2_}by fitting the values for E_{2^{(n)_}and}scrE^{(n)}with the expression (8). The fits of the two sets of data yield the same extrapolated estimates, within statistical error bars. By using the single hole energy values, we estimate the magnitude of the binding energy of the two hole bound state in the infinite lattice to be greatly reduced from its value on a 4 times 4 lattice.

To illustrate this point more clearly, we plot in Fig. 5 the estimates of DELTA / t obtained as the difference delta E_{2^{(n)}-}2 delta E_{1^{(n)_,}at}J / t = 0.4 on an 8 times 8 lattice; the dotted line is obtained from the fitting curves for both E_{2^{(n)_}and}E_{1^{(n)_.}The}solid line represents the exact value of the binding energy for a 4 times 4 lattice. Despite the size of the error bars, the results shown in Fig. 5 indicate rather clearly that the magnitude of the binding energy is considerably reduced, namely, from a value of 0.349t on a 4 times 4 lattice to about 0.1t on an 8 times 8 lattice.

Our results for the two hole energy delta E_{2_}at J/t = 0.4 indicate that delta E_{2_}does not change significantly in the infinite lattice from its value on the 8 times 8 lattice; thus, because the single hole energy also remains unchanged beyond 8 times 8 ,^{8}we estimate the binding energy DELTA to be (- 0 +- 0 ) t in the infinite lattice. The extrapolated estimates for delta E_{2_,}delta E_{1_,}and DELTA for different values of J/t and on the various lattices studied are reported in Table I.

The main contribution to the significant size dependence found for DELTA arises from the single hole ground state energy values. A possible explanation for the different size dependence of the results for one and two holes may come from a long range effect of a single hole on the antiferromagnetic background; the addition of a second hole may result in an opposite contribution and therefore in a substantial cancellation of such an effect. As an example, consider the possible long range planar distortion of the antiferromagnetic moment of the spin background caused by the motion of a single hole.^{16}Far away from the hole, the distortion delta m^{dag}is proportional to k vec cdot r vec / r, where k vec is the hole momentum; in a two hole state with total momentum equal to zero, the distortions caused by the two holes have opposite signs and the net effect is a dipolarlike field decayings as app 1 / r^{3}at large distances, as opposed to app 1 / r for a single hole.

Further insight into the two hole ground state properties is provided by the calculation of the two hole distribution function C(r), defined as

C (r) = 1 over N ^ sum from { i < j } ^ langle h hat_{i}h hat_{j}rangle delta ( | R vec_{j}- R vec_{i}| - r ) ,   (11)

where h hat_{i}= 1 - n hat_{i_}is the hole number operator associated with the ith site, which is positioned at R vec_{i_.}In Figs. 6 and 7, we show the results for C (r) on a 4 times 4 and on an 8 times 8 lattice at J / t = 1; different symbols refer to the values of C(r) computed after different numbers of GFMC iterations. The dotted line refers to the data from the zeroth iteration, whereas the dashed line refers to the data from the last iteration performed. As we see, the holes are initially located at distances 1, sqrt 2, and 2 from one another; that is, they are either located on nearest neighboring or next nearest neighboring sites. As the algorithm is iterated, the weight is partially shifted to larger distances, but 1 and sqrt 2 remain largely dominant. Since the values of C(r) change very little in the last three iterations sown for both lattices, we conclude that the values from the last iteration give a correct representation of the ground state two hole distribution. We compute the rms hole separation and obtain R_{roman}rms = 1 +- 0.2_ at the last iteration for the 4 times 4 lattice; this value is in agreement with the exact one.^{4}On the 8 times 8 lattice, the rms separation is equal to R_{roman}rms = 1 +- 0.03_.

Figures 8 and 9 show the same calculation for the two lattices at J/t = 0.4. Here, as the algorithm is iterated, the weight is significantly redistributed to the larger distances, and little evidence of convergence can be found. The value of R_{roman}rms_ obtained from the last iteration on the 4 times 4 lattice is R_{roman}rms = 1 +- 0.03_, lower than the exact one. We attribute this difference to the need of iterating the GFMC algorithm longer in order to recover fully the large distance features of the function C(r). The fact that the same number of iteration yields better evidence of convergence for the energy than for the two hole distribution can be explained as being from the energy not being very sensitive to the form of C(r) at large distances. From the results for C(r) on an 8 times 8 lattice, we can see that there is an evident tendency of the two holes to ``spread'' at distances larger than the ones at which they are positioned in the initial state and larger than the ones available on a 4 times 4 lattice. This is consistent with the decrease in the binding energy found as the lattice size is increased and is a further indication of important finite size effects affecting the calculation on a 4 times 4 lattice. We estimated the rms hole separation R_{roman}rms_ by fitting the values of R_{roman}rms^{(n)_}obtained at different GFMC iterations with an exponential. On a 4 times 4 lattice, we find R_{roman}rms = 1 +- 0.04_, in agreement with the exact value within error bars. On lattices of larger size, the asymptotic value of R_{roman}rms_ given by the fit is significantly larger than the one from our last iteration, which indicates the need for a larger number of iterations. For instance, on a 10 times 10 lattice, R_{roman}rms^{inf}= 2 +- 0.04_, whereas the estimate from the last (48th) iteration is only equal to 1 +- 0.03. From the fit we estimate the number of iterations needed to reach convergence to be of the order of 200 on the 10 times 10 lattice, whereas we could only perform 60 iterations because of the rapid increase of the size of the error bars. The extrapolated value of R_{roman}rms_ is 2 +- 0.04 on both the 8 times 8 and 10 times 10 lattices.

The rapid decrease of DELTA with increasing t/J for an infinite lattice is in marked contrast with the results on a 4 times 4 lattice, where DELTA seems to vary mildly with t/J (in fact, DELTA app - 0.8J for 0 <= t/J <= 5 .^{4}Our results suggest that hole binding will no longer occur for J less than a critical value J_{c}< 0.4 t_; such a critical value is known to exist in two dimensions for a d-wave bound state in the continuum case.^{18}We can estimate the value of J_{c_}by fitting the data in Table I by means of an extrapolation formula; it should be noted that in region J < 0.4t the binding energy is very small and our GFMC calculation cannot resolve the energy difference delta E_{2}- 2 delta E_{1_}with sufficient accuracy. We can obtain an extrapolation formula by modeling the complicated many body problem of two mobile holes in a quantum antiferromagnet with a simple effective Hamiltonian describing the motion of two interacting quasiparticles in a square lattice with N sites and periodic boundary conditions:

H hat_{eff}= mark ^ sum from { i,j < i , sigma } ^ t prime ( R vec_{i}, R vec_{j}) ( b hat_{{}i sigma }^{dag}b hat_{{}j sigma } + H.c.)  

lineup + sum from { i,j < i } ^ J prime ( R vec_{i}, R vec_{j}) n hat_{i}n hat_{j} ,   (12)

where b hat^{dag}is a ``dressed'' quasiparticle creation operator and n hat is a number operator; the sum runs over all pairs of lattice sites, as well as over all spin projections sigma; R vec_{i_}is the position vector of the ith site. The quasiparticles move in a band E ( k vec ), which is assumed to have a minimum at k vec = ( +- pi / 2 , +- pi / 2 ). We assume that the effective interaction J prime ( R vec_{i}, R vec_{j})_ is attractive and short ranged, i.e., J prime ( R vec_{i}, R vec_{j}) = - J prime_ for nearest neighbors and zero otherwise; the effective Hamiltonian (12) neglects effects of quantum spin fluctuations. The eigenvalue problem can be solved exactly (see also Ref. 19); if we seek a solution corresponding to a state with zero total momentum and with the spatial symmetry of a d wave, we find the following condition for the existence of a solution:

1 = - 2 J prime ( I_{xx}- I_{xy}) ,   (13)

which determines the ground state energy E. We have introduced the quantities

I_{xx}= 1 over N ^ sum from k vec ^ { cos^{2}( k_{x}) } over { E - 2 E ( k vec ) } ,  

I_{xy}= 1 over N ^ sum from k ^ { cos ( k_{x}) cos ( k_{y}) } over { E - 2E ( k vec ) } ,  N (14)

and we have used E ( k vec ) = E ( - k vec ). In the vicinity of k vec_{min}= ( pi / 2 , pi / 2 )_, where E ( k vec ) attains its minimum, we can set E ( k vec ) approx E_{min}+ t_{1^{sprime}q}sub 1^{2}+ t_{2^{sprime}q}sub 2^{2_,}where q vec = k vec - k vec_{min_}and q_{1_}and q_{2_}are the components of q vec along the high symmetry direction (0,0) -> ( pi , pi ) and (0 , pi ) -> ( pi , 0 ), respectively. We introduce the binding energy DELTA prime = E - 2 E_{min_;}in the limit where DELTA prime is small, it is straightforward to check that, upon neglecting terms of order DELTA^{{}prime 2 } and higher, (13) can be expressed as

1 over {J prime } = 1 over { J_{c^{sprime}}}[ 1 - gamma prime DELTA prime ^ ln ( DELTA prime / curlep prime ) ] ,   (15)

where we have introduced the constants J_{c^{sprime_,}gamma}prime, and curlep prime. We have expressed all energy scales in units of t_{1^{sprime_,}since}t_{1^{sprime}>>}t_{2^{sprime_}for}a single hole in the 2D t \(hy J model,^{2}and neglected terms of order DELTA^{{}prime 2 } and higher (a linear term in DELTA prime has been conveniently absorbed in the cutoff energy scale curlep prime in the logarithm). If we assume that the function t/J ieq f ( DELTA / t ) for the problem of two holes in the 2D t \(hy J model will have an expansion analogous to (17) around DELTA / t = 0 (i.e., around J = J_{c})_, we write

1 over J = 1 over { J_{c}} ( 1 - gamma DELTA over t ln ( DELTA / curlep ) right ) ,   (16)

where we have omitted terms of order ( DELTA / t )^{2}and ( DELTA / t ) ln ( DELTA / curlep )^{2.}A simple fit of our extrapolated values of DELTA (Table I) yields J_{c}app 0.27t_.

The main conclusion that can be drawn from the results of the GFMC calculation for two holes outlined in this paper is that since DELTA (J/t < 0.4 ) < 0.1t, the actual value of J/t in the copper oxide superconductors (taking t app 0.5 eV) cannot be much different than J/t = 0.4 if this model is relevant for their superconductivity. There is an open question, namely, what happens if a finite fraction of holes is introduced on the infinite square lattice? The existence of a two hole bound state is known to be related to the occurrence of a many body pairing instability, to the extent that the system can be approximated by a continuum system where holes interact via a two body potential.^{18}Another possibility is that phase separation may take place for some range^{13,14}or at any value of J/t .^{12}

ACKNOWLEDGMENTS

This work was supported by the Office of Naval Research under Contract No. N00014 93-1-0189 and by the Supercomputer Computations Research Institute, which is partially funded by the U.S. Department of Energy under Contract No. DE FC05 85ER 250000.

FIGURE AND TABLE CAPTIONS




View Page ImagesPDF (693 kB), or ASCII text (41 kB)

---------------


Previous article | Next article | Up to abstract

Home | Retrieve | Search | Browse | Feedback | Help

E-mail: prola@aps.org