Extending J Walking to Quantum Systems: Applications to Atomic Clusters

The J-walking (or jump-walking) method is extended to quantum systems by incorporating it into the Fourier path integral Monte Carlo methodology. J walking can greatly reduce systematic errors due to quasiergodicity, or the incomplete sampling of configuration space in Monte Carlo simulations. As in the classical case, quantum J walking uses a jumping scheme ,to overcome configurational barriers. It couples the usual Metropolis sampling to a distribution generated at a higher temperature where the sampling is sufficiently ergodic. The J-walker distributions used in quantum J walking can be either quantum or classical, with classical distributions having the advantage of lower storage requirements, but the disadvantage of being slightly more computationally intensive and having a more limited useful temperature range. The basic techniques are illustrated first on a simple one-dimensional double well potential based on a quartic polynomial. The suitability of J walking for typical multidimensional quantum Monte Carlo systems is then shown by applying the method to a multiparticle cluster system consisting of rare gas atoms bound by pairwise Lennard-Jones potentials. Different degrees of quantum behavior are considered by examining both argon and neon clusters. Remarkable improvements in the convergence rate for the cluster energy and heat capacity, analogous to those found in classical systems, are found for temperatures near the cluster transition regions.


INTRODUCTION
We recently presented a method called J walking (for jump walking)' that was shown to greatly reduce systematic errors occurring in random walks in classical systems that arise because of quasiergodicity,' the incomplete sampling of configuration space. We demonstrated the technique's utility on small atomic clusters bound by Lennard-Jones potentials. Recently, Tsai and Jordan3 applied J walking together with histogram methods4 to a Monte Carlo simulation of the water octomer, and Strozak, Lopez, and Freeman applied J walking to a simulation of argon clusters absorbed on graphite.5 In this work, we extend the method to quantum systems, incorporating it into the Fourier path integral (FPI) Monte Carlo methodol-~gy.~ Path integral methods based on Feynman's formulation of quantum statistical mechanics' provide a prescription for the computer simulation of quantum many-body systems. In the FPI scheme, the quantum paths associated with each coordinate are represented by Fourier series expansions, leading to a generalized configuration space that has been extended to higher dimensions by the inclusion of the auxiliary degrees of freedom from the Fourier coefficients. This enhanced multidimensional space is amenable to the usual methods associated with the Metropolis, Rosenbluth, Rosenbluth, Teller, and Teller algorithm.* As in the classical case, quasiergodicity arises in quantum systems where the sample space contains two or more regions having a very low transition probability between them, resulting in bottlenecks that effectively confine the sampling to only one of the regions. ' The dichotomy of time scales characterizing the walks produces rapid motion within a confined region and very slow movement between the regions." This results in systematic errors that arise because of the finite length of the walks. They diminish with increased walk length, disappearing in the limit of infinitely long walks. A prototypical example is the double well potential where the wells are separated by a large barrier. For sufficiently low temperatures, the random walker is unable to cross over the barrier within the duration of the walk and hence never samples from the other well. J walking addresses this problem by coupling the usual small scale Metropolis moves with occasional large scale jumps that move the Metropolis walker to different regions of configuration space, in essence "jumping" over the barriers rather than walking over them. This is accomplished in classical systems by using the Boltzmann distribution obtained at a temperature sufficiently high for the sampling to be ergodic. The distribution generated by the high temperature walker (J walker) becomes the sampling distribution for the low temperature walker's attempted jumps. Because the distribution's peaks correspond to the potential minima, the J-walker's motion remains biased about the minima, greatly increasing the likelihood an attempted jump would be accepted. An analogous scheme can be developed for quantum systems simply by sampling from the distribution generated by a high temperature quantum walker. However, because of their larger number of param-eters, FPI methods offer a greater degree of flexibility than available in the Metropolis simulation of classical systems, allowing several possible variations on the theme. Since the configuration space is often dominated by the coordinate subspace, an alternative scheme is for the low temperature quantum walker to sample from a classical J-walker distribution. Another possibility is to use a quantum J walker with the same temperature, but having a lower mass, since quasiergodicity becomes less important in highly quantum systems where larger zero-point energies and tunneling effectively lower the configurational barriers.
cle cluster systems consisting of either Ne or Ar atoms bound by pairwise Lennard-Jones potentials and compare it to the usual FPI methods. Finally, in Sec. V we summarize our findings.

A. Fourier path integral Monte Carlo methods
Path integral formulations of quantum statistical mechanics have become increasingly useful in recent years as continuing developments in computer technology have provided the computational power necessary for their application to realistic systems. There have been two major numerical algorithms developed to evaluate path integrals, discretized path integral methods," and Fourier path integral methods,6 although at a fundamental level, the two prescriptions are essentially the same. " We describe here the application of J walking to the Fourier methods only and note in passing that, in principle, J walking can also be extended to discretized methods.

I. Action integrals
Two complementary implementations for generating the classical J-walker Boltzmann distributions were originally presented. The first ran the J walker in tandem with the low temperature walker. The low temperature walker occasionally attempted jumps to the current J-walker position by simply using the J-walker coordinates for its trial move. The second scheme ran the J walker beforehand and periodically stored the configurations in an external array. Subsequent jump attempts could be made by accessing the stored configurations via randomly generated indices. The tandem walker scheme required that the J walker be moved an extra number of steps whenever a jump was attempted in order to reduce correlations between the two walkers that caused systematic errors. The extra number of steps needed depended on the temperature difference between the high and low temperature walkers, increasing the computation time greatly as the difference became larger. This scheme is therefore more suited for parallel computation. The external distribution scheme, on the other hand, had only a modest overhead in computation (mostly the time required to generate the distribution), but demanded large storage facilities to handle the distribution arrays. Because of the additional computational overhead inherent in the FPI method (typically up to an order of magnitude), we did not consider the tandem walker method to be feasible for the simulation of quantum syswhere tems on scalar computers and so we limited our investigation to externally stored J-walker distributions.
Recent reviews of FPI methods give a comprehensive description of their development and application.6 We review briefly here the major concepts that are instrumental to the application of J walking. For notational simplicity, we limit the discussion to one-dimensional systems; the extension to multidimensional systems is straightforward and is described in the reviews.
Path integral methods enable the calculation of the quantum density matrix, in coordinate representation We begin in Sec. II with a brief review of the FPI method and the essential concepts that pertain to J walking. We then develop two quantum J walking methods based on the use of quantum and classical distributions, and contrast them with the more familiar classical Jwalking method. As in our original study of classical systems, we have tested the quantum J-walker method on a number of simple systems. For heuristic purposes, we begin in Sec. III with the simple model of the prototypical one-dimensional double well potential defined by a quartic polynomial. This potential consists of one well of fixed depth containing the global minimum and a second well of variable depth separated by a fixed barrier. The variable depth of the second well allows us to examine in a systematic manner quasiergodic behavior as a function of well depth as well as a function of temperature. A solution of high accuracy can be obtained easily using standard harmonic oscillator basis set expansion techniques so that errors arising from quasiergodicity can be quantified. In Sec. IV, we apply the J-walker method to quantum multiparti- where (6) wk,,,(w@> =pf#exp -1:: @p262k) 1 is the width of the Gaussian-like Fourier coefficient distribution and pfP is the free particle density The Gaussian width ak provides a natural length scale, growing as the temperature is decreased; like the thermal deBroglie wavelength, it provides a measure of the extent of quantum contributions in the system. The quantum density given by Eq. (5) can be used in the canonical ensemble to evaluate the expectation value of any operator depending only on the coordinates.t3 For example, the expectation value of the potential energy is Freeman,and Doll: Extending J walking to quantum systems 5715 The quantum average is obtained in the limit that all the Fourier coefficients are included *im (f?kmBx+(ne Cax-+m For notational simplicity, we will omit the k,,, subscript from subsequent averages. Equation ( 13) has the form of a classical Monte Carlo average, differing only in the increased dimensionality of the integrations due to the introduction of the Fourier coefficients ak as additional integration variables. This addition of auxiliary degrees of freedom in the quantum expressions is characteristic of path integral treatments in general.
(8) 2. Kinetic energy estimators -P j-i Ux(u> ldu] (9) is the action weight factor. Since the quantum paths begin and end at the same point for properties dependent only on the coordinates, the free particle density reduces to a constant The calculation of the kinetic energy expectation value is complicated by the nondiagonality of the kinetic operator in coordinate representation. Several kinetic energy estimators have been developed, two of which are the T method14 (for temperature differentiation) and the H method15 (for Hamiltonian operating on the density matrix). The T method evaluates the internal energy U in its entirety using the statistical mechanical expression (15) (10) and so we denote the weight factor simply as w(x,a;P). The effects of the higher-order Fourier coefficients on the expectation value decrease with increasing k, and in practice, the infinite sums are truncated at some suitably selected maximum value of k denoted k,,, giving an approximate path where Q is the standard canonical partition function Q(T,V,iV>= dxp(x;P). I (16) Using JZq. (5) gives (17) x,(u) =x+ c (Ik sin(k?ru). k=l (11) The number of Fourier coefficients required for a given calculation depends on the system being studied, the temperature, and the level of accuracy desired. For notational convenience, we represent the average of some function along an approximate path as for the systems we consider, the T-method kinetic energy can be written explicitly as The expectation value for the potential then becomes (nL= s dx da wkmuha;P> V(x) .f dx da w&(x,a$) ' (13) where Another kinetic energy estimator closely related to the Tmethod estimator can be obtained by noting that in the limit km,, -+ 03, the term (( qx(u)]),-V(x)) in Eq. (18) (where the average is taken over all the paths) vanishes.16 The expression resulting from the neglect of this term at finite k,,, can also be obtained directly from the mass differentiation of Q and hence is denoted the M method (for mass differentiation) <Th=;(y-( z 3)).

Partial averaging
The H-method estimator can be obtained direc$y as the expectation value of the kinetic operator T= -(#/ 2m)V2, Since the number of Fourier coefficients required in the FPI calculations increases rapidly with decreasing temperature or as the system becomes highly quantum, it is useful to employ techniques that improve the convergence of the calculations with respect to k,,. One such method is partial averaging, " which takes into account the effects of the truncated Fourier coefficients of order higher than k,, that are neglected by the direct Fourier method. The system potential is replaced with an appropriate effective potential obtained as the Gaussian transform of the bare potential (T)H= S dx $~(x,x';P) 1 x=x) S dx p(x,x';P) ' (20) which gives (T)H=$- where Both the T and M methods have the advantage of being computationally simple to evaluate, but suffer the disadvantage of having a variance that grows linearly with k,,.
The H method is more computationally intensive, but its variance growth is usually much weaker, even decreasing with increasing k,, for some systems.

Heat capacity
As with kinetic energy estimators, there are several estimators for the constant volume heat capacity C, The estimators we used are analogous to the T-method kinetic estimator-they are derived by temperature differentiation of the energy (24) and hence are called T-method heat capacities. Inserting the T-method energy expression ( 17) into Eq. (24) gives the T-method heat capacity (Cv)T=kB [y-2( 12 $) +Bi(O--oi)], (25) where (26) while using the sum of the potential energy and the Hmethod kinetic energy given in Eq. (21) gives a hybrid TH-method heat capacity This leads to an analogous partial average action weight factor wPA(wd) =pfp@)exp vd[x(u) 1). .
For classical systems, the usual Metropolis method uses a random walker to sample the configuration space by making moves from an initial coordinate xi to a final coordinate xf with a probability of acceptance where (38) is the acceptance ratio, p(x) =27'exp[-@V(x)] is the Boltzmann distribution with Z the standard configuration integral, and T(x' Ix) is the transition matrix or sampling distribution.
The sampling distribution is usually generated from uniform deviates 6 over a finite stepsize range A to give'* T(x'lx)= I l/A for Ix'-xl <A/2 0 otherwise and (4-o) Attempted moves are generated according to x'=x+ (g -;)A with the maximum stepsize A/2 usually adjusted to give acceptance probabilities of approximately 50%. The required size decreases with increasing p. This temperature dependence can lead to quasiergodic behavior whenever the step size becomes too small relative to the potential barrier heights and widths; the walker becomes effectively trapped within a region of configuration space and is unable to sample the whole space in a representative manner. (39) where w(x,a;p) is given by Eq. ( 14), or by Eq. (30) if partial averaging is used. Since the coordinate moves are usually more important than the Fourier coefficient moves, the usual implementation of the transition matrix is to move a randomly selected Fourier coefficient simultaneously with each attempted coordinate move. Usually the selection scheme for the Fourier coefficients weights the lower-order coefficients more highly than the higher-order coefficients. While the coordinate moves use the typical Metropolis scheme of uniform deviates scaled by a suitable box size A,, the Gaussian-like distribution of the Fourier coefficients allows attempted coefficient moves to be sampled either from a Gaussian distribution of width ok or from uniform deviates scaled by box sizes A+. Both schemes work well, although sampling from a Gaussian distribution results in more rejected moves at lower temperatures. The J-walker method is based on the observation that Quantum J walking then couples this scheme with atthe Boltzmann distributions governing the Metropolis tempted jumps to a higher temperature distribution of walker are largely dependent on the form of the potential (xJ,aJ) configurations in a manner analogous to the claswith distribution maxima corresponding to potential min-sical case. Whenever a jump is attempted, the transition ima. The widths of the distribution peaks are inversely matrix used is simply the action weight factor for the dependent on B with higher temperatures resulting in higher temperature Frantz, Freeman, and Doll: Extending J walking to quantum systems 5717 wider distributions. In essence, higher temperature walkers are less constrained since their larger stepsizes allow them to overcome the potential barriers more effectively. For a Metropolis walk of a given length, there is a threshold temperature such that walks undertaken above the temperature are ergodic. J walking then occasionally replaces the usual Metropolis moves with an attempted jump to the position occupied by such a higher temperature walker (J walker). This is equivalent to replacing the usual Metropolis transition matrix given by Eq. (39) with the Boltzmann distribution at the higher temperature whenever a jump is attempted where PJ is the jump attempt probability. This gives an acceptance ratio of where fiJ is the J-walker temperature parameter. In the high temperature limit pJ -+O, the acceptance ratio reduces to the standard Metropolis expression given in Eq. (40).
Because the Boltzmann distribution broadens as the temperature increases, J walking in this limit reduces to simple jumping with a large stepsize AJ, and so the likelihood of a jump being accepted decreases sharply. In the limit /3,+fl, qJ(xJlx) + 1 since the low temperature walker is now effectively sampling from its own distribution.
For quantum systems, the extra degrees of freedom associated with the Fourier coefficients result in an acceptance ratio The J subscript on temperature-dependent parameters such as ukJ indicates they have had fi replaced with fl& Inserting Eq. (44) together with Eqs. ( 10) and ( 14) into Eq. (43) gives q+h,aJI x,a) =eAsQ, where AS', is the change in the action in moving from (x,a) to (xJ,aJ). For the direct FPI formulation, ~q,=hSe+BJC(f[x,(u);aZ,(Ir)l). , In the high temperature limit a,+O, the Fourier coefficient distribution width okJ vanishes and the acceptance ratio reduces to the Metropolis expression (with large step sizes implicit ), analogous to the classical case, so that the jump acceptance becomes negligibly small. In the limit #IJ-p, OkJ -+ CQ, and o$(u) +2(u), giving qJp(xJ,aJIx,a) + 1, so all attempted jumps are accepted since the low temperature walker is again sampling from its own distribution. The inverse dependence between /3 and the mass m as evident from the ratio B/m appearing in the expressions for parameters such as pfp, ak, and a(u) suggests alternative formulations of J walking where the temperature is held constant and the J walker mass decreased, or where both the J-walker temperature and mass are varied. We have not investigated these formulations, but only note in passing that they might be considered in those situations where J walking at higher temperatures is problematic.
The disadvantage of using quantum J-walker distributions is that the storage requirements are increased by a factor k,,. Since quasiergodicity results primarily from the large barriers associated with coordinate space, it seems plausible to use a distribution that is mostly dependent on the coordinates. Hence an alternative to using the quantum distribution for J walking is to use a classical J-walker distribution obtained from the Boltzmann distribution for coordinate moves and to then sample the Fourier moves from a normalized Gaussian distribution. The composite transition matrix for attempted jumps is then The prime on okJ signifies that the coordinate and Fourier coefficient J-walker temperatures need not be the same, and in fact the resultant expression for the change in the action can be greatly simplified by setting the Fourier coefficient J-walker temperature equal to that of the low temperature walker. This results in the cancellation of the Fourier coefficient contributions, giving for the direct Fourier path integral formulation qJc(xhaJIx,a) =flc, The limiting behavior of the acceptance ratio is not as straightforward as in the purely classical and purely quantum cases discussed previously. In the high temperature limit BJ-O, the J-walker effective potential ( V[X~(U)]),, reduces to the potential V(x,). Again because the distribution width increases with increasing temperature, the configuration x, will not likely correspond to a minimum and the probability of the jump being accepted is small. In the limit flJ+/?, the change in the action AS, (or tiCpA) will not be zero as was the case for the quantum distribution, and consequently the jump acceptance will not be unity. Physically, the low temperature walker is not sampling from its own distribution if it is sampling from a classical distribution at the same temperature. For higher temperatures, the quantum and classical distributions converge since the effective potentials approach the bare potentials and the partial average contributions become small, but for lower temperatures, these differences become appreciable. So in the limit PJ=8-0, q&wA x,a) -4 1, but in the low temperature limit flJ =B+CO, qJc(x,,,aJIx,a) -+ 0. Thus J walking from a classical distribution is ultimately limited in its use at low temperatures (or low particle masses).

III. QUASIERGODICITY IN THE QUANTUM DOUBLE WELL POTENTIAL
We illustrate the effects of quasiergodicity on quantum systems by considering a one-dimensional double well potential defined by the quartic polynomial This function has a minimum located at x = b with V(b) =0, a barrier of height h located at x=0, and a second minimum located at x= --a. By restricting a to the range O<a<b, the x=b well is the global minimum and the potential varies from a single well and barrier for a=0 to a symmetric double well for a= b. We can parametrize the degree of asymmetry in the double well for a given b and h by recasting a in terms of the ratio of the well depths For a given b and y, a can be obtained easily by iteration. This potential is a generalization of the double well potential we examined in our study of J walking on classical systems,' where the barrier height and global minimum were fixed at h = b = 1. For a particle of given mass m, the degree of quantum mechanical behavior in terms of spacings between quantized energy levels can be varied by adjusting the "box size" b and barrier height h. Alternatively, for a given potential, the degree of quantum behavior can be adjusted by varying the system mass. In order to compare the quantum results with the classical results obtained earlier, we have chosen the latter scheme and have used the previous classical potential with h and b fixed at unity and a variable. We have arbitrarily selected m = 1,O. 1, and 0.01 amu to model low, moderate, and high degrees of quantum behavior, respectively. Figure 1 shows the potentials for some representative values of y.

A. Exact results
The problems arising from quasiergodicity in these potentials can be seen by comparing their average energies obtained by FPI methods with the numerically exact solutions obtained directly from the eigenenergies. For systems at a defined temperature and volume, the average of a Hermitian operator P^ in energy representation is and using standard numerical techniques on the resulting matrix equations. A suitable oscillator frequency can be obtained in terms of the barrier height h with o = $%%. Since FPI methods can provide averages for both the potential and kinetic energies, we calculated expectation values for both the potential and kinetic operators (in amu) As a check, the eigenenergies were also obtained directly from the differential equations using common two point boundary value numerical methods." Figure 2 shows the energy levels for the y=O.8 potential for each of the three masses, while Fig. 3 shows the average internal energy as a function of p for the same potential, as well as the average classical energy for comparison.

FPI results
Using the FPI formulation, the average potential energy for a given potential at some temperature defined by /3 is given by Eq. ( 13 ) , while kinetic energies can be obtained with Eqs. (32)-( 34). Applying partial averaging to the double well potential given by Eq. (52) gives the effective potential correction term for Eq. ( 3 1) in terms of the polynomial coefficients Likewise, the T-and M-method correction terms in Eqs.
(32) and (33) are Partial averaging provided significant improvements in convergence for this potential and all results reported subsequently were obtained with partial averaging invoked. Hence for the sake of notational simplicity, the partial average subscripts will be omitted. tialized FPI Metropolis walks, each consisting of lo4 warm-up moves followed by lo4 moves with data accumulation. Each attempted move consisted of a coordinate move together with a Fourier coefficient move (sampled from a Gaussian distribution of width ak>, with the coefficient index selected from a Gaussian distribution of width k,,,/2 to give higher weight to the lower-order coefficients. The action integrals were evaluated using Simpson's rule with n quad = 16 quadrature points. In our study of quasiergodic behavior in classical systems we noted that the direction of the error resulting from the incomplete sampling depended on the walk initialization and that this property could be used as an indicator of quasiergodicity. For quasiergodic walks originating at the global minimum, the average value of the potential energy (v) is too low since the higher energy wells are insufficiently sampled, while walks originating from higher wells result in values that are too high. This is the case for quantum systems as well. Figure 4 shows the results obtained for (V) and for the H-method kinetic energy (T), for quantum particles with masses of m = 1,O. 1, and 0.01 amu. The plots show the variation with P for a fixed potential having y=O.9, and the variation with y for a fixed temperature corresponding to fi= 10. As in the classical case, there is a threshold temperature that marks the onset of quasiergodicity for walks of a given length. It occurs at the corresponding j3 point where the (V) curves representing global initialization diverge from the (V) curves representing random initialization. The threshold temperature depends slightly on the mass, being lowest for the m =O.Ol system. This can be attributed to the effective lowering of the energy barrier by the higher quantum zero point energy with the lower mass and by the increased importance of tunneling. The effects are particularly noticeable in the fI=lO, m=f, km?=1 . 7=0.9, m=l, km==1 1 ,,' I . , . , , , , 030 @+J 1.8, In each case, the smooth curves are the numerically exact solutions, while the Metropolis curves were obtained from 100 independently initialized walks consisting of 10" warm-up moves followed by lo4 moves with data accumulation, all with partial averaging invoked. Each attempted move consisted of a coordinate move sampled from a uniform distribution, together with one Fourier coefficient move sampled from a Gaussian distribution. The Fourier coefficient moved was also selected randomly from a Gaussian distribution. The number of Fourier coefficients needed to achieve convergence in the kinetic energy at p= 10 is given by k,,. Some representative single standard deviation error bars have been included. The potential energy curves exhibit the bifurcation characteristic of quasiergodicity in classical systems. In each plot, the lower branch corresponds to walks initialized at the global minimum x= 1, while the upper branch corresponds to walks randomly initialized. The kinetic energy curves exhibit no bifurcation. The H-method kinetic energies show no evidence of quasiergodicity and are in good agreement with the numerically exact values throughout. Similar results were obtained for the T method. Figure 5 compares the H-and T-method kinetic energies as functions of p for various k,,, for the y=O.9 potential for m=O.Ol. As expected, the convergence in k,, is significantly better for the H-method kinetic estimator than for the T method. While the T method has a lower variance at low k,,,, its variance increases much more rapidly with increasing k,, than does the H-method variance, so that its variance is comparable to the H method's in the convergence limit. The larger variance at higher k,,, can make quasiergodicity less of a problem for highly quantum systems since the longer walks required to reduce the variance in the kinetic energy will also increase the likelihood of the walker passing through the bottlenecks in configuration space.

C. J-walking results
The ability of J walking to eliminate the systematic errors arising from quasiergodicity in the quantum double well potential can be seen in Fig. 6, which shows results obtained using J walking from classical and quantum distributions. The plots depict the internal energies as functions of y for /3=10 for masses m=l and m=O.l. The jump attempt probability was P,=O.l. Comparison with the analogous FPI Metropolis results depicted in Fig. 4 shows that the large systematic errors in the average potential energy have been eliminated. The results were obtained in a similar manner as those in Fig. 4. Although we have depicted only four combinations, we ran calculations for k---1, 2, 4, 8, 16, and 32 for all three masses for J walking from both quantum and classical distributions. All combinations gave similarly good results, differing only with the increased noise at higher values of k,,,. The effects of various J-walking parameters on the quantum results are very similar to the effects on the analogous classical results we obtained in our earlier work,' and so we refer the reader to that paper for a detailed discussion of the effects of varying the jump attempt frequency and distribution parameters, and instead concentrate here on the differences between J walking from quantum distributions and classical distributions.
The quantum distributions consisted of lo5 configurations sampled every ten moves from individual high temperature FPI Metropolis walks at pJ=3. Each configuration consisted of the coordinate and its associated k,,, Fourier coefficients. An attempted jump simply required selecting one of the configurations randomly and evaluating the acceptance ratio given by Eq. (45). The computational overhead for the evaluation of the acceptance ratio was kept very low by also storing the Fourier coefficient term Xka2,/2o$, th e path averaged effective potential ( VIdmaJI and the partial average term mJo4;o2J(um, f or each configuration in the distribution. Consequently, the extra computational cost required to implement J walking from a quantum distribution was due mostly to the generation of the distribution from a long FPI Metropolis walk (the length being the product of distribution size and the sampling frequency). The classical distributions similarly consisted of lo5 configurations sampled every ten moves from individual high temperature classical Metropolis walks at flJ=3. However, each classical configuration consisted only of the coordinate, reducing the storage requirements by a factor of k-m compared to the quantum distributions. An attempted jump using a classical distribution consisted of randomly selecting a configuration from the coordinate distribution and generating each of the Fourier coefficients uk randomly from a Gaussian distribution of width ok. The attempted jump would then be accepted or rejected according to the acceptance ratio given by Eq. (49). The evaluation of the acceptance ratio was made more efficient by storing the potential energy for each configuration, but because all the Fourier coefficients were generated with each jump attempt, the effective potential and partial average terms also had to be generated. This resulted in slightly longer walk times compared to J walking from quantum distributions. However, since the classical distributions could be generated much more quickly than the quantum distributions, the use of classical distributions was more efficient overall than the use of quantum distributions.
Although useful for heuristic purposes, onedimensional potentials are of limited use as tests since methods that give goods results in one-dimensional systems can fail badly when applied to multidimensional systems. For example, simply using step sizes larger than the barrier width gives good results for the classical onedimensional double well potential since a relatively large  Fig. 4, the smooth curves represent the numerically exact solutions, while the J-walker curves were obtained from 100 randomly initialized walks consisting of lo4 warm-up moves followed by lo4 moves with data accumulation; partial averaging was used throughout. The J-walker distributions used consisted of 10' configurations sampled every ten moves from high temperature Metropolis walks at fi,=3. Classical distributions contained configurations obtained from individual classical Metropolis walks, while quantum distributions contained configurations and their associated k,, Fourier coefficients generated from individual Metropolis FPI walks. The jump attempt frequency was P,=O.l. number of attempted steps will land close to the other minimum, resulting in higher acceptance rates. The likelihood of landing near a minimum in a multidimensional space with a large random step is quite small though, and so the method is useless for the higher-dimensional spaces typically encountered in most Metropolis simulations. The great advantage of J walking is that the locales near the various potential minima are visited representatively by the high temperature walker and their corresponding conflgurations stored for future use. The major limitation in classical systems arises from the acceptance ratio depending exponentially on the difference between the temperature parameters [as is evident in Eq. (42)], so that the likelihood of a jump being accepted decreases for larger temperature differences. For quantum systems, the dependence of the jump acceptance probability on the temperature is more complicated, with J walking from classical distributions inherently limited at lower temperatures and J walking from quantum distributions practically limited at low temperatures because of the increased storage requirements with increased k,,. These limitations did not arise in the double well potential, and both quantum J-walking schemes gave equally good results, although J walking from classical distributions was still the preferred method for this system because of the lower storage requirements and shorter computational times required to generate the distributions. The differences between the use of classical and quantum distributions are more evident in multidimensional systems and a more thorough comparison of the two is given in the next section where we apply the methods to atomic clusters.

IV. J WALKING IN QUANTUM CLUSTERS
In our initial study of classical J walking,' we tested the method on rare gas clusters bound by pairwise additive LennardJones potentials. Evidence of quasiergodic behavior had been noted in previous Monte Carlo calcula-tions20'21 and molecular dynamics studies,22 and had been postulated to account for differences between Monte Carlo and molecular dynamics results for Ar13.23 The clusters' small sizes and consequent modest storage requirements and computational times made them ideal candidates for a realistic test of J walking. We found that quasiergodicity was indeed evident in the constant volume heat capacity (and to a lesser extent in the internal energy) over a significant part of the temperature range encompassing the transition from "solid-like" behavior to "liquid-like" behavior,24*25 and we showed that J walking was very successful in eliminating the problem.
The large majority of theoretical studies of the equilibrium and dynamical behavior of clusters to date have been classical. While classical mechanical methods are appropriate for many cluster studies, there are many interesting problems where quantum mechanical effects are significant enough to render classical mechanical studies inadequatea The low temperatures typically encountered in rare gas cluster simulations suggest that quantum effects can be significant in the transition regions. Recent studies26s27 on small argon and neon clusters found moderate quantum effects in the cluster melting region for argon and large quantum effects for neon. In particular, the large zeropoint fluctuations found for neon make the validity of even a qualitative classical description of the dynamical behavior of small neon clusters suspect. Contrast this with the bulk melting temperature of neon (23.48 K), which is quite close to the value obtained from applying the law of corresponding states to the melting temperature of bulk argon (83.95 K). The recent introduction of quantum simulation methods such as path integral Monte Carlo combined with the increased computational performance of modern workstation computers has made the study of many-body quantum mechanical systems feasible. Thus clusters are also ideal for testing quantum J-walking methods.

A. Metropolis simulations of classical Ne7 and Ar,
The comparison of quantum and classical Metropolis calculations for small clusters can help provide insight into the limitations of the classical studies, as well as indicate possible quasiergodic behavior in the quantum systems and so we begin with the classical results. The clusters were modeled by the usual pairwise additive Lennard-Jones potential v= 2 VLJ(@, icj (61) with E= 119.4 K and a=3405 A for argon, and e=35.60 K and a=2.749 A for neon." Small clusters are known to become unstable beyond a threshold temperature TB that varies with the cluster size n.28V29 For the Lennard-Jones potential under free volume conditions, the average energy vanishes in the limit of an infinite number of configurations. Consequently, the choice of boundary conditions can have a pronounced effect on some of the properties of small clusters.24 We have followed Lee, Barker, and Abraham3' and confined the clusters by a perfectly reflecting constraining potential of radius R, centered on the cluster's center of mass. The constraining radius was set to R,=4o for all runs. 5724 Frantz, Freeman, and Doll: Extending J walking to quantum systems The classical internal energy and heat capacity for an n-atom cluster are given by Curves for ( 17) and (C,) as functions of T were generated for argon over the temperature range from l<T<60 K with a mesh size of AT= 1 K using the usual Metropolis methods. The lowest energy configuration, the pentagonal bipyramid,3' was used to initialize the T= 1 K walk; subsequent temperatures used the final configuration from the preceding walk as their initial configuration. Each data point consisted of lo5 warm-up passes followed by IO' passes with data accumulation. Because classical Lennard-Jones clusters obey the principle of corresponding states, results for Ne7 were obtained by simply scaling the Ar, results.
B. J walking for classical Ne, and Ar, As we found in our earlier study of classical Art3, a single J-walker distribution generated at a temperature of FIG. 7. A comparison of J-walking and Metropolis results for classical Ar, and Ne,. The open circles were obtained from Metropolis walks consisting of 10' warm-up passes followed by 10' passes with data accumulation. Single standard deviation error bars have been included for those points where the standard deviation is larger than the symbol size. The curves were generated using J walking, with jumps attempted with a frequency P,=O.l from external J-walker distributions that each contained lo6 configurations sampled every 100 passes. The AI7 distributions were generated in stages at T~=50, 38, 24, and 14 K. The lower temperature distributions were obtained using J walking, with jumps attempted from the preceding distribution. Each J-walking point consists of 104 warm-up passes followed by 10' passes with data accumulation. The plot at left is the internal energy and the plot at right is the constant volume heat capacity; reduced units are used with lP = U/e and C$ = Cv/kB T=50 K is not useful at very low temperatures since the jump acceptance decreases sharply as the temperature difference between the low-temperature walker and the Jwalker distribution increases. Hence, to cover the entire temperature range, the J-walker distributions were generated in stages at T=50, 38, 24, and 14 K. The T=50 K distribution consisted of lo6 configurations sampled every 100 passes from a single Metropolis walk. This distribution was used for J walking over the temperature range 38<T<50 with a mesh size of AT=0.2 K. The J walks were similar to the Metropolis walks, except that jumps were attempted with a frequency of PJ=O. 1. The distribution was then used to generate a T=38 K distribution of the same size, which was then used for the next temperature range and then to generate the T=24 K distribution, etc. Figure 7 compares the Metropolis and J-walking results for the classical energy and heat capacity as functions of the temperature for the two systems. There is very good agreement between the J-walking and Metropolis results except in the low temperature side of the C, peak where the Metropolis results are very noisy and unconverged, indicating problems due to quasiergodicity. The shoulder in the heat capacity curve occurring at about 10-20 K for Ar7 corresponds to the transition region where Berry and co-workers32 found evidence of solid-like and liquid-like clusters coexisting. We graphically displayed samples of the Metropolis configurations in this region and found them to consist primarily of the four stable isomers and the compact transition forms reported by Berry. Despite the large configurational barriers reported in this region, we did not observe problems due to quasiergodicity. This is in sharp contrast to our earlier study of Ar13, which showed pronounced effects of quasiergodicity in the solid-liquid coexistence region. Quasiergodic effects are evident in the slightly higher temperature range from 20 to 40 K. When we graphically displayed samples of the Metropolis configurations from this region, we also noticed the compact isomers and transition forms appearing in the shoulder region, but in addition, we also saw a large percentage of loose random forms that could no longer be related to recognizable isomers. We also saw occasional dissociations consisting of a compact Ar6 core and a nearby lone atom. This was not the case in the high temperature side of the C, peak. Here, both Metropolis and J walking are well converged and in agreement with each other. This temperature range corresponds to the cluster dissociation region and many of the Metropolis configurations that we examined appeared consistent with a dense gas enclosed in a spherical container. The low temperature side of the C, peak appears then to be another coexistence region containing a wide variety of configurations. The discrepancy between Metropolis and J-walking results in this region also occurred in the analogous quantum curves and we present a more detailed examination later when we discuss the quantum results.

C. Metropolis FPI for quantum Ne, and Ar,
The FPI treatment of Lennard-Jones clusters has been detailed before6 and so we merely recap the pertinent points. The computational overhead relative to analogous classical Metropolis simulations increases rapidly as the number of Fourier coefficients k,,.,a, needed for convergence grows. There is a slight increase due to the additional dimensionality associated with the extra Fourier coefficients, but since a typical FPI move consists of a coordinate vector move together with one Fourier coefficient vector move (with a bias favoring the lower order coefficients), the computational cost is only moderate. The major cost is found in the calculation of the action integrals. Using Simpson's rule to evaluate the integrals, the total computation time has a roughly linear dependence on the number of quadrature points nquad used. Since the paths corresponding to higher k,, become increasingly irregular," nquad needs to be increased as the number of Fourier coefficients is increased. We determined the minimum number of quadrature points needed for a given k,, by increasing lts,& until the changes in the average energy and heat capacity were within the standard deviation. The values of rzquad that were found to be sufficient for each k,,, were k,,l, 2, 4, 8, 16, 32, 64, nq,,&-+ 4, 8, 8, 16, 32, 64. (64) In our preliminary investigations, we also found the number of Fourier coefficients needed for convergence was reduced substantially by using partial averaging within the gradient approximation,6 and so all the results reported here are with partial averaging invoked. Again for notational simplicity, we will leave out the partial average subscripts. Applying partial averaging to the second-order Taylor series expansion of the Lennard-Jones potential gives for the partial average terms for Eqs. (3 I)-( 36), Although we found the H-method kinetic energy estimator to be superior to the T-method estimator for the double well potential, we primarily used the T-method estimator for the cluster simulations. The computational requirements for the H method (especially for the heat capacity) were substantially greater than those for the T method, and since values of k,,, of 1 or 2 were found to be adequate for the partial averaged T-method kinetic energy and heat capacity over much of the temperature range investigated for each cluster, the faster convergence of the H method with k,,,, was not an issue. Another technical difficulty with the H method that disfavored its use at higher temperatures concerns the constraining potential. Unlike the T-method case, the evaluation of the H-method integrals requires a finite constraining potential, such as a high power polynomial V,(r) = [(r--R,,. )/Rd2', where R,,, is the center of mass of the cluster. The evaluation of such a potential adds even more computational overhead. We checked some of our T-method results by comparing them with their H-method counterparts and found excellent agreement. We obtained FPI Metropolis results for the potential, kinetic, and total energies, and for the heat capacity for Ar, over a temperature range from 1 to 60 K on a mesh size of AT= 1 K, and for Ne, over a range of 0.5-20 K on a mesh size of AT=1 for higher temperatures and AT=05 for lower temperatures. Each result was obtained from walks having a total length of lo7 passes, initialized from previously warmed configurations. We ran simulations for each point for k,,,,,= 1,2,4,... until the total energy and heat capacity converged, up to a maximum of k,,=32, with the number of quadrature points for each k,,, given by Eq. (64). Attempted moves for each pass were made by randomly displacing an atom's coordinate vector together with one of its Fourier coefficient vectors, with a bias in favor of the lower order coefficients. Both coordinate and Fourier coefficient displacements were sampled from uniform distributions with box sizes fixed to give approximately 50% acceptance ratios for the combined moves. aThe J-walker distributions were each generated from a single walk with configurations stored every 100 passes. Those distributions generated using J walking have their corresponding J-walker temperature listed in the T, column (quantum distributions generated with J walking used classical J-walker distributions). The columns labeled Size contain the total number of configurations stored. Classical configurations consisted of coordinates and the cluster potential energy, while quantum configurations consisted of coordinates, Fourier coefficients, and various cluster energies.
D. FPI J walking for quantum Ne, and Ar, We ran quantum J-walking calculations using both quantum and classical distributions for both Ar, and Ne,. The two methods gave statistically identical results for each cluster over the common temperature ranges and they both gave good agreement with FPI Metropolis over most of the temperature range, differing only in the low temperature side of the Cy peak, as in the classical case. Classical Ar, distributions were generated in stages at temperatures of T,=50,38,24, and 14 K as described in Sec. IV A, with classical J walking used to generate the lower temperature distributions from the preceding ones. Quantum Ar, distributions with increasing k,,, were generated at temperatures of T,=60, 40, 25, 15, and 9 K. Similarly, classical The data for each cluster were obtained in the same manner as was the Metropolis FPI data, except that the total walk length was only lo6 passes instead of lo7 passes; we found the J-walking data to typically have a much lower variance than the Metropolis data, despite the shorter walk length. The jump attempt probability was the same as in the classical J walking PJ=O. 1. For J walking from classical distributions, an attempted jump consisted of randomly selecting a cluster coordinate configuration from the distribution and randomly generating each of the Fourier coefficients for each coordinate from a Gaussian distribution. This allowed us to run simulations for each point with a complete km,, series for full comparison with the Metropolis results. For J walking from quantum distributions, an attempted jump consisted of randomly selecting a cluster configuration containing both the coordinates and their attendant k,,, Fourier coefficients. This required the low temperature walker to also have the same value of km, and so the convergence in k,, could not be compared to Metropolis and J walking using classical distributions without generating separate quantum distributions for each k,,,. Because the Metropolis FPI and the J walking from classical distributions showed similar convergence and because generating quantum distributions was very time consuming, we simply took the results obtained using the classical distributions to determine the converged value of k,, needed for each quantum distribution. We The open symbols were obtained from FPI Metropolis walks consisting of 10' passes initialized from previously warmed configurations. Single standard deviation error bars have been included for those points where the standard deviation is larger than the symbol size. The curves were generated using J walking from externally stored J-walker distributions with jumps attempted with a frequency P,=O.l. The distribution parameters are listed in Table I. The J-walking walk length was only lo6 passes. The plot at left is for Ne7 using J walking from quantum distributions; the plot at right is for Ar, using J walking from classical distributions. For J walking using classical distributions, the values of km were the same as for the corresponding FPI Metropolis values, but for J walking using quantum distributions, the values were fixed by the distribution k,, values (listed in Table I).
also checked the sensitivity of J-walking results on the distributions by extending the temperature ranges to overlap with other distributions. For example, the Ar7 T,= 50 classical distribution had a useful range of 35<T<60 K, while the T,=38 K distribution had a useful range of 20<T<45 K. Data obtained using each distribution agreed over the common range 35<T<45 K. A comparison of J-walking and Metropolis results for the total energy, kinetic energy, and heat capacity is shown in Figs. 8,9,and 10,respectively. Because Ne, shows larger quantum effects than Ar7, we have shown Ne, re- O'.. sults obtained using quantum distributions and Ar, results obtained using classical distributions. The results for Ne, obtained using classical distributions and for Ar7 obtained using quantum distributions were similar. The plots show very good agreement between J walking and Metropolis over most of the temperature range. As in the cases of classical Ar7 and Ne,, discrepancies occur in the temperature range corresponding to the low temperature side of the heat capacity peak. The Metropolis data in this region is very noisy for both the total energy and the heat capacity, and generally lies above the J-walker curves. In an attempt to account for these discrepancies, we checked the convergence for both classical and quantum Ar7 at 35 K for total walk lengths of lo4 to lo8 passes, using a classical distribution at T,=38 K for each. The results are summarized in Table II. For both the classical and quantum Jwalking results, the heat capacity and energy both show the l/,/N decrease in the standard deviation expected for asymptotically convergent stochastic processes, while the Metropolis results are clearly not converged. Quantum Jwalking results obtained using a quantum distribution at Tp=40 K with km,= 1 gave good agreement with the results obtained using the classical distribution. These results also agreed with J-walking results obtained using the T, =50 K classical distribution. A similar analysis of quantum Ne, at T= 10 K gave similar results-the Metropolis results were noisy and unconverged, while the J-walker results from both classical and quantum distributions at different temperatures were properly converged and gave good agreement with each other.
The T-method kinetic energies obtained using Metropolis did not show the same convergence difficulties, indicating the problem is isolated mostly in the potential energy. There is a slight discrepancy between Ne7 Metropolis and J-walking results at about 9-l 1 K, indicating perhaps some effects of quasiergodicity in the kinetic energy as well, but there could also be effects arising from incomplete convergence in k,,,. Different properties converge with k,,, at different rates, and the kinetic energy is more slowly con- Metropolis vergent than the total energy.6 The J-walker curve was fixed at k,,= 2 because the quantum distribution had been generated at that value (which was adequate for the total energy and heat capacity over that temperature range). The k,,,= 4 J-walking results obtained using a classical distribution gave good agreement with the Metropolis results. Figure 11 compares the jump acceptance rates for J walking using quantum and classical distributions. For J walking from quantum distributions, the jump acceptance for each distribution decreases from 100% (for T= T,), approaching 0% as the temperature difference becomes  Table I. large. The steepness of the decrease is greater at lower temperatures, reflecting the narrowing of the distributions with decreasing temperature. The useful range of the Tq = 15 K distribution ( k,, = 1) is from 15 down to about 10 K, while the useful range of the Tq= 6 K distribution (k,,, =8) is only from 6 to 4 K. In addition, the number of Fourier coefficients required for convergence increases rapidly at low temperatures, further limiting the useful temperature range of a low quantum J-walker distribution since the low temperature walker km,, is fixed at the same value as the J-walker distribution. Because the distribution size and the time required to generate the distribution grow rapidly with increasing k,,, J walking using quantum distributions becomes practically unfeasible at very low temperatures. For J walking from classical distributions, the maximum jump acceptance for a given distribution was less than lOO%, even when T= To since the classical configurations and generated Fourier coefficients are different than the corresponding quantum configurations. As the temperature becomes lower and the system becomes more quantum, the difference becomes larger. The T,= 15 K classical distribution had a useful temperature range of about 20-I 1 K, with a maximum jump acceptance of about 62%, but the T,= 11 K distribution had a useful range of 12-8 K with a maximum jump acceptance of only 7%. Hence classical distributions have fundamental limitations at low temperatures where the Fourier coefficients become strongly coupled to the coordinates. Classical distributions have the advantage of allowing variable km,. Although only curves for converged values of km, are shown in Fig.  11, the jump acceptance did not vary much with k,,. Figure 12 compares the quantum and classical energies and heat capacities for Ne, and Ar,. The curves are comprised mostly of J-walking data, with Metropolis results used for those temperature regions where J-walking data had not been obtained. The quantum neon Cy curve is shifted substantially to lower temperatures, with the peak much reduced relative to the classical curve, and the U curve is much higher than the classical energy throughout. This indicates that quantum effects for Ne, are considerable over the entire temperature range, including the cluster dissociation region. The shoulder in the classical C, curve corresponding to the liquid-solid coexistence region is completely absent in the quantum curve. This is a consequence of the quantum zero-point motion and tunneling, which effectively raise the minima and lower the saddles of the multidimensional potential surface and allow for easier isomerization.26 These results are consistent with classical and quantum quench studies,27 which showed significant differences between classical and quantum quenched isomer population distributions. The shoulder is also largely diminished in the quantum Ar, Cy curve, and there is a small, but significant difference between the classical and quantum argon curves for higher temperatures up to the cluster dissociation region. Clearly, it is necessary to incorporate quantum effects for a proper understanding of dynamical and equilibrium behavior of small argon and neon clusters.

V. DISCUSSION
Path integral methods form an important set of tools for extending Monte Carlo methods to quantum manybody systems. Because of the greater computational requirements incurred by quantum simulations relative to similar classical Metropolis simulations, ensuring proper convergence and developing methods that increase the rate of convergence are especially important. Two sources of slow convergence are inherent in FPI Metropolis simulations. The first is due to the increased number of Fourier coefficients needed for highly quantum systems at low temperatures. This increases the effective dimensionality of the simulation, as well as the computational overhead associated with the numerical calculation of the action integrals. The effect is purely quantum; classical simulations converge more quickly as the temperature becomes very low, but quantum simulations become increasingly difficult. The second source of slow convergence is due primarily to bottlenecks in configurational space and is common to classical simulations as well. While the effects of quasiergodicity in quantum systems appear to be diminished somewhat in comparison to analogous classical systems, they can still be quite formidable. Quantum J walking is a very useful technique for substantially reducing errors arising from quasiergodicity and increasing the rate of convergence, thus extending the power and scope of FPI methods. The method is related closely to classical J walking and shares many of its features. J walking is essentially a variant of the usual Metropolis algorithm, and so is very easy to incorporate into existing code. Quantum J walking is slightly more complicated to implement than classical J walking, but it offers more flexibility. Two variations-J walking using fully quantum distributions and J walking using classical coordinate distributions together with Gaussian distributed Fourier coefficients-were shown to be successful over temperature ranges where both quasiergodicity was evident and quantum effects were substantial. Both variations are completely compatible with various FPI techniques such as partial averaging, as well as the many kinetic energy and heat capacity estimators that have been developed. Although J walking showed better convergence than FPI Metropolis even at lower temperatures where quasiergodicity was not an issue, its limitations at very low temperatures preclude its general use in this region for overcoming the slow convergence that is a consequence of including large numbers of Fourier coefficients; the method is best suited for handling quasiergodicity in configurational space (these limitations are inherent in scalar computers-this would not be the case with parallel computers).
Because of the close kinship between quantum and classical J walking, many of the caveats we originally discussed remain valid. Care must be taken in properly selecting the J-walking temperature. It must be high enough to ensure a fully ergodic distribution, but low enough to provide sufficient jump acceptance rates. For the multidimensional systems we have examined to date, both criteria cannot be met with a single distribution and multistage distributions generated at several temperatures over the required range are necessary. Quasiergodicity in the lower temperature distributions can be eliminated with J walking, using the previously generated distributions for sampling. The walks used to generate the distributions need to be sufficiently lengthy to obtain representative samples of the configuration and Fourier coefficient spaces. This is especially important if the J-walker distribution is generated at a temperature where quasiergodicity is substantial. The results obtained using such distributions can be checked for consistency by generating additional distributions at the same temperature or at nearby temperatures. Quantum J-walking results can be checked further by us-ing both quantum and classical distributions for sampling.
To keep distribution sizes within the limitations imposed by finite computer resources,33 configurations can be sampled periodically, say every ten or 100 passes. This also helps reduce correlations in the distributions. Since most computers have much larger disk storage than memory storage, the large distributions can be stored as a collection of several files, each of which can then be loaded randomly into the computer's primary memory in a periodic manner as required. Writing the files in a parallel manner while the distribution is being generated also helps reduce the correlations within each file. Optical storage devices with capacities of hundreds of megabytes are well suited for J walking. Fully read-w&able devices are commercially available at prices roughly comparable to fixed magnetic disk drives and with access times only slightly slower.34 Although each cartridge can typically store only a few hundred megabytes per side, very large distributions can be stored on several cartridges, greatly increasing the effective storage capacity of the system. At present, we still do not know the necessary criteria for determining the minimum size required for a given distribution to ensure that it is a representative sample of the configuration space. Consequently, we have tended to make the distributions as large as possible to minimize systematic errors due to unrepresentative distributions.35 classical distributions. For Ne7, J walking using quantum distributions gave good results for the solid-liquid coexistence region (corresponding to the shoulder in the Cy curve), a temperature range not accessible using classical distributions, but this was also a region where standard FPI Metropolis gave sufficiently good results.

ACKNOWLEDGMENTS
One of us (DDF) wishes to thank the University of Lethbridge and the Department of Chemistry for their generous support of this research. One of us (DLF) acknowledges the Donors of the Petroleum Research Fund of the American Chemical Society for partial support of this work. We also thank Dr. Gustav0 Lopez, Dr. David Leitner, and Dr. Steve Rick for stimulating discussions.
For the systems we examined, J walking using classical distributions was preferable to the use of quantum distributions. These distributions could be generated much more quickly than their quantum counterparts, and since they consumed far less storage per configuration, we could include more configurations in each distribution. The classical distributions can of course be used for J-walking studies of the corresponding classical systems, allowing for comparison between the classical and quantum systems. J walking from classical distributions also has the advantage of being capable of handling arbitrary Jr,,, so that k, convergence can be checked easily for each temperature. J walking from quantum distributions fixes the lowtemperature walker's k,, value to the distribution's value, requiring that the converged value of km over the useful temperature range of the distribution be ascertained before the distribution is generated. Classical distributions are more limited than quantum distributions for handling very low temperatures. For example, for Ne,, classical distributions were useful down to about 8 K, while quantum distributions were useful down to 4 K. However, much of the quasiergodic behavior evident in multiparticle systems such as clusters is a consequence of bottlenecks in the coordinate subspace and becomes evident only at higher energies. Consequently, the corresponding temperatures are usually high enough to allow the use of classical distributions. Even in the case of Ne,, where quantum effects were quite large over much of the temperature range, J walking using classical distributions was just as effective as J walking using quantum distributions in eliminating quasiergodicity. The use of quantum distributions is generally preferable for handling very low temperatures or for providing an independent check of the results obtained from