Quantum Melting and Isotope Effects from Diffusion Monte Carlo

Subscriber access provided by UNIV OF NEWCASTLE

Article

Quantum melting and isotope effects from diffusion Monte Carlo studies of p-H2 clusters. Quantum Melting and Isotope Effects From Diffusion Monte Carlo Studies of -H Clusters. p

2

Joel D. Mallory, and Vladimir A. Mandelshtam J. Phys. Chem. A, Just Accepted Manuscript • DOI: 10.1021/acs.jpca.7b06649 • Publication Date (Web): 25 Jul 2017 Downloaded from http://pubs.acs.org on July 30, 2017

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a free service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are accessible to all readers and citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

The Journal of Physical Chemistry A is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 10

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Quantum Melting and Isotope Effects from Diffusion Monte Carlo Studies of p-H2 Clusters. Joel D. Mallory and Vladimir A. Mandelshtam∗ Department of Chemistry, University of California, Irvine, California 92697, USA ABSTRACT We present a rigorous characterization of the ground state structures of p-H2 clusters and their isotopologues using diffusion Monte Carlo combined with the inherent structures analysis. For the N = 19 cluster we explore the effect of “quantum melting” by quantifying the contributions of local minima to the ground state as a function of continuously varying particle mass. Doubling the cluster size leads to an enormous increase of its complexity: the ground state of (p-H2 )38 is highly delocalized over a large number of minima representing all the funnels of the potential energy surface. The ground state of (o-D2 )38 is also delocalized, but over a smaller subset of minima, which exclusively belong to the same disordered motif.

INTRODUCTION

Small p-H2 clusters have inspired significant interest in the field of condensed-matter physics due to their possibility for superfluidity at low temperatures.1–21 The related physical characteristics that have been studied intensely but that have proven challenging to pin down definitively even on a qualitative level were the structural motifs and phases of p-H2 clusters22–29 . Moreover, the simplest Lennard-Jones (LJ) potential was often used to model the p-H2 clusters, its isotopologues, as well as clusters of rare gas atoms,30–36 with the de Boer quantum delocalization parameter Λ being used to distinguish the LJ systems with different degrees of quantum strength. An extensive body of literature is devoted to elucidating the structural properties of these systems in their ground states using methodologies such as path integral ground state (PIGS)27–29,37–39 or diffusion Monte Carlo (DMC).5,7,18,19,26,40–43 In addition, quantum effects on equilibrium properties of small LJ clusters have been extensively studied as a function of the de Boer parameter Λ at finite temperatures with path integral Monte Carlo (PIMC)44–48 or with the harmonic superposition method,31,32 or using the variational Gaussian wavepackets.33–36 The previous ground state studies involving PIGS or DMC focused primarily on ascertaining the phases of pure p-H2 or o-D2 systems as a function of cluster size using either the Silvera-Goldman49 (SG) or the Buck50 potential energy surfaces (PES). These works relied heavily on the analysis of chemical potentials (i.e., the energy difference, µ(N ) = E0 (N ) − E0 (N − 1), between the clusters of consecutive sizes) to address the structural properties of these systems.26–29,40–42 For example, the cluster size N at which µ(N ) has a peak would be associated with a“magic number” implying that the corresponding cluster is especially stable. In addition, the

∗

[email protected]

physical state of a p-H2 system was often inferred from the shape of its radial density profile ρ(r)26–29 , from Lindemann criteria,28,29 or even from ground state energies determined by running importance-sampled DMC simulations using “liquid” or “solid” trial wavefunctions.26 However, none of the aforementioned quantities standing alone or reported together could provide the necessary structural information that must be robust enough to enable a correct characterization of the ground state wavefunction. Most importantly, the structural assessment can be considered reliable only upon assuming that the results are converged. However, the results of the present work will make it obvious that these systems are significantly more challenging than what has been assumed so far, inevitably requiring sophisticated sampling techniques and/or enormous computational resources due to the need to sample millions of local energy minima of the PES that all contribute to the system’s ground state. It should be noted though that Ref. 51 tried to make a comparative assessment of the convergence of PIGS and DMC in an attempt to resolve the discrepancies between the earlier DMC26,40,41 and PIGS27,28 results, which for relatively large clusters (e.g., with N > 26) disagreed even qualitatively. For example, the DMC chemical potential µ(N ) turned out to be a smooth function of the cluster size N , while µ(N ) computed by PIGS showed pronounced structure. Based on calculations for (p-H2 )48 Ref. 51 concluded that the DMC results were poorly converged, in particular, due to the severe bias with respect to the random walker population NW . Indeed, the apparent convergence of a DMC simulation for relatively small populations NW , especially when importance sampling is implemented (as was always the case in all the previous DMC studies of these systems), can often be misleading such that the use of the so called “solid” trial wavefunction, localized in any a priori chosen minimum of the PES26 would be unjustified. Moreover, although the authors of Refs. 27–29 posit that the PIGS results were well converged, this claim was based exclusively on the numerical behavior of the computed ground state energy as a function of the time step and the projection time

ACS Paragon Plus Environment


Page 2 of 10 2

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

rather than on a definitive knowledge of a large number of configurations contributing appreciably to the ground state wavefunction itself, which is delocalized over thousands (and perhaps even millions) of local minima. In order to obtain an accurate qualitative picture of the ground state so as to characterize it as “solid-like” or “liquid-like,” one must examine its structural properties by either directly exploring the configurations arising during the simulation (i.e., the configurations contributing to the ground state wavefunction) or by computing the relevant observables such as various spatial correlation functions or order parameters. Incidentally, both DMC and PIGS studies18,26–29,40,42 often deduced the phase of the ground state from analysis of the structure of the radial density profile ρ(r). We note though that two different configurations belonging to the same structural motif (e.g., anti-Mackay) would likely have noticeably different ρ(r) (especially near the origin r = 0) due to the sensitivity of this observable to the position of the cluster center of mass. Our present results (see below) indicate that the ground state wavefunctions for both p-H2 and o-D2 clusters are delocalized over many more than one minimum, which should result in density profiles that are difficult to interpret. Our ultimate objective is to characterize the structural properties of the ground state wavefunctions of pH2 clusters and their isotopologues. We employ DMC as originally proposed by Anderson,52,53 i.e., without importance sampling. This well-established implementation of DMC has been applied extensively by others and by us to determine the ground state properties of a wide variety of transient or “floppy” systems with weak 55 56 57 interactions,54 namely, CH+ , H+ , H+ , small HCl 5 5 7 58 59–65 clusters , and small water clusters. This variant of DMC uses a population of NW random walkers with equivalent weights that samples the configuration space bounded by a PES and collectively represents the wavefunction of the many-body system at projection time τ . At each time step of length ∆τ , some random walkers are replicated and some are killed to keep the random walker population virtually constant, and hence, to stabilize the average estimate of the ground state energy hEref i. At sufficiently long τ the ensemble of random walkers approaches a stationary distribution and the instantaneous energy Eref (τ ) fluctuates in a predictable way about its average value hEref i. In the limit of τ → ∞, ∆τ → 0, and NW → ∞, the random walker distribution converges to the exact ground state wavefunction with hEref i = E0 , i.e., the exact ground state energy. While we acknowledge that importance sampling has been used many times in DMC studies of (p-H2 )N and related systems5,18,19,26,40–43,51 and is an appealing approach for reducing the statistical errors in DMC calculations,26,41,66 the trial wavefunction introduces another source of bias into the simulations, especially if it does not overlap appreciably with the true ground state wavefunction, which for a floppy system as complex as a (p-H2 )N cluster is hard to avoid.

We use the SG PES,49 which is the most popular model utilized in previous studies.26–28,41 However, we believe that using other available models38,39,50,67 would unlikely affect our results qualitatively. (In fact, our preliminary calculations involving the LJ19 cluster resulted in virtually the same quantitative results.) For each SGN system, besides the usual calculation of the ground state energy, we utilize the descendant weighting (DW) method54,68,69 to compute the pair correlation function g(r) as a correctly defined observable. The ultimate aim of the DW method is to obtain a second copy of the wavefunction such that observables can be correctly averaged over the true density |Ψ(r)|2 as opposed to simply over the wavefunction Ψ(r) alone (i.e. the ensemble of random walkers). To this end, the progeny of each random walker are counted over a pre-determined projection time into the future τDW to determine the respective weights of the walkers. Random walkers that give birth to a large number of progeny have large weights while those that are more likely to be killed have negligible or zero weights. In addition, in the spirit of the inherent structures method,70 during each DMC simulation the random walker population would occasionally be quenched using the conjugate gradient method. The quenched configurations were assigned by comparing them with the minimized structures stored in the isomer library. Thus, we computed the corresponding isomer fractions, which quantify the contribution of the minima of the PES to the ground state wavefunction Ψ. Unfortunately, perhaps because only a small fraction of the random walker population accumulates very large weights for relatively long τDW (i.e. the DW projection time), we were unable to succeed in obtaining unambiguous (i.e., converged with respect to τDW ) results for the isomer fractions using the aforementioned DW method. Consequently, we considered the most natural quantity for the isomer fraction within the DMC framework, i.e. Ψ-averaged weights fk , computed by simply counting the number of quenched configurations arriving at each minimum of the PES.

RESULTS N=19

We start with the relatively simple, albeit still physically interesting, case of SG19 . This system is often identified as a so-called “magic number” cluster due to the double icosahedral (D5h ) structure of its global energy minimum26–29 which makes it relatively more stable than other clusters with similar sizes. Yet, we believe that SG19 represents a generic case for small SGN clusters undergoing the “quantum melting” transition. In order to better characterize this effect we performed a series of 64 calculations by regularly changing the particle mass in the range m ∈ [1; 15] amu, where m ≈ 2, 4, and 6 amu would correspond to p-H2 , o-D2 , and p-T2 , respectively. However, we report our results as a function of the


Page 3 of 10

The Journal of Physical Chemistry 3

FIG. 1. DMC ground state energies per particle E0 /N for SG19 as a function of the time step ∆τ for fixed NW = 2.5 × 106 (left column) and as a function of the reciprocal population size 1/NW for fixed ∆τ = 200 au (right column) using three different values of the effective quantum parameter λ = m−1/2 . All simulations were carried out to a projection time of τmax = 2 × 107 au. The statistical errors are on the order of 0.01 K per particle or lower.

E0/N (K)

-58.92

E0/N (K)

-58.80

λ=0.40

-58.96

-58.85

-59.00

-58.90

-59.04

-58.95

-40.46

-40.05

λ=0.55

-40.53

-40.20

-40.60

-40.35 -40.50

-40.67

-25.95 E0/N (K)

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

-25.85

λ=0.70

-26.10

-25.88

-26.25

-25.90 -25.93

-26.40 0

200

400 600 ∆τ (au)

800

1000

more relevant effective quantum parameter λ = m−1/2 . Note here that in her studies using PIMC calculations Chakravarty44 demonstrated the effect of quantum melting for LJ13 and LJ19 . In her studies she relied on the behavior of the Lindemann index and the equilibrium occupation of the global minimum as a function of the de Boer parameter Λ. Also note Refs. 35 and 36, in which the whole (albeit approximate) N − Λ phase diagram of LJN clusters was produced using the variational Gaussian wavepacket method. In order to establish the behavior and extent of the systematic error in our DMC simulations, we performed bias studies with the projection time τmax = 2 × 107 au using different magnitudes of the time step ∆τ and the random walker population NW for three representative values of λ. Fig. 1 demonstrates that the time step bias is small for ∆τ = 200 au (i.e. on the order of ∼ 0.01 K

0.0

0.5

1.0 6 10 /NW

1.5

2.0

per particle for the largest value of λ) upon comparison to the value extrapolated to the ∆τ → 0 limit. (The extrapolation was done by fitting the E(∆τ ) data points by a fifth-order polynomial with the constraint dE/dτ = 0 at ∆τ = 0.) Notably, the ∆τ bias also becomes apparently larger with increasing λ, as for a fixed time step ∆τ the error in the split-operator approximation, as employed in DMC, increases with the increase of the kinetic energy term. As seen in the right column of Fig. 1 the random walker population bias for the ground state energy is on the order of ∼ 0.1 K per particle for the largest population considered, i.e., NW = 2.5 × 106 . Although not very small, as apparent from the figure, it is under control, i.e, the E(NW ) curve can be extrapolated reliably to the NW → ∞ limit. Moreover, the bias of this size is unlikely a cause for concern for the following structural analysis of these systems. In addition, we also determined that



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 ACS Paragon Plus Environment

Page 4 of 10

Page 5 of 10


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

P TABLE 1. DMC ground state energies E0 and the entropy S0 = − k fk lnfk for different values of the projection time τmax and the random walker population NW for the N = 19 and N = 38 systems studied in this work. All DMC runs were performed with a time step of ∆τ = 200 au, and the statistical errors are given in parentheses. DMC or PIGS data from Refs. 26–28,41 are reported in the last column. System τmax (au) NW S0 /N E0 /N (error) (K) E0 /N [ref.] (p-H2 )19 2 × 107 2.5 × 106 0.40 -25.61 (0.008) -25.57626 , -25.59541 , -25.0327 , -25.7328 (o-D2 )19 2 × 107 2.5 × 106 0.08 -46.47 (0.01) -45.7928 6 8 (p-H2 )38 1.5 × 10 1 × 10 0.40 -34.47 (0.03) 8.5 × 105 2 × 108 0.40 -34.41 (0.04) -34.53126 , -34.53241 , -34.8228 6 × 105 3 × 108 0.40 -34.58 (0.05) (o-D2 )38 1.4 × 106 1 × 108 0.21 -58.73 (0.10) 5.5 × 105 2 × 108 0.19 -58.52 (0.07) -59.0428 5 8 6 × 10 3 × 10 0.20 -58.80 (0.15)

N=38

We now turn the presentation to the case of the N = 38 system, which has been the focus of numerous publications in the past,31–34,71–81 although most of the previous work considered the classical LJ38 cluster as even without quantum effects it already presented a great challenge for numerical simulations. Also, the two pair potentials (i.e., LJ and SG) are very similar, resulting in nearly the same properties. Due to its double-funnel potential energy landscape,32,71–74 corresponding to the octahedral (Oh ) global and icosahedral (C5v ) local minima, but with the C5v minimum stabilized by the entropic effects, this system can undergo a solid-solid transition induced by either thermal or quantum fluctuations.32 Moreover, the icosahedral motif appears in two forms, namely Mackay and anti-Mackay34,78 that differ by the structure of the overlayer. The lowest in energy anti-Mackay minimum of the PES is higher than the lowest Mackay minimum, yet the former has higher vibrational entropy, and can in turn be stabilized by either thermal or quantum fluctuations. Consequently, as demonstrated in Refs. 33,35, and 36 neither non-icosahedral nor Mackay icosahedral motifs are energetically stable for Lennard Jones clusters at sufficiently high values of the de Boer quantum delocalization parameter Λ. In particular, the ground state of the weakly quantum Ne38 (Λ ≈ 0.095) already has an anti-Mackay structure.33 Therefore, we expect the more quantum (p-H2 )38 (Λ ≈ 0.28) and (o-D2 )38 (Λ ≈ 0.20) to have disordered (or liquid-like) structures. While the systematic error in a DMC calculation is decreased with increasing random walker population NW , the statistical error for a fixed time step ∆τ is roughly proportional to (NW τmax )−1/2 . Therefore, given the total computer time, in order to minimize the total error, which includes both the statistical and systematic errors, it is preferable to utilize the largest population NW that can be fit into the fast computer memory. Thus, for both (p-H2 )38 and (o-D2 )38 we performed three independent DMC calculations using very large random walker populations NW = 1 × 108 , 2 × 108 , 3 × 108 . Each single calculation on a 64-core processor took on the order of 21 days to reach the projection times reported in the table

for the N = 38 clusters. The results for the ground state energies E0 and entropies S0 are summarized in Table. 1 along with the NW values used in the DMC simulations. We note that the projection times τmax for the N = 38 clusters are more than an order of magnitude smaller than those for N = 19 and that running the simulations to longer τmax would further decrease the statistical error. However, the large random walker populations keep the statistical errors for all of the energies on the order of ∼ 0.1 K per particle or lower. All results for E0 agree within the statistical errors and are consistent with the importance-sampled DMC26,41 and PIGS27,28 estimates reported previously, but we do not assume that our E0 values, particularly those for the N = 38 clusters, are more or less accurate than those reported by others. Rather, we emphasize that the value of the ground state energy E0 alone and its associated statistical error for any particular cluster may be the least interesting quantity and is a poor means by which to determine the convergence of the simulation. The fact of the apparent convergence (or lack thereof) for E0 cannot guarantee the convergence (or non-convergence) for other properties, particularly, the structural properties defined by the ground state wavefunction itself (e.g., g(r) distributions and isomer fractions). In fact, in none of the previous publications18,26–29,41 a correct structural assessment of either (p-H2 )38 or (o-D2 )38 has been reported. This could be the consequence of not using an adequate structural analysis of the ground state wavefunction or poor numerical convergence (e.g., associated with possible sampling problems when using path integral-based approaches) and/or bias introduced by the importance sampling (as was always implemented in the previous DMC calculations). To this end Fig. 5 shows the isomer fractions fk versus the orientational bond order parameter Q6 for the two systems. For (p-H2 )38 the random walker population, and thus the ground state wavefunction, is delocalized over a very large number of minima. The broad range of the bond order parameters, Q6 ∈ [0; 0.6], indicates that essentially all the structural motifs are represented by the ground state of (p-H2 )38 , including the octahedral, Mackay, and anti-Mackay motifs. In contrast, the situation becomes somewhat more




Page 6 of 10

Page 7 of 10


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

ACKNOWLEDGMENTS

grateful to Craig Martens for his useful comments and suggestions.

This work was supported by the National Science Foundation (NSF) Grant No. CHE-1566334. We are

REFERENCES

[1] P. Sindzingre, D. M. Ceperley, and M. L. Klein, “Superfluidity in clusters of p-H2 molecules,” Phys. Rev. Lett. 67, 1871 (1991). [2] D. Scharf, M. L. Klein, and G. J. Martyna, “Pathintegral Monte Carlo studies of para-hydrogen clusters,” J. Chem. Phys. 97, 3590 (1992). [3] D. Scharf, G. J. Martyna, and M. L. Klein, “Isotope effect on the melting of para-hydrogen and orthodeuterium clusters,” Chem. Phys. Lett. 197, 231 (1992). [4] M. A. McMahon, R. N. Barnett, and K. B. Whaley, “Rotational excitations of quantum liquid clusters: He7 and (H2 )7 ,” J. Chem. Phys. 99, 8816 (1993). [5] M. A. McMahon and K. B. Whaley, “Variational and diffusion Monte Carlo studies of (H2 )N clusters,” Chem. Phys. 182, 119 (1994). [6] K. B. Whaley, “Structure and dynamics of quantum clusters,” Int. Rev. Phys. Chem. 13, 41 (1994). [7] G. Tejeda, J. Fern´ andez, S. Montero, D. Blume, and J. Toennies, “Raman spectroscopy of small para-H2 clusters formed in cryogenic free jets,” Phys. Rev. Lett. 92, 223401 (2004). [8] E. Rabani and J. Jortner, “Spatial delocalization in paraH2 clusters,” J. Phys. Chem. B 110, 18893 (2006). [9] F. Mezzacapo and M. Boninsegni, “Superfluidity and quantum melting of p-H2 clusters,” Phys. Rev. Lett. 97, 045301 (2006). [10] S. Khairallah, M. Sevryuk, D. Ceperley, and J. Toennies, “Interplay between magic number stabilities and superfluidity of small parahydrogen clusters,” Phys. Rev. Lett. 98, 183401 (2007). [11] F. Mezzacapo and M. Boninsegni, “Structure, superfluidity, and quantum melting of hydrogen clusters,” Phys. Rev. A 75, 033201 (2007). [12] F. Mezzacapo and M. Boninsegni, “Superfluidity of isotopically doped parahydrogen clusters,” Phys. Rev. A 76, 021201 (2007). [13] F. Mezzacapo and M. Boninsegni, “Local superfluidity of parahydrogen clusters,” Phys. Rev. Lett. 100, 145301 (2008). [14] M. Boninsegni, “Quantum statistics and the momentum distribution of liquid parahydrogen,” Phys. Rev. B 79, 174203 (2009). [15] F. Mezzacapo and M. Boninsegni, “Classical and quantum physics of hydrogen clusters,” J. Phys.: Condensed Matter 21, 164205 (2009). [16] F. Mezzacapo and M. Boninsegni, “Parahydrogen clusters: Numerical estimates and physical effects,” J. Phys.: Conference Series 150, 032059 (2009). [17] F. Mezzacapo and M. Boninsegni, “On the possible supersolid character of parahydrogen clusters,” J. Phys. Chem. A 115, 6831 (2011). [18] S. Warnecke, M. Sevryuk, D. Ceperley, J. Toennies, R. Guardiola, and J. Navarro, “The structure of para-

hydrogen clusters,” Eur. Phys. J. D 56, 353 (2010). [19] O. Osychenko, R. Rota, and J. Boronat, “Superfluidity of metastable glassy bulk para-hydrogen at low temperature,” Phys. Rev. B 85, 224513 (2012). [20] J. Navarro and R. Guardiola, “Thermal effects on small para-hydrogen clusters,” Int. J. Quantum Chem. 111, 463 (2011). [21] T. Zeng and P.-N. Roy, “Microscopic molecular superfluid response: theory and simulations,” Rep. Prog. Phys. 77, 046601 (2014). [22] C. Chakravarty, “Fourier path integral simulations of para-H2 and ortho-D2 clusters,” Mol. Phys. 84, 845 (1995). [23] C. Chakravarty, “Structure of binary quantum clusters,” Phys. Rev. Lett. 75, 1727 (1995). [24] C. Chakravarty, “Isothermal-isobaric ensemble simulations of melting in quantum solids,” Phys. Rev. B 59, 3590 (1999). [25] J. Deckman and V. A. Mandelshtam, “The ground state estimation by global optimization of an effective potential. Application to binary para-H2 /ortho-D2 molecular clusters,” J. Phys. Chem. A 114, 9820 (2010). [26] E. Sola and J. Boronat, “Solidification of small p-H2 clusters at zero temperature,” J. Phys. Chem. A 115, 7071 (2011). [27] J. E. Cuervo and P.-N. Roy, “Path integral ground state study of finite-size systems: Application to small (parahydrogen)N (N = 2 − 20) clusters,” J. Chem. Phys. 125, 124314 (2006). [28] J. E. Cuervo and P.-N. Roy, “On the solid- and liquidlike nature of quantum clusters in their ground state,” J. Chem. Phys. 128, 224509 (2008). [29] M. Schmidt, S. Constable, C. Ing, and P.-N. Roy, “Inclusion of trial functions in the Langevin equation path integral ground state method: Application to parahydrogen clusters and their isotopologues,” J. Chem. Phys. 140, 234101 (2014). [30] C. Chakravarty, “Melting of neon clusters: Path integral Monte Carlo simulations,” J. Chem. Phys. 102, 956 (1995). [31] F. Calvo, J. Doye, and D. Wales, “Quantum partition functions from classical distributions: Application to rare-gas clusters,” J. Chem. Phys. 114, 7312 (2001). [32] J. P. Doye and F. Calvo, “Entropic effects on the structure of Lennard-Jones clusters,” J. Chem. Phys. 116, 8307 (2002). [33] C. Predescu, P. A. Frantsuzov, and V. A. Mandelshtam, “Thermodynamics and equilibrium structure of Ne38 cluster: Quantum mechanics versus classical,” J. Chem. Phys. 122, 154305 (2005). [34] P. A. Frantsuzov, D. Meluzzi, and V. A. Mandelshtam, “Structural transformations and melting in neon clusters: Quantum versus classical mechanics,” Phys. Rev. Lett.



Page 8 of 10 8

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

96, 113401 (2006). [35] J. Deckman, P. A. Frantsuzov, and V. A. Mandelshtam, “Quantum transitions in Lennard-Jones clusters,” Phys. Rev. E 77, 052102 (2008). [36] J. Deckman and V. A. Mandelshtam, “Effects of quantum delocalization on structural changes in Lennard-Jones clusters,” J. Phys. Chem. A 113, 7394 (2009). [37] J. E. Cuervo and P.-N. Roy, “Weakly bound complexes trapped in quantum matrices: Structure, energetics, and isomer coexistence in (para-H2 )N (ortho-D2 )3 clusters,” J. Chem. Phys. 131, 114302 (2009). [38] N. Faruk, M. Schmidt, H. Li, R. J. Le Roy, and P.N. Roy, “First-principles prediction of the Raman shifts in parahydrogen clusters,” J. Chem. Phys. 141, 014310 (2014). [39] M. Schmidt, J. M. Fernandez, N. Faruk, M. Nooijen, R. J. Le Roy, J. H. Morilla, G. Tejeda, S. Montero, and P.N. Roy, “Raman vibrational shifts of small clusters of hydrogen isotopologues,” J. Phys. Chem. A 119, 12551– 12561 (2015). [40] R. Guardiola and J. Navarro, “Structure of small clusters of parahydrogen molecules,” Phys. Rev. A 74, 025201 (2006). [41] R. Guardiola and J. Navarro, “A diffusion Monte Carlo study of small para-hydrogen clusters,” Cent. Eur. J. Phys. 6, 33 (2008). [42] J. Navarro and R. Guardiola, “Energetics and structure of small para-hydrogen clusters,” J. Low Temp. Phys. 148, 857 (2007). [43] F. Operetto and F. Pederiva, “Diffusion Monte Carlo study of the equation of state of solid para-H2 ,” Phys. Rev. B 73, 184124 (2006). [44] C. Chakravarty, “Quantum delocalization and cluster melting,” J. Chem. Phys. 103, 10663 (1995). [45] C. Chakravarty and R. Ramaswamy, “Instantaneous normal mode spectra of quantum clusters,” J. Chem. Phys. 106, 5564 (1997). [46] C. Chakravarty, “Path integral simulations of quantum Lennard-Jones solids,” J. Chem. Phys. 116, 8938 (2002). [47] C. Chakravarty, “Energy landscapes of quantum Lennard-Jones solids,” J. Phys. Chem. A 115, 7028 (2011). [48] M. B. Sevryuk, J. P. Toennies, and D. M. Ceperley, “Why are para-hydrogen clusters superfluid? A quantum theorem of corresponding states study,” J. Chem. Phys. 133, 064505 (2010). [49] I. F. Silvera and V. V. Goldman, “The isotropic intermolecular potential for H2 and D2 in the solid and gas phases,” J. Chem. Phys. 69, 4209 (1978). [50] U. Buck, F. Huisken, A. Kohlhase, D. Otten, and J. Schaefer, “State resolved rotational excitation in D2 +H2 collisions,” J. Chem. Phys. 78, 4439 (1983). [51] M. Boninsegni and S. Moroni, “Population size bias in diffusion Monte Carlo,” Phys. Rev. E 86, 056712 (2012). [52] J. B. Anderson, “A random-walk simulation of the Schr¨ odinger equation: H+ 3 ,” J. Chem. Phys. 63, 1499 (1975). [53] J. B. Anderson, “Quantum chemistry by random walk. 1 ′ 3 + 1 + 1 H 2 P, H+ 3 D3h A1 , H2 Σu , H4 Σg , Be S,” J. Chem. Phys. 65, 4121 (1976). [54] A. B. McCoy, “Diffusion Monte Carlo approaches for investigating the structure and vibrational spectra of fluxional systems,” Int. Rev. Phys. Chem. 25, 77 (2006).

[55] L. M. Johnson and A. B. McCoy, “Evolution of structure in CH+ 5 and its deuterated analogues,” J. Phys. Chem. A 110, 8213 (2006). [56] P. H. Acioli, Z. Xie, B. J. Braams, and J. M. Bowman, “Vibrational ground state properties of H+ 5 and its isotopomers from diffusion Monte Carlo calculations,” J. Chem. Phys. 128, 104318 (2008). [57] P. Barrag´ an, R. P´erez de Tudela, C. Qu, R. Prosmiti, and J. M. Bowman, “Full-dimensional quantum calculations of the dissociation energy, zero-point, and 10 K + properties of H+ 7 /D7 clusters using an ab initio potential energy surface,” J. Chem. Phys. 139, 024308 (2013). [58] J. S. Mancini, A. K. Samanta, J. M. Bowman, and H. Reisler, “Experiment and theory elucidate the multichannel predissociation dynamics of the HCl trimer: Breaking up is hard to do,” J. Phys. Chem. A 118, 8402 (2014). [59] A. Shank, Y. Wang, A. Kaledin, B. J. Braams, and J. M. Bowman, “Accurate ab initio and ”hybrid” potential energy surfaces, intramolecular vibrational energies, and classical ir spectrum of the water dimer,” J. Chem. Phys. 130, 144314 (2009). [60] G. Czak´ o, Y. Wang, and J. M. Bowman, “Communication: Quasiclassical trajectory calculations of correlated product-state distributions for the dissociation of (H2 O)2 and (D2 O)2 ,” J. Chem. Phys. 135, 151102 (2011). [61] Y. Wang and J. M. Bowman, “Communication: Rigorous calculation of dissociation energies (D0 ) of the water trimer, (H2 O)3 and (D2 O)3 ,” J. Chem. Phys. 135, 131101 (2011). [62] Y. Wang, V. Babin, J. M. Bowman, and F. Paesani, “The water hexamer: Cage, prism, or both. Full dimensional quantum simulations say both,” J. Am. Chem. Soc. 134, 11116 (2012). [63] J. D. Mallory, S. E. Brown, and V. A. Mandelshtam, “Assessing the performance of the diffusion Monte Carlo method as applied to the water monomer, dimer, and hexamer,” J. Phys. Chem. A 119, 6504 (2015). [64] J. D. Mallory and V. A. Mandelshtam, “Binding energies from diffusion Monte Carlo for the MB-pol H2 O and D2 O dimer: A comparison to experimental values,” J. Chem. Phys. 143, 144303 (2015). [65] J. D. Mallory and V. A. Mandelshtam, “Diffusion Monte Carlo studies of MB-pol (H2 O)2−6 and (D2 O)2−6 clusters: Structures and binding energies,” J. Chem. Phys. 145, 064308 (2016). [66] R. Assaraf and M. Caffarel, “Zero-variance zero-bias principle for observables in quantum Monte Carlo: Application to forces,” J. Chem. Phys. 119, 10536 (2003). [67] R. J. Hinde, “A six-dimensional H2 –H2 potential energy surface for bound state spectroscopy,” J. Chem. Phys. 128, 154308 (2008). [68] M. A. Suhm and R. O. Watts, “Quantum Monte Carlo studies of vibrational states in molecules and clusters,” Phys. Rep. 204, 293 (1991). [69] G. L. Warren and R. J. Hinde, “Population size bias in descendant-weighted diffusion quantum Monte Carlo simulations,” Phys. Rev. E 73, 056706 (2006). [70] F. H. Stillinger and T. A. Weber, “Inherent structure in water,” J. Phys. Chem. 87, 2833 (1983). [71] J. P. Doye, D. J. Wales, and M. A. Miller, “Thermodynamics and the global optimization of Lennard-Jones clusters,” J. Chem. Phys. 109, 8143 (1998).


Page 9 of 10


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

[72] J. P. Doye and D. J. Wales, “Thermodynamics of global optimization,” Phys. Rev. Lett. 80, 1357 (1998). [73] J. P. Doye, M. A. Miller, and D. J. Wales, “The doublefunnel energy landscape of the 38-atom Lennard-Jones cluster,” J. Chem. Phys. 110, 6896 (1999). [74] J. P. Doye, M. A. Miller, and D. J. Wales, “Evolution of the potential energy surface with size for Lennard-Jones clusters,” J. Chem. Phys. 111, 8417 (1999). [75] F. Calvo, J. Neirotti, D. L. Freeman, and J. Doll, “Phase changes in 38-atom Lennard-Jones clusters. II. A parallel tempering study of equilibrium and dynamic properties in the molecular dynamics and microcanonical ensembles,” J. Chem. Phys. 112, 10350 (2000). [76] J. Neirotti, F. Calvo, D. L. Freeman, and J. Doll, “Phase changes in 38-atom Lennard-Jones clusters. I. A parallel tempering study in the canonical ensemble,” J. Chem. Phys. 112, 10340 (2000). [77] P. A. Frantsuzov and V. A. Mandelshtam, “Sizetemperature phase diagram for small Lennard-Jones clus-

ters,” Phys. Rev. E 72, 037102 (2005). [78] V. A. Mandelshtam and P. A. Frantsuzov, “Multiple structural transformations in Lennard-Jones clusters: Generic versus size-specific behavior,” J. Chem. Phys. 124, 204511 (2006). [79] L. Zhan, J. Z. Chen, and W.-K. Liu, “Determination of structural transitions of atomic clusters from local and global bond orientational order parameters,” J. Chem. Phys. 127, 141101 (2007). [80] R. M. Sehgal, D. Maroudas, and D. M. Ford, “Phase behavior of the 38-atom Lennard-Jones cluster,” J. Chem. Phys. 140, 104312 (2014). [81] H. M. Cezar, G. G. Rondina, and J. L. Da Silva, “Parallel tempering Monte Carlo combined with clustering Euclidean metric analysis to study the thermodynamic stability of Lennard-Jones nanoclusters,” J. Chem. Phys. 146, 064114 (2017). [82] D. J. Wales, “Decoding heat capacity features from the energy landscape,” Phys. Rev. E 95, 030105(R) (2017).




Page 10 of 10

Quantum Melting and Isotope Effects from Diffusion Monte Carlo

Recommend Documents