Coating Polycyclic Aromatic Hydrocarbon Cations with Helium Clusters

Dec 8, 2014 - The classical and quantum structures of cationic polycyclic aromatic hydrocarbon molecules (benzene, pyrene, coronene, and circumcoronen...
0 downloads 12 Views 4MB Size
Article pubs.acs.org/JPCA

Coating Polycyclic Aromatic Hydrocarbon Cations with Helium Clusters: Snowballs and Slush Florent Calvo* University of Grenoble Alpes, LIPHY, F-38000 Grenoble, France and CNRS, LIPHY, F-38000 Grenoble, France S Supporting Information *

ABSTRACT: The classical and quantum structures of cationic polycyclic aromatic hydrocarbon molecules (benzene, pyrene, coronene, and circumcoronene) coated by helium atoms have been theoretically investigated using a variety of computational methods. Classical shell filling, as determined from global optimization, is generally found to proceed by epitaxial additions on the graphitic surfaces before peripheral closure. From the quantum mechanical perspective provided by path-integral molecular dynamics simulations, vibrational delocalization is found to generally decrease the size of this first solvation shell, but also to give rise to a variety of situations in which the helium atoms are more or less localized depending on their environment, with strong finite size effects depending both on the hydrocarbon cation and the number of coating helium atoms. While the graphitic planes tend to bind helium sufficiently to give rise to snowball precursors, the peripheral regions are less dense and more delocalized, not as liquid as the outer layers but within an intermediate slushy character. first shell. One particularly suitable experimental method to determine the size of the snowball is undoubtly mass spectrometry, in which the relative stabilities of neighboring cluster sizes can be inferred from the abundance distributions, anomalies obtained under nonequilibrium conditions being often related to special behavior such as shell or subshell closing, as particularly relevant in the case of snowballs. This approach has been pursued for rare gas impurities12−15 and, more recently, for a great variety of atomic16 and molecular17,18 ions by the Scheier group. The interpretation of experimental measurements requires dedicated computations that account, as much as possible, for the highly quantum nature of helium and, as accurately as possible, for the detailed interactions with the impurity. Alkali cations have been particularly scrutinized,18−20 and revealed significant variations in the size of the snowball depending on the cation. Calculations26−29 have notably found an icosahedral shell of 12 atoms around Na+, similarly to earlier findings on embedded argon and kryption cations,12−15 also found for lead.30 However, these results have been recently questioned by experiments which reevaluated the snowball size to be closer to 9 atoms instead,18 suggesting subtle effects that were missing in the previous quantum Monte Carlo simulations.31 Besides alkali ions and dimers,18−20 snowball sizes have been experimentally

1. INTRODUCTION Helium droplets are chemically inert solvents that can reach extremely low temperatures enabling high resolution spectroscopy.1 The fascinating properties of helium droplets also make them a natural cryogenic laboratory for studying chemical reactions in cold environments such as the interstellar medium,2,3 extending the possibilities offered by heavier rare gases in which the reaction takes place at the surface.4 One particular feature of helium compounds is the emergence of superfluidity at sufficiently low temperatures. Superfluidity occurs both for 4He and 3He but at different temperatures, and in small droplets the λ transition was predicted to depend on size,5 as commonly observed in conventional finite size phase transitions.6 The superfluid phase in helium droplets has been evidenced indirectly from spectroscopic7−10 or transport11 measurements. Both types of experiments rely on atomic or molecular dopants or impurities such as electron bubbles, and probe the different response of those impurities under isolated or helium environments, possibly of either isotope. Ionic impurities can have a marked influence on the behavior of helium droplets, due to the stronger binding they exert in their vicinity owing to polarization forces. The ability of ions to drag tens of helium atoms with them has been demonstrated in early experiments on bulk helium,1 and the phenomenon responsible for this localization has since been denoted as electrostriction. One manifestation of electrostriction in helium droplets doped with ionic impurities is the formation of socalled snowballs, referring to a first solvation shell undergoing stronger binding, greater confinement associated with the loss of superfluidity, fluidlike behavior being recovered above the © XXXX American Chemical Society

Special Issue: Jean-Michel Mestdagh Festschrift Received: October 28, 2014 Revised: December 6, 2014

A

dx.doi.org/10.1021/jp510799h | J. Phys. Chem. A XXXX, XXX, XXX−XXX

The Journal of Physical Chemistry A

Article

determined for protonated hydrogen,21 halogens,22 and larger molecular ions such as fullerenes23 and even fullerene clusters.24 Obviously, larger molecules support larger solvation shells and this has been shown in calculations of fullerenes such as C20 (ref 25), or in C+60 where the interpretation of the measured snowball size of 60 helium atoms appeared nontrivial,23 calculations32,33 indicating strong size effects, and even reentrant freezing similar to the phenomenology observed for bulk helium coating graphitic substrates.34 One particular case of molecular dopant is aromatic hydrocarbons, which have been heavily studied in the limit of low35−42 and high11,43−46 helium coverage by means of electronic-vibrational spectroscopy. Detailed experiments35,38,40−42 and calculations40,47−51 concurred to show that the interaction between the helium atoms and the molecular substrate is highly anisotropic, especially in electronic excited states, giving rise to equally anisotropic delocalization and confinement manifested by specific spectral features such as so-called zero-phonon lines and zero-phonon wings. Polycylic aromatic hydrocarbons (PAHs) have recently gained a lot of attention from the physical chemistry community because of their relevance as combustion products52 and in astrochemistry as possible carriers of specific infrared bands in the emission spectra of interstellar clouds.53 Such rather large molecular compounds could be synthesized upon electron irradiation of smaller aromatics in helium droplets.54 The structural similarity between large PAHs and graphene would suggest similar helium adsorption properties; however, quantum Monte Carlo simulations indicate that in extended monolayers the delocalized fluid phase competes with the localized √3 × √3 phase.55,56 As revealed by similar methods, finite 4He clusters could be adsorbed on both sides of graphene and also exhibit localized structures commensurate with the √3 × √3 lattice, provided that the substrate corrugation is properly accounted for.57 The finite extent of the graphitic clusters, the presence of peripheral hydrogens, and the ionic character required by mass spectrometry experiments are additional factors that make helium solvation around larger polyaromatics a more complicated problem than that around fullerenes. The present work aims to address this issue by means of dedicated molecular modeling and to evaluate the extent of delocalization in helium clusters coating three archetypal cationic PAHs, namely pyrene (C16H+10), coronene (C24H+12), + ), supplemented by benzene and circumcoronene (C54H18 (C6H+6 ). Contrary to the vast majority of earlier theoretical works, we focus here on size effects rather than investigating in details a few selected sizes. Our objectives in this work are 4fold: (i) from a qualitative point of view we wish to explore classically and quantum mechanically the stable solvation structures accross the first shell; (ii) more quantitatively, we wish to evaluate the onset of snowballing, or the size above which additional atoms form extra liquid layers; (iii) analyzing our atomistic simulations in terms of specific indices characterizing spatial delocalization and relative rigidity should provide a more detailed picture about the collective behavior of the helium solvent around the highly anisotropic PAH molecules; and (iv) finally we propose to correlate those indices to more global probes of the system’s rigidity based on an energy landscape perspective considering the equilibrium distributions of inherent structures.58 The computational framework chosen to address the above issues is that of path-integral molecular dynamics (PIMD) for

spinless 4He particles, that is, neglecting bosonic exchange statistics. This approach is thus limited to the normal fluid state (T > 0.4 K); however, it is practical enough to allow clusters containing hundreds of helium atoms to be simulated efficiently. The results obtained using the PIMD approach confirm earlier conclusions47−51 about the anisotropic character of delocalization, atoms adsorbed over the aromatic center being more confined than peripheral atoms. Upon increasing the number of solvent atoms, both faces of the aromatic hydrocarbon are coated before the outer regions close to the molecular plane are also filled. Eventually a highly anisotropic snowball is formed with liquid layers on both sides of the PAH, the belt region in the neighborhood of hydrogens showing intermediate delocalization akin to slush rather than a fully developed liquid. In larger droplets, confinement is again favored in the immediate vicinity of aromatic units, the distinction between the slush and outer liquid phases becoming less clear. The next section describes the computational tools implemented to tackle the problem of PAH coating by helium clusters, including the potential energy surfaces, structural and quantum dynamical simulation methods, as well as the various indicators that were used to characterize the different phases and extent of confinement. The results are presented and discussed in section 3, emphasizing the strongly nonmonotonic effects found as the size of the hydrocarbon or the number of helium atoms are varied independently from one another, but also some common trends such as snowballing or the presence of peripheral slush phases. Concluding remarks finally close the paper in section 4.

2. METHODS The present computational investigation relies on a combination of atomistic tools that were successfully used in a previous study of C+60Hen clusters.32 We briefly describe those tools and highlight the novelty and additional features with respect to this earlier effort. In the following, the cationic aromatic dopant CmH+p is taken at the minimum geometry as predicted by standard density-functional theory calculations employing the B3LYP functional and the 6-311++G** basis set, and the atomic charges that reproduce best the electrostatic potential around the molecule are evaluated.59 Those charges and the detailed geometries are all documented as Supporting Information. It is important to notice here that the cationic nature of the dopant has at least two consequences on its interactions with helium atoms. First, it enhances the overall binding due to polarization forces. In addition, it slightly alters the molecular geometry due to Jahn−Teller distortion. In practice though, these latter effects were found to be rather minor especially in polycyclic aromatic molecules where the graphitic structure is essentially preserved. The optimized geometry was then placed in the z = 0 plane and assumed to be fixed there and rigid, the molecule being frozen in its ground vibrational state. The problem we wish to solve is therefore the optimal location of helium atoms in a finite number, in the classical and quantum mechanical descriptions of nuclear motion. We use the simplified notation of PAH+Hen to denote the system of the polyaromatic hydrocarbon cation coated with n 4He atoms. 2.1. Potential Energy Surfaces. We denote by R = {ri} the set of Cartesian coordinates of the n helium atoms, and by R′ = {ri′} the coordinates of the m + p carbon and hydrogen atoms of the dopant. Because of the relatively heavy B

dx.doi.org/10.1021/jp510799h | J. Phys. Chem. A XXXX, XXX, XXX−XXX

The Journal of Physical Chemistry A

Article

structures for the subsequent molecular dynamics simulations rather than trying to build realistic configurations by hand. We have used the standard optimization method of basinhopping64 to determine the putative global minima of PAH+Hen for n ranging from 1 to 60 (benzene), 90 (pyrene), 120 (coronene), and 160 (circumcoronene), all sizes being considered in these ranges. For each system, five series of 104 local minimizations were performed and the Monte Carlo displacements were accepted based on the Metropolis probability at T = 5 K. 2.3. Path-Integral Molecular Dynamics. To capture the essential effects of zero-point delocalization, we have carried out path-integral molecular dynamics simulations of the nuclear wave function of selected PAH+Hen clusters, using a Trotter number of P = 64 at the temperature T = 1 K where exchange effects can be ignored for 4He. It is not the purpose of the present article to describe in length the PIMD method, recent reviews on the subject being available.65,66 Briefly, the PIMD equations of motion are propagated in normal mode coordinates and thermostatting is applied independently to each of the P replicas using the Nosé−Hoover method with the recommended thermostat mass Q = 3n/P/kBT, where kB is the Boltzmann constant.65 The contraction scheme of Markland and Manolopoulos67 was applied to the polarization energy, this contribution being evaluated at the centroids only. The PIMD simulations used a time step of 1 fs, and for each system studied averages and properties were accumulated from 1 ns long trajectories after a period of 200 ps discarded for equilibration. Finite temperature methods based on path integrals are more straightforwardly implemented than variational and diffusion Monte Carlo simulations, which require carefully chosen trial wave functions in order to be efficient. They are less accurate in predicting ground state energies, and in the case of PIMD they do not account for bosonic statistics that could further contribute to delocalization.29 In carrying PIMD simulations at the moderately high temperature of 1 K, we somehow compensate for this underestimation of quantum mechanical delocalization by enhancing the thermal component. 2.4. Delocalization Indicators. Various properties were considered to analyze the structures and solvation patterns obtained from our simulations. The classical binding energies of the global minima naturally give insight into the dominant contributions to binding in the system. The harmonic contribution to the zero-point energy, as given by diagonalizing the dynamical matrix at the minimum, further provides some information about the role of nuclear delocalization, although the harmonic approximation is expected to fail once highly delocalized regions are formed.68,69 To evaluate the amount of spreading of the nuclear wave function, and neglecting the Jahn-Teller distortions of the cationic dopants, we initially considered 2D density plots in cylindrical coordinates (r, z), z being the axis perpendicular to the aromatic molecule and r being the distance to this axis. However, because of usual the lack of axial symmetry for the atoms in direct contact with the dopant, those plots turned out to be often not insightful. The 3D density plots were chosen instead. In addition we considered two separate indicators to measure the overall degrees of delocalization and rigidity in the system. The extent of localization of each atom i was simply measured from the fluctuations of the individual positions of the various replicas in the PIMD simulations,

computational cost of path-integral simulations, it is necessary to use a practical representation of the interactions between helium atoms and the frozen hydrocarbon molecule. The potential energy surface V(R) is thus mostly additive, with pairwise contributions between all mobile atoms, and completed with a many-body polarization contribution reflecting the spreading of the charge on the molecule. We thus write V as a sum of He−He and He−PAH contributions: V (R) = VHe(R) + VPAH(R, R′) VHe(R) =



VJA(rij)

i , j ∈ Hen; i < j

VPAH(R, R′) =

∑ ∑ i ∈ PAH j ∈ Hen

(i) V pair (rij) +



Vpol(ri, R′)

i ∈ Hen

In addition to the Janzen−Aziz potential VJA between helium atoms,60 the total potential contains helium−PAH interactions VPAH(R, R′) taken as pairwise terms V(i) pair between all atoms involved, and a sum of polarization contributions Vpol(ri, R′) acting on each helium atom i but originating from the partial charges distributed on the PAH atoms. The interactions between neutral atoms and molecules being treated as pairwise are assumed as transferable among the various PAHs. The pair interaction between helium and carbon atoms was borrowed from the work on helium coating fullerenes,32 in which the pair interaction was carefully parametrized against electronic structure calculations61 and will not be repeated here. For the helium−hydrogen interaction, and in the absence of reference data of comparable quality as the aforecited helium−carbon interaction, we used a simple Lennard-Jones form with parameters already used by Heidenreich and Jortner for anthracene48 and taken from an earlier work by Lim62 on the collisional energy transfer between helium and toluene. Finally, for the polarization contribution we only consider the dominant contribution from the electric fields Ei created at the sites of each helium atom as α Vpol(ri, R′) = − He Ei2 (1) 2 where αHe = 0.205 Å3 is the atomic polarizability of helium.63 The main approximation is thus the neglect of three-body forces between induced dipoles, but this would significantly increase the computational cost without any expected gain on the overall solvation processes, at least as long as quantitative details are not sought.31 Also, we note that a more refined electrostatic description with charges on bonds or higher-order multipoles was not deemed necessary owing to the relatively weak contribution of such corrections, especially in larger PAHs where the charges are significant only in peripheral regions. 2.2. Classical Approach: Global Optimization. While the most rigorous approach to study solvation of helium around molecular dopants is probably quantum (variational and diffusion) Monte Carlo, it is also computationally intensive and not practical for large dopants and a broad range of cluster sizes. Before even considering nuclear quantum effects, it is useful to tackle the solvation issue on a purely classical basis, by locating the lowest-energy structures of the coated compounds. This strategy not only provides a first estimate of the size of the first solvation shell, but also gives some insight into the growth mechanisms. Furthermore, it can be justified a posteriori if the helium atoms appear confined rather than delocalized. More pragmatically it is also best to evaluate suitable starting C

dx.doi.org/10.1021/jp510799h | J. Phys. Chem. A XXXX, XXX, XXX−XXX

The Journal of Physical Chemistry A σi =

1 P

Article

P

∑ [⟨(r(i k))2 ⟩ − ⟨r(i k)⟩2 ]1/2

(2)

k=1

r(k) i

where denotes the Cartesian positions of replica k of atom i. A global measure of delocalization ⟨σ⟩ can then be introduced as the average over σi. Complementary to these geometrical quantities we also considered global indicators measuring the amount of rigidity in the dynamics, as the root-mean-square bond length fluctuation index ⟨δ⟩ also known as the Lindemann index, evaluated in the PIMD simulations from the centroids as ⟨δ⟩ =

2 n(n − 1)

1 δi = n−1

∑ j≠i



⟨rij2⟩ − ⟨rij⟩2 ⟨rij⟩

i 0.2). However, it was found to be unsufficient here owing to many intermediate cases in which parts of the systems are rigid and in equilibrium with less ordered parts (vide infra). Similarly to the index σi we have thus also used the atom-resolved indices δi, from which ⟨δ⟩ is simply the average, to go beyond this globally averaged information. Besides σ and δ, the energy landscape perspective was also adopted in order to correlate the extent of delocalization or fluidlike character to specific underlying local minima. This approach has long been fruitful in the study of liquids58 and glasses71−73 and in the present case should provide a lot of insight into the extent of spreading of the nuclear wavefuction. From the PIMD sample, periodic quenches of the centroid configurations were performed and the corresponding local minima were recorded. The excitation energies and the distributions of these minima, also called inherent structures, were then obtained.

Figure 1. Classical (black symbols) and quantum-corrected (red symbols) binding energies of cationic polycyclic aromatic hydrocarbons as a function of the number N of helium atoms coating them. The arrows locate the sizes of the snowballs or first solvation shells: (a) benzene; (b) pyrene; (c) coronene; (d) circumcoronene.

structures, and this result also holds in the quantum mechanical case. Such distortions will therefore be ignored from now on by refering to the molecular symmetries as (approximately) sixfold. In general, small numbers of helium atoms tend to cover the graphitic surface initially near the center, with additional atoms alternating on the two faces of the PAH. At this stage the energy is thus minimized by placing the largest possible number of atoms close to the center of the PAH, maximizing the helium-PAH interaction at the expense of He−He interactions. This effect can even be significant if helium atoms are brought closer to peripheral regions. The coronene cation prefers for instance to have five helium atoms distributed as 3 + 2 on the two facets, a repartition that is slightly lower (by 1.6 cm−1) than the most favorable 4 + 1 repartition, and much lower (by 42.4 cm−1) than having all five atoms on one side. Epitaxial covering on the honeycomb lattice is favored especially for the 6-fold symmetric molecules, with a particularly stable 7-atom hexagonal motif on both faces reached at n = 14 and depicted in Figure 2 for benzene. The √3 × √3 commensurate solid found in quantum descriptions of helium monolayers on bulk graphitic substrates55,56,74 is also occasionally seen, as in the case of coronene+He11 also shown in Figure 2. However, the stabilization of this less dense phase is expected to occur precisely when quantum effects are included, and the classical structures usually exhibit the same honeycomb lattice as the dopant. Once the two faces of the hydrocarbon molecule are fully covered, additional atoms tend to reside closer to the molecular plane, in the peripheral region, making weak bridges between the two planar helium layers already formed. The distance between these layers (about 6.5 Å) is too large for He−He

3. RESULTS AND DISCUSSION The classical picture of the global minima is insightful as it provides structural information about the growth mechanisms of the helium droplet over the hydrocarbon cations. The binding energies of the PAH+Hen clusters have been calculated in this classical approximation for the four dopant molecules and all sizes ranging to the aforementioned maximum values, and in the global minimum configuration of helium atoms. Their variations with n are represented in Figure 1. For the four dopants the binding energies show comparable behaviors with a stronger slope at small sizes and more gentle variations in larger clusters, the transition between the two regimes being less and less marked as the hydrocarbon molecule enlarges. These two regimes are associated with stronger binding at low sizes concomitant with direct helium−PAH interactions, and weaker helium−helium contributions once a first complete layer has been formed around the dopant. Classical snowball sizes can be evaluated already at this stage simply by looking at the global minima, for which a selection is represented in Figure 2. These minima also illustrate the classical growth mechanisms for the four hydrocarbons. As a first comment, it is important to notice that the small Jahn-Teller distortion in the cationic molecules do not have any detectable consequence on the stable D

dx.doi.org/10.1021/jp510799h | J. Phys. Chem. A XXXX, XXX, XXX−XXX

The Journal of Physical Chemistry A

Article

doing so we neglect here the possible structural transitions associated with a high-lying isomer becoming the most stable structure upon including the zero-point energy contribution.75 The harmonic correction significantly shifts the binding energy for all dopants, as shown in Figure 1. In larger clusters, it even leads to a decrease (in magnitude) of the binding energy with increasing number of atoms, an unphysical result that conveys the approximate character of the harmonic correction to zeropoint energies in the regime dominated by helium−helium interactions.68,69 In particular, the size at which the binding energy starts decreasing is comparable to the classical snowball size, providing another empirical way to define nC. Interestingly, we notice that the variations of the binding energy harmonically corrected for the zero-point contribution are somewhat irregular, especially in the cases of benzene and pyrene although this is also true but less immediately visible for the other dopants. These complex variations are probably the result of the aforementioned ZPE-induced structural transitions, and would likely be washed out had the zero-point energies been included for all minima in our global optimization searches. The path-integral molecular dynamics approach is much more rigorous in capturing nuclear quantum effects, and it was applied to a large number of PAH+Hen clusters across the ranges in which the global minima were determined. For all hydrocarbons, all sizes n = 1−20 were studied. All other sizes n = 21−60 were studied for the benzene cation, while for larger dopants only even sizes were addressed in the range n = 20−30 and sizes of multiples of 5 above this limit, plus some extra sizes. Although the PIMD approach is suitable to unravel the general thermodynamical behavior and efficient enough to study large clusters, we have not used it to determine accurate binding energies, as the method is not as competitive as diffusion Monte Carlo. We show in Figure 3 and Figure 4 the variations with increasing numbers of helium atoms of the globally averaged Lindemann index ⟨δ⟩ and the geometrical radius ⟨σ⟩ measuring the spreading of the nuclear wave function, respectively. The fluctuations of those averaged properties, which measure the heterogeneity in the dynamics and delocalization among individual atoms, are depicted by error bars. Both properties turn out to be generally correlated to one another, and exhibit strong nonmonotonic size effects for all systems, especially at small sizes n < 30. The smallest systems, for n = 1 and 2, are essentially rigid, but already for three atoms or more significant degrees of fluxionality are found. In the large cluster limit, that is above the classical snowball size nC, also a common behavior is found with a steady increase in both indicators that conveys the increasing proportion of liquidlike atoms in such clusters. The systematically high fluctuations of both δ and σ in this large size limit indicate that the dynamics and localization are strongly heterogeneous, consistent with the picture of a snowball with rigidlike regions close to the dopant in equilibrium with liquidlike outer regions. But the rich variations below this upper limit remain to be interpreted. We first notice that these variations are qualitatively similar to the case of helium adsorbed on C+60, with small clusters exhibiting localization, intermediate clusters being fluidlike due to structural disorder, and reentrant freezing at n = 60 corresponding to the snowball, eventually covered by additional liquid layers above this size.23,32,33 The phenomenology of helium clusters around PAH cations thus bears resemblance to the fullerenes, but with notable differences and a marked dependence on the hydrocarbon dopant.

Figure 2. Global minima of selected cationic aromatic hydrocarbons coated with helium atoms. For each hydrocarbon, three sizes are shown below, at, and above the classical snowball size.

bonds, and those additional atoms are located off the plane by about 1.2 Å, alternating laterally in a zigzag fashion. For all dopants, there is a maximum size nC of helium atoms above which additional atoms can no longer be accommodated in the belt region, and a second layer is nucleated instead on top of one face. The classical snowball size nC is evaluated from our calculations to be 34, 56, 70, and 112 for the four increasingly large dopants. Those values have been highlighted in Figure 1. They are associated with the transition between the two binding regimes, an effect that is particularly clear in the cases of benzene and pyrene with a marked change in slope. When nC is represented against the number of aromatic cycles, the number of total or carbon atoms does not reveal any immediate correlation, which suggests that circumcoronene is still relatively small in the scaling regime of finite size effects. As further atoms are added to the clusters at size nC, both layers on each side of the dopant are themselves covered rather homogeneously, with comparable numbers of atoms over each face. Examples of global minima having a partially filled second layer are depicted in Figure 2 for the largest sizes considered in the global optimization study, namely 60, 90, 118, and 160. The classical treatment for nuclei is questionable for helium, except perhaps at small sizes where the atoms are most strongly bound to the dopant. We first discuss the role of nuclear quantum effects on the binding energies, applying here a simple harmonic approximation to the zero-point energy from the normal-mode frequencies obtained after diagonalizing the dynamical matrix at the global minimum configuration. In E

dx.doi.org/10.1021/jp510799h | J. Phys. Chem. A XXXX, XXX, XXX−XXX

The Journal of Physical Chemistry A

Article

Figure 3. Average Lindemann fluctuation index ⟨δ⟩ of helium clusters coating cationic aromatic hydrocarbons as a function of their number N, as obtained from path-integral molecular dynamics simulations at 1 K. The error bars denote the fluctuations among the individual contributions to the average: (a) benzene; (b) pyrene; (c) coronene; (d) circumcoronene.

Figure 4. Average spatial extension ⟨σ⟩ of helium clusters coating cationic aromatic hydrocarbons as a function of their number N, as obtained from path-integral molecular dynamics simulations at 1 K. The error bars denote the fluctuations among the individual contributions to the average: (a) benzene; (b) pyrene; (c) coronene; (d) circumcoronene.

At this stage the global indicators only provide definitive clues about the localized state of the helium clusters in the lowand high-value limits. When ⟨δ⟩ < 0.1, which approximately corresponds to ⟨σ⟩ < 1.5 Å, the entire system is rigidlike and localized. Examples include benzene+He29, pyrene+He25, coronene+He50, and circumcoronene+He38. Conversely, if ⟨δ⟩ > 0.15 (and ⟨σ⟩ > 2.5 Å) then the entire cluster is fluidlike and delocalized. This occurs preferentially for benzene and circumcoronene, for example, for benzene+He16 and circumcoronene+He18, but we cannot rule out that the other two dopants may also exhibit this stronger delocalization for sizes that we have not simulated by PIMD. In between those limits the clusters do not display a fully uniform collective behavior, as shown by particularly large fluctuations around the average values. Atom-resolved indicators are necessary in order to better characterize their dynamics and degree of delocalization or confinement and the heterogeneity of this behavior. Before discussing the phenomenology associated with those strong finite size effects, it is instructive to relate the previous findings to the inherent structures sampled by the PIMD centroids wave function. For each simulated cluster, a global indicator was again constructed from the set {α} of local minima with energies {Eα} and their occurences {gα} in the sample. The statistical fluctuations δE of this distribution were calculated and used as a measure of the spreading of the wave function in the energy landscape. δE is also a measure of nonrigidity, as a strictly positive value is necessarily associated with multiple inherent structures indicating a behavior that is at

least partly liquidlike. However, δE is unable to distinguish between permutational isomers, and delocalization involving such minima may not be visible even though they have clearer signatures on ⟨δ⟩ > 0.1. The variations of δE with increasing numbers of coating helium atoms are represented in Figure 5 for the four cationic dopants. The energy indicator also displays nonmonotonic size effects for all hydrocarbons, with a notable maximum near the snowball size and a smooth increase above this limit of a few wavenumbers with each added helium atoms. Finite excitations may occur even at small sizes, especially for the largest dopants, and the presence of several minima is interpreted as the consequence of some mobility of the small helium clusters on each side of the hydrocarbon. It is worth noticing that a mild correlation between the geometrical and dynamical indicators is found with δE, at least in the globally localized cases. This is best seen on the rigidlike structures of benzene+He22, coronene+He50, or circumcoronene+He70, for instance. To some extent the same observation is made for the globally delocalized cases with the examples of benzene+He24, pyrene+He8, or circumcoronene+He50. However, there are cases such as benzene+He32 where δE is high with no sign of fluidlike character. Those unexpected situations actually reflect some structural transitions as mentioned above in relation with zero-point energy corrections, and will be further discussed in what follows. F

dx.doi.org/10.1021/jp510799h | J. Phys. Chem. A XXXX, XXX, XXX−XXX

The Journal of Physical Chemistry A

Article

Figure 6. Atom-resolved dynamical properties of cationic benzene coated with selected numbers N of helium atoms: (upper panel) spatial extension versus Lindemann index; (middle panel) distribution of inherent structure energies; (lower panel) density plots. (a) N = 6; (b) N = 24; (c) N = 32. Figure 5. Fluctuations of the energy minima found by locally minimizing the PIMD sample of configurations, as a function of the number of helium atoms coating the cationic aromatic hydrocarbons: (a) benzene; (b) pyrene; (c) coronene; (d) circumcoronene.

vacancy. There are three inherent structures differing in the location of these five helium atoms, but the global dynamics is heavily liquidlike, only the single atom on the opposite side remaining in a confined vibrating state. The cluster at n = 24 provides an example of a globally disordered system, the liquidlike character being again caused by the presence of several vacancies in the belt peripheral region close to the dopant plane. For this system all atoms are able to interchange between the layers on either side through this peripheral region, resulting in the nearly uniform liquid layer seen on the density plot. At n = 32 all atoms appear well localized but the cluster occupies mostly a minimum lying 17 cm−1 above the classical global minimum. The spontaneous evolution of the cluster toward this mestastable minimum is related to its greater stability in terms of zero-point contribution. At the harmonic level, the metastable isomer becomes lower than the global minimum by more than 53 cm−1, which confirms the results of the PIMD trajectories. This system illustrates a spontaneous structural transition induced by zero-point effects. Our discussion of helium clusters coating hydrocarbon cations can now shift to the first polycyclic dopant, pyrene. The three sizes chosen for illustration are n = 8, 35, and 80, and the same properties previously considered for benzene are represented in Figure 7. The smallest cluster provides an example of localized coating on the √3 × √3 commensurate lattice. Although the atoms are confined and only visit one inherent structure, their wave function spreads over nearly 1.8 Å, significantly more than the isolated localized atom in benzene+He6. This increased spreading is consistent with the location of the helium atoms at the edge of the dopant. At n = 35 the starting structure has complete helium layers of 14 atoms each on each side of the dopant; the 7 extra atoms lie in the peripheral region. The PIMD trajectories predict a rather mixed phase in which these atoms are distinctively delocalized

In larger clusters above the classical snowball size, the increasingly fluidlike character is associated with the formation of a liquid layer, whose contribution to ⟨δ⟩ and ⟨σ⟩ is in proportion to the number of actual delocalized atoms. Figure 5 shows that those structures are also associated with increasingly high-lying isomers. The slope of the binding energy in Figure 1 in this regime is significantly higher than a few wavenumbers per atom, and this difference originates from the large overestimation of binding energies in the classical limit owing to the neglect of zero-point corrections. We are now in position to shed light onto the size-dependent behavior of PAH+Hen clusters, by considering atom-resolved indicators δi and σi instead of the averaged quantities, and considering also the distributions of inherent structures rather than their global energy fluctuations δE. For the benzene dopant the sizes n = 6, 24, and 32 were selected. The correlation between the dynamical and structural indicators, the distribution of inherent structures and a pictorial representation70 of the three-dimensional helium density at equilibrium have been represented in Figure 6. The hexamer is an archetypal case of vacancy-mediated disorder.32,56 The most stable configuration has one atom above the benzene plane along the approximate 6-fold axis, while the five others stand on the opposite side, once again along the axis. The four remaining atoms would be classically located on four sites of the hexagonal lattice, but nuclear quantum effects blur this picture by allowing the four atoms to occupy all the 6 nearly equivalent sites, as well as some configurations commensurate with the √3 × √3 structure. In fact, even the central atom fluctuates a lot, a phenomenon that is absent in the pentamer with a single G

dx.doi.org/10.1021/jp510799h | J. Phys. Chem. A XXXX, XXX, XXX−XXX

The Journal of Physical Chemistry A

Article

Figure 7. Atom-resolved dynamical properties of cationic pyrene coated with selected numbers N of helium atoms: (upper panel) spatial extension versus Lindemann index; (middle panel) distribution of inherent structure energies; (lower panel) density plots. (a) N = 8; (b) N = 35; (c) N = 80.

Figure 8. Atom-resolved dynamical properties of cationic coronene coated with selected numbers N of helium atoms: (upper panel) spatial extension versus Lindemann index; (middle panel) distribution of inherent structure energies; (lower panel) density plots. (a) N = 10; (b) N = 32; (c) N = 70.

and liquidlike, while the atoms standing next to the pyrene planes remain localized on the honeycomb lattice. This coexistence is well seen on the correlation between the dynamical and geometrical indicators δi and σi and is associated with multiple inherent structures differing only in the location of the peripheral atoms and all lying about 13 cm−1 above the classical global minimum. The larger cluster at n = 80, above the classical snowball size, has an even greater delocalized character, but still contains localized subparts on both sides of the flat dopant with σ < 1 Å and δ ≈ 0.1. Here the delocalization is again intermediate for atoms in the peripheral region close to the molecular plane (δ ≈ 0.15−0.2 and σ ≈ 2−3 Å), but the atoms above the localized layers are clearly liquidlike with δ > 0.2 and σ > 3 Å. This liquidlike character is manifested by numerous inherent structures extending more than 50 cm−1 above the classical global minimum. The density plot nicely shows the coexistence between the localized first layer and the liquidlike second layer, as well as the intermediate slushy phase closer to the molecular plane. Coronene-coated clusters show comparable behaviors, but yet other phenomenology for selected sizes. Figure 8 shows the atom-resolved properties obtained for n = 10, 32, and 70, the latter size being chosen exactly at the classical snowball size. The 10-atom cluster shows an interesting lowest-energy structure that is preserved upon including harmonic ZPE corrections, and which consists of a rhombic tetramer occupying sites of the √3 × √3 lattice on one side, the six remaining atoms lying epitaxially on the honeycomb lattice with one vacancy. The PIMD simulations predict that all atoms except the central one on the 6-fold symmetry axis are delocalized: In the 6-atom set, the vacancy easily moves between sites, whereas the opposite tetramer is not as constrained as its counterpart on pyrene (in the case of the octamer) and is free to rotate on the coronene surface. The delocalization found on both layers is thus not strictly

characteristic of homogeneous liquidlike phases, as this would entail larger values of δi, but is still associated with several inherent structures that only differ in the relative locations of the opposite layers. The coronene+He32 cluster shows uncomplete layers in epitaxial contact with the honeycomb layer, without any atom in the peripheral belt. The nuclear wave function at 1 K is again that of a mostly localized system (the 14 atoms that have maximum contact with the carbon atoms for which σ ∼ 1 Å) surrounded by more disordered rings in which the vacancies switch between equivalent sites. In the cluster at n = 70 = nC, the initial single-shell structure expands, and about four atoms on each side are extruded and float above this expanded structure. For this system, the equilibrium phase at 1 K again consists of localized and rigidlike first layers parallel to the molecule, delocalized and liquidlike atoms beyond these layers, and an intermediate belt of atoms in a slushy phase. The expansion relative to the classical structure is a natural consequence of vibrational delocalization, and shows here that the classical snowball is not accurately appropriate to describe the quantum mechanical system. Simulations performed at larger sizes confirm the number of atoms remaining in the first shell, with fluctuations amounting to one or two atoms at the maximum, and we are thus able to reevaluate the size of the (quantum mechanical) snowball as nQ = 62 instead of 70. The same procedure leads to nQ = 31, 50, and 97 for the benzene, pyrene, and circumcoronene cations, which approximately represents a drop by more than 10% with respect to nC. The three sizes chosen to illustrate the consequences of coating the circumcoronene cation are n = 18, 80, and 120, and the corresponding properties are represented in Figure 9. At n = 18 the stable structure shows two layers of nine atoms each in epitaxy with the honeycomb lattice. The 7-atom central regions are both occupied, and the four remaining atoms are on H

dx.doi.org/10.1021/jp510799h | J. Phys. Chem. A XXXX, XXX, XXX−XXX

The Journal of Physical Chemistry A

Article

Figure 10. Atom-resolved spatial extension versus average distance from the dopant center for 500 helium atoms coating cationic circumcoronene, as obtained from path-integral molecular dynamics simulations at 1 K. The inset shows the corresponding density plot accumulated along the trajectory.

essentially axially symmetric ellipsoid around the dopant, in agreement with the high symmetry of circumcoronene that is preserved upon ionization. The correlation plot indicates that the most localized atoms are also the closest to the dopant. What this scatter plot further shows is the continuity between the solidlike localization of these atoms in contact with the graphitic cation and the outermost, liquidlike atoms. This continuity proceeds thanks to the peripheral region, which is not sufficiently dense to allow for confinement and can exchange atoms with the more distant liquid atoms. As anticipated, the analysis of the PIMD simulations offers a drastically different picture than the classical global optimization study, which generally fails, even at small sizes, except for some particularly stable, magic-like systems lacking vacancies or defects and that accommodate ideally to the dopant structure in a commensurate way. The strongly size-dependent extent of delocalization revealed by the global indicators hides a broad variety of situations in which delocalization is promoted by rotational or translational disorder. The less dense belt regions closer to the molecular plane are more prone to delocalization, and only for benzene coated by 32 helium atoms did we find a case were the complete solvation shell is fully localized once a structural transition from the global minimum had occurred. Interestingly, benzene coated with slightly more than 32 helium atoms tends to form a 31-atom shell with floating excess atoms. This result may originate from unsufficiently converged simulations, but we believe it is another genuine finite size effect. Repeating indeed the simulations for benzene+He32 after adding one 33-atom outside the belt region leads to the extrusion of one atom from the 32-atom shell, confirming the quantum snowball size to be 31 in large enough droplets, although the special structural stability at n = 32 may be uncovered in mass spectrometry experiments. In general, the slushy phase affecting the peripheral region is best characterized once the snowball has just been formed, and before the complete solvation by additional layers as in large droplets where the distinction with the liquidlike region becomes less clear. It is reminiscent of partially melted phases in the classical thermodynamics of clusters, as occurs in the context of surface melting of van der Waals clusters76−78 or confined ion clouds79 but also for alkali halide clusters.80

Figure 9. Atom-resolved dynamical properties of cationic circumcoronene coated with selected numbers N of helium atoms: (upper panel) spatial extension versus Lindemann index; (middle panel) distribution of inherent structure energies; (lower panel) density plots. (a) N = 18; (b) N = 80; (c) N = 120.

peripheral sites. According to the PIMD trajectories, all atoms undergo significant diffusional motion on the graphitic surface, hopping between the honeycomb lattice sites. This is again a vacancy-induced type of delocalization not found in the symmetric 14-atom cluster and associated with excitations in the energy landscape that are significant. The circumcoronene+He80 cluster is comparable to the pyrene+He35 system, in the sense that its global minimum has fully occupied layers parallel to the dopant and a few extra atoms in the belt separating them. Those extra atoms remain in the belt once nuclear expansion has taken place, and they are more mobile and delocalized than the atoms from the flat layers, allowing the system to explore minima lying up to 7 cm−1 above the global minimum despite contributing significantly to increasing δ or σ. The onset of the classical snowball for the circumcoronene cation is marked by 120 helium atoms, but as was the case for coronene+He70 some atoms overflow and leave the flat first layers that expand as a response to vibrational delocalization. More than 20 atoms are thus left floating in a liquidlike state coexisting with two other rigidlike and slushy phases close to the dopant plane and in the belt around it, respectively. The nearly continuous variations exhibited by the δi and σi indicators confirm the heterogeneity in the dynamics and localization properties of the helium clusters. This analysis is further supported by considering a much larger system having n = 500 helium atoms, still for the circumcoronene cationic dopant. We have represented in Figure 10 the correlation between the atom-resolved delocalization index σi and the average distance ⟨ri⟩ to the center of the dopant, chosen for the additional information it carries with respect to the Lindemann index. For this mediumsize droplet, the vast majority of the atoms belong to a fluidlike phase, but some more localized regions are still identified from the low values of σi below 2 Å. The density plot reveals an I

dx.doi.org/10.1021/jp510799h | J. Phys. Chem. A XXXX, XXX, XXX−XXX

The Journal of Physical Chemistry A

Article

Is is finally relevant to compare the results obtained for the present clusters coating finite planar polyaromatics to previous investigations dealing with the extended graphene substrate.55−57 Such a comparison appears particularly relevant in the case of larger PAHs such as coronene or circumcoronene, for which the partial charges on the central carbon atoms are vanishingly small. The near degeneracy between the localized and fluidlike phases of the 4He bulk monolayer adsorbed on graphene found by Gordillo and Boronat55 is consistent with our finding of such competing phases for finite clusters. Another consistent result is the stabilization of some clusters in the commensurate √3 × √3 lattice, although it seems that epitaxial localization on the honeycomb sites can be favored as well. We speculate that the finite extent of the dopant provides more room for the helium atoms to rearrange perpendicular to it, stabilizing the epitaxial arrangement at least for systems not exhibiting major defects and a symmetry compatible with that of the dopant, for example, for n = 14 in the cases of benzene, coronene, and circumcoronene or n = 8 for pyrene. In their quantum Monte Carlo study at zero and finite temperatures, Vranješ Markić and co-workers57 emphasized the key roles of helium-substrate interactions and in particular corrugation on the equilibrium properties of helium clusters, and also noticed significant size effects. Our results confirm those size effects, and notably show that even 2D clusters below the snowball limit can exhibit fluidlike or localized character depending on the precise number of atoms. In comparison, what the present work specifically highlights is the importance of the boundaries on the extent of delocalization, especially once the two faces of the molecular dopant have both been covered.

liquid, and is washed out in large droplets where the atoms continuously exchange with the outer liquid. Before this situation is reached, medium-size clusters display a peculiar equilibrium state with well-localized regions on both faces of the dopant, liquidlike atoms floating away from them, and a slushy belt of atoms near the molecular plane. It is not clear how the slushy phase would impact the relative stability of the various clusters, and in particular whether the snowball sizes predicted in the present work would be observed in mass spectrometry experiments. Higher relative abundances could well take place before the less stable peripheral atoms start to fill in. Evaluating the ground state energies of the various clusters in the size range leading to the snowball and beyond would be a natural extension of the present work, using more accurate but also more demanding computational methods such as diffusion or path-integral ground state quantum Monte Carlo. More systematic characterization of the collective and individual behaviors of the atoms in the slush phase could be attempted as well, for example, through the evaluation of diffusion constants from ring-polymer molecular dynamics.81 It would also be useful to clarify the boundaries between the various phases, for instance by computing suitable spatial correlation functions. Besides molecular geometry, one further cause for interaction anisotropy lies in the nonuniform repartition of the excess charge on these cationic molecules, central carbon atoms in larger compounds carrying much smaller charges on their way to graphene. It would naturally be interesting to repeat the present calculations on neutral PAHs to see how the results are affected by the global charge. In addition, the case of anionic dopants would also be relevant in experiments, especially in the case of small molecules where the charge effects would be strongest. Previous studies already emphasized the crucial importance of the interactions on the extent of delocalization and confinement,31,57 and although they did not discuss the possible contribution of polarization forces82 it would not be surprising that the snowballs for anions differ from the results in cations, even though the phenomenology of slushy phases remains. Such a result would incidentally pave the way toward some possible experimental verification of such fine effects and contribute to validating ab initio theories or explicit analytical potentials.

4. CONCLUDING REMARKS The localization of helium atoms around ionic impurities or extended substrates is one of the most intriguing features exhibited by this element. As they combine both properties, cationic molecular dopants embedded in helium droplets are natural candidates to produce the so-called snowballs in which the nearby atoms lose their liquidlike character. The present work was aimed at investigating theoretically such particular phases in the case of polycyclic aromatic hydrocarbon solutes ranging from benzene to circumcoronene. Using a combination of classical and semiclassical simulation methods, we could locate the completion of the first solvation shell, or snowball, to occur between 34 and 112 atoms in the classical approximation of nuclear motion, these numbers being decreased by more than 10% once quantum effects are included. The strong anisotropy in the helium−PAH interaction gives rise to complex finite-size effects affecting the collective behavior of the helium cluster around the dopant. Cases of strong localization were found generally related with particularly stable arrangements commensurate with the dopant structure, either on the honeycomb lattice or on the √3 × √3 sites. The presence of vacancies or excess defects usually breaks those rigidlike structures either locally (on the edges) or more globally, as it does in the case of fullerene dopants.32 Rotational delocalization was occasionally found with the central atom in the dopant symmetry axis remaining localized, similarly to earlier computational findings in the case of neutral benzene.47,49 Once the two faces of the dopants are covered, the peripheral region becomes increasingly filled, but this belt is less dense and more prone to delocalization. This slushy phase exhibits properties intermediate between those of a solid and a



ASSOCIATED CONTENT

S Supporting Information *

Optimized geometries and partial atomic charges of the cationic aromatic hydrocarbons, as obtained from DFT calculations. This material is available free of charge via the Internet at http://pubs.acs.org.



AUTHOR INFORMATION

Corresponding Author

*E-mail: fl[email protected]. Notes

The authors declare no competing financial interest.



ACKNOWLEDGMENTS The author wishes to acknowledge generous computational resources from the regional Pôle Scientifique de Modélisation Numérique in Lyon. J

dx.doi.org/10.1021/jp510799h | J. Phys. Chem. A XXXX, XXX, XXX−XXX

The Journal of Physical Chemistry A



Article

Energetics and Dynamics of Helium Adsorbed on Isolated Fullerene Ions. Phys. Rev. Lett. 2012, 108, 076101. (24) Echt, O.; Kaiser, A.; Zöttl, S.; Mauracher, A.; Denifl, S.; Scheier, P. Adsorption of Polar and Non-polar Molecules on Isolated Cationic C60, C70 and Their Aggregates. ChemPlusChem. 2013, 78, 910−920. (25) Hernandez, E. S.; Cole, M. W.; Boninsegni, M. Adsorption of Atoms and Fluids on Spherical Surfaces. Phys. Rev. B 2003, 68, 125418. (26) Galli, D. E.; Buzzacchi, M.; Reatto, L. Pure and Alkali-IonDoped Droplets of 4He. J. Chem. Phys. 2001, 115, 10239. (27) Marinetti, F.; Coccia, E.; Bodo, E.; Gianturco, F. A.; Yurtsever, E.; Yurtsever, M.; Yildirim, E. Bosonic Helium Clusters Doped by Alkali Metal Cations: Interaction Forces and Analysis of Their Most Stable Structures. Theor. Chem. Acc. 2007, 118, 53−65. (28) Coccia, E.; Bodo, E.; Marinetti, F.; Gianturco, F. A.; Yildirim, E.; Yurtsever, E.; Yurtsever, M. Bosonic Helium Droplets with Cationic Impurities: Onset of Electrostriction and Snowball Effects from Quantum Calculations. J. Chem. Phys. 2007, 126, 124319. (29) Galli, D. E.; Ceperley, D. M.; Reatto, L. Path Integral Monte Carlo Study of 4He Clusters Doped with Alkali and Alkali-Earth Ions. J. Phys. Chem. A 2011, 115, 7300−7309. (30) Slavicek, P.; Lewerenz, M. Quantum Solvation and Coordination: Lead Ions Inside Small Helium Droplets. Phys. Chem. Chem. Phys. 2010, 12, 1152−1161. (31) Issaoui, N.; Abdessalem, K.; Habli, H.; Ghalla, H.; Oujia, B.; Calvo, F. Theoretical Investigation of the Relative Stability of Na+Hen (n = 2−24) Clusters: Many-Body versus Delocalization Effects. J. Chem. Phys. 2014, 141, 174316. (32) Calvo, F. Size-Induced Melting and Reentrant Freezing in Fullerene-Doped Helium Clusters. Phys. Rev. B 2012, 85, 060502(R). (33) Shin, H.; Kwon, Y. Commensurate−Incommensurate Transition of 4He Adsorbed on a Single C60 Molecule. J. Chem. Phys. 2012, 136, 064514. (34) Crowell, P. A.; Reppy, J. D. Reentrant Superfluidity in 4He Films Adsorbed on Graphite. Phys. Rev. Lett. 1993, 70, 3291−3294. (35) Smalley, R. E.; Wharton, L.; Levy, D. H.; Chandler, D. W. Fluorescence Excitation Spectrum of s-tetrazine Cooled in a Supersonic Free Jet: Van der Waals Complexes and Isotopic Species. J. Chem. Phys. 1978, 68, 2487−2491. (36) Zwier, T. S.; Carrasquillo, E.; Levy, D. H. The Spectroscopy and Single Vibronic Level Fluorescence Quantum Yields of Jet Cooled Trans-Stilbene and its Van der Waals Complexes. J. Chem. Phys. 1983, 78, 5493−5505. (37) Semmes, D. H.; Baskin, J. S.; Zewail, A. H. Picosecond RealTime Studies of Mode-Specific Vibrational Predissociation. J. Chem. Phys. 1990, 92, 3359−3376. (38) Yu, L.; Williamson, J.; Foster, S. C.; Miller, T. A. High Resolution Laser Spectroscopy of Free Radical−Inert Gas Complexes: C5H5·He, C5H5·He2, C5H5·Ne, and CH3−C5H4·He2. J. Chem. Phys. 1992, 97, 5273−5282. (39) Coutant, B.; Bréchignac, P. Anomalous Complex Shift of LowFrequency Out-of-Plane Vibrations in Aniline-M Van der Waals Complexes (M = He, Ne, Ar). J. Chem. Phys. 1994, 100, 7087−7092. (40) Bach, A.; Leutwyler, S.; Sabo, D.; Bačić, Z. Very Large Amplitude Intermolecular Vibrations and Wave Function Delocalization in 2,3-Dimethylnaphthalene·He Van der Waals Complex. J. Chem. Phys. 1997, 107, 8781−8793. (41) Hartmann, M.; Lindinger, A.; Toennies, J. P.; Vilesov, A. F. Laser-Induced Fluorescence Spectroscopy of Van der Waals Complexes of Tetracene−ArN (N ≤ 5) and Pentacene−Ar within Ultracold Liquid He Droplets. Chem. Phys. 1998, 239, 139−149. (42) Even, U.; Jortner, J.; Noy, D.; Lavie, N.; Cossart-Magos, C. Cooling of Large Molecules below 1 K and He Clusters Formation. J. Chem. Phys. 2000, 112, 8068−8071. (43) Lehnig, R.; Slenczka, A. Microsolvation of Phthalocyanines in Superfluid Helium Droplets. ChemPhysChem 2004, 5, 1014−1019. (44) Schmied, R.; Carçabal, P.; Doktor, A. M.; Lonij, V. P. A.; Lehmann, K. K.; Scoles, G. UV Spectra of Benzene Isotopomers and Dimers in Helium Nanodroplets. J. Chem. Phys. 2004, 121, 2701− 2710.

REFERENCES

(1) Atkins, K. R. Ions in Liquid Helium. Phys. Rev. 1959, 116, 1339− 1343. (2) Toennies, J. P.; Vilesov, A. F. Superfluid Helium Droplets: A Uniquely Cold Nanomatrix for Molecules and Molecular Complexes. Angew. Chem., Int. Ed. 2004, 43, 2622−2648. (3) Krasniokutski, S. A.; Huisken, F. Ultra-Low-Temperature Reactions of Mg Atoms with O2 Molecules in Helium Droplets. J. Phys. Chem. A 2010, 114, 7292−7300. (4) Mestdagh, J.-M.; Gaveau, M.-A.; Gée, C.; Sublemontier, O.; Visticot, J.-P. Cluster Isolated Chemical Reactions. Int. Rev. Phys. Chem. 1997, 16, 215−247. (5) Sindzingre, P.; Klein, M. L.; Ceperley, D. M. Path-Integral Monte Carlo Study of Low-Temperature 4He Clusters. Phys. Rev. Lett. 1989, 63, 1601−1604. (6) Binder, K. Finite Size Effects on Phase Transitions. Ferroelectrics 1987, 73, 43−67. (7) Toennies, J. P.; Vilesov, A. F. Spectroscopy of Atoms and Molecules in Liquid Helium. Annu. Rev. Phys. Chem. 1998, 49, 1−41. (8) Nauta, K.; Miller, R. E. Stark Spectroscopy of Polar Molecules Solvated in Liquid Helium Droplets. Phys. Rev. Lett. 1999, 82, 4480− 4483. (9) Callegari, C.; Conjusteau, A.; Reinhard, I.; Lehmann, K. K.; Scoles, G.; Dalfovo, F. Superfluid Hydrodynamic Model for the Enhanced Moments of Inertia of Molecules in Liquid 4He. Phys. Rev. Lett. 1999, 83, 5058−5061. (10) Stienkemeier, F.; Lehmann, K. K. Spectroscopy and Dynamics in Helium Nanodroplets. J. Phys. B.: At., Mol. Opt. Phys. 2006, 39, R127−R166. (11) Henne, U.; Toennies, J. P. Electron Capture by Large Helium Droplets. J. Chem. Phys. 1998, 108, 9327−9338. (12) Kojima, T. M.; Kobayashi, N.; Kaneko, Y. Helium Cluster Ions RgHe+x (Rg = Ne, Ar and Kr, x ≤ 14) Formed in a Drift Tube Cooled by Liquid Helium. Z. Phys. D 1992, 22, 645−650. (13) Callicoatt, B. E.; Forde, K.; Ruchti, T.; Jung, L. L.; Janda, K. C.; Halberstadt, N. Capture and Ionization of Argon within Liquid Helium Droplets. J. Chem. Phys. 1998, 108, 9371−9382. (14) Kim, J. H.; Peterka, D. S.; Wang, C. C.; Neumark, D. M. Photoionization of Helium Nanodroplets Doped with Rare Gas Atoms. J. Chem. Phys. 2006, 124, 214301. (15) Bartl, P.; Leidlmair, C.; Denifl, S.; Scheier, P.; Echt, O. On the Size and Structure of Helium Snowballs Formed around Charged Atoms and Clusters of Noble Gases. J. Phys. Chem. A 2014, 118, 8050−8059. (16) Schöbel, H.; Bartl, P.; Leidlmair, C.; Denifl, S.; Echt, O.; Märk, T.; Scheier, P. High-Resolution Mass Spectrometric Study of Pure Helium Droplets, and Droplets Doped with Krypton. Eur. Phys. J. D 2011, 63, 209−214. (17) Ferreira da Silva, F.; Bartl, P.; Denifl, S.; Echt, O.; Märk, T. D.; Scheier, P. Argon Clusters Embedded in Helium Nanodroplets. Phys. Chem. Chem. Phys. 2009, 11, 9791−9797. (18) An der Lan, L.; Bartl, P.; Leidlmair, C.; Jochum, R.; Denifl, S.; Echt, O.; Scheier, P. Solvation of Na+, K+ and Their Dimers in Helium. Chem.Eur. J. 2012, 18, 4411−4418. (19) Müller, S.; Mudrich, M.; Stienkemeier, F. Alkali-Helium Snowball Complexes Formed on Helium Nanodroplets. J. Chem. Phys. 2009, 131, 044319. (20) Theisen, M.; Lackner, F.; Ernst, W. E. Cs Atoms on Helium Nanodroplets and the Immersion of Cs+ into the Nanodroplet. J. Chem. Phys. 2011, 135, 074306. (21) Bartl, P.; Leidlmair, C.; Denifl, S.; Scheier, P.; Echt, O. Cationic Complexes of Hydrogen with Helium. ChemPhysChem 2012, 14, 227− 232. (22) Ferreira da Silva, F.; Waldburger, P.; Jaksch, S.; Mauracher, A.; Denifl, S.; Echt, O.; Märk, T. D.; Scheier, P. On the Size of Ions Solvated in Helium Clusters. Chem.Eur. J. 2009, 15, 7101−7108. (23) Leidlmair, C.; Wang, Y.; Bartl, P.; Schöbel, H.; Denifl, S.; Probst, M.; Alcamı ́, M.; Martı ́n, F.; Zettergen, H.; Hansen, K.; et al. Structures, K

dx.doi.org/10.1021/jp510799h | J. Phys. Chem. A XXXX, XXX, XXX−XXX

The Journal of Physical Chemistry A

Article

(45) Birer, Ö .; Moreschini, P.; Lehmann, K. K. Electronic Spectroscopy of Benzo[g,h,i]perylene and Coronene inside Helium Nanodroplets. Phys. Chem. Chem. Phys. 2008, 10, 1648−1657. (46) Pentlehner, D.; Slenczka, A. Microsolvation of Anthracene inside Superfluid Helium Nanodroplets. Mol. Phys. 2012, 15−16, 1933−1940. (47) Kwon, Y.; Whaley, K. B. Localization of Helium at an Aromatic Molecule in Superfluid Helium Clusters. J. Chem. Phys. 2000, 114, 3163−3169. (48) Heidenreich, A.; Even, U.; Jortner, J. Nonrigidity, Delocalization, Spatial Confinement and Electronic-Vibrational Spectroscopy of Anthracene−Helium Clusters. J. Chem. Phys. 2001, 115, 10175− 10185. (49) Huang, P.; Whaley, K. B. Localized Helium Excitations in 4 HeN−Benzene Clusters. Phys. Rev. B 2003, 67, 155419. (50) Whitley, H. D.; Huang, P.; Kwon, Y.; Whaley, K. B. Multiple Solvation Configurations around Phthalocyanine in Helium Droplets. J. Chem. Phys. 2005, 123, 054307. (51) Whitley, H. D.; DuBois, J. L.; Whaley, K. B. Spectral Shifts and Helium Configurations in 4HeN−Tetracene Clusters. J. Chem. Phys. 2009, 131, 124514. (52) Fetzer, J. Polycyclic Aromatic Hydrocarbons and Combustion, in Handbook of Combustion; Wiley-VCH Verlag (2010). (53) Léger, A.; Puget, J.-L. Identification of the Unidentified IR emission Feature of Interstellar Dust. Astron. Astrophys. 1984, 137, L5−L8. (54) Gschliesser, D.; Ferreira da Silva, F.; Bartl, P.; Denifl, S.; Scheier, P. Synthesis of Polycyclic Aromatic Hydrocarbons in He Nanodroplets. J. Phys.: Conf. Ser. 2012, 388, 052053. (55) Gordillo, M. C.; Boronat, J. 4He on a Single Graphene Sheet. Phys. Rev. Lett. 2009, 102, 085303. (56) Gordillo, M. C.; Cazorla, C.; Boronat, J. Supersolidity in Quantum Films Adsorbed on Graphene and Graphite. Phys. Rev. B 2010, 83, 121406(R). (57) Vranješ Markić, L.; Stipanović, P.; Bešlić, I.; Zillich, R. E. 4He Clusters Adsorbed in Graphene. Phys. Rev. B 2013, 88, 125416. (58) Stillinger, F. H.; Weber, T. A. Hidden Structure in Liquids. Phys. Rev. A 1982, 25, 978−989. (59) Frisch, M. J.; Trucks, G. W.; Schlegel, H. B.; Scuseria, G. E.; Robb, M. A.; Cheeseman, J. R.; Scalmani, G.; Barone, V.; Mennucci, B.; Petersson, G. A.; Nakatsuji, H.; Caricato, M.; Li, X.; Hratchian, H. P.; Izmaylov, A. F.; Bloino, J.; Zheng, G.; Sonnenberg, J. L.; Hada, M.; Ehara, M.; Toyota, K.; Fukuda, R.; Hasegawa, J.; Ishida, M.; Nakajima, T.; Honda, Y.; Kitao, O.; Nakai, H.; Vreven, T.; Montgomery, J. A., Jr.; Peralta, J. E.; Ogliaro, F.; Bearpark, M.; Heyd, J. J.; Brothers, E.; Kudin, K. N.; Staroverov, V. N.; Kobayashi, R.; Normand, J.; Raghavachari, K.; Rendell, A.; Burant, J. C.; Iyengar, S. S.; Tomasi, J.; Cossi, M.; Rega, N.; Millam, J. M.; Klene, M.; Knox, J. E.; Cross, J. B.; Bakken, V.; Adamo, C.; Jaramillo, J.; Gomperts, R.; Stratmann, R. E.; Yazyev, O.; Austin, A. J.; Cammi, R.; Pomelli, C.; Ochterski, J. W.; Martin, R. L.; Morokuma, K.; Zakrzewski, V. G.; Voth, G. A.; Salvador, P.; Dannenberg, J. J.; Dapprich, S.; Daniels, A. D.; Farkas, O.; Foresman, J. B.; Ortiz, J. V.; Cioslowski, J.; Fox, D. J. Gaussian 09, revision D.01; Gaussian Inc.: Wallingford, CT, 2014. (60) Janzen, A. R.; Aziz, R. A. An Accurate Potential Energy Curve for Helium Based on ab Initio Calculations. J. Chem. Phys. 1997, 107, 914−919. (61) Varandas, A. J. C. Helium−Fullerene Pair Interactions: An ab Initio Study by Perturbation Theory and Coupled Cluster Methods. Int. J. Quantum Chem. 2011, 111, 416−429. (62) Lim, K. F. Quasiclassical Trajectory Study of Collisional Energy Transfer in Toluene Systems. II. Helium Bath Gas: Energy and Temperature Dependences, and Angular Momentum Transfer. J. Chem. Phys. 1994, 101, 8756−8787. (63) Masili, M.; Starace, A. F. Static and Dynamic Dipole Polarizability of the Helium Atom Using Wave Functions Involving Logarithmic Terms. Phys. Rev. A 2003, 68, 012508.

(64) Wales, D. J.; Doye, J. P. K. Global Optimization by BasinHopping and the Lowest Energy Structures of Lennard-Jones Clusters Containing up to 110 Atoms. J. Phys. Chem. A 1997, 101, 5111−5116. (65) Pérez, A.; Tuckerman, M. E.; Müser, M. H. A Comparative Study of the Centroid and Ring-Polymer Molecular Dynamics Methods for Approximating Quantum Time Correlation Functions from Path Integrals. J. Chem. Phys. 2009, 130, 184105. (66) Witt, A.; Ivanov, S. D.; Forbert, H.; Marx, D. On the Applicability of Centroid and Ring Polymer Path Integral Molecular Dynamics for Vibrational Spectroscopy. J. Chem. Phys. 2009, 130, 194510. (67) Markland, T. E.; Manolopoulos, D. E. An Efficient Ring Polymer Contraction Scheme for Imaginary Time Path Integral Simulations. J. Chem. Phys. 2008, 129, 024105. (68) Oleksy, K.; Karlický, F.; Kalus, R. Structures and Energetics of Helium Cluster Cations: Equilibrium Geometries Revisited through the Genetic Algorithm Approach. J. Chem. Phys. 2010, 133, 164314. (69) Calvo, F.; Naumkin, F. Y.; Wales, D. J. Interplay between Charge and Vibrational Delocalization in Cationic Helium Clusters. J. Chem. Phys. 2011, 135, 124308. (70) The 3D density plots have been constructed using the Mayavi2 suite of Python scripts Ramachandran, P.; Varoquaux, G. M. Mayavi: 3D Visualization of Scientific Data. IEEE Comput. Sci. Eng. 2011, 13, 40−51. (71) Angelani, L.; Parisi, G.; Ruocco, G.; Viliani, G. Connected Network of Minima as a Model Glass: Long Time Dynamics. Phys. Rev. Lett. 1998, 81, 4648−4651. (72) Chowdhary, J.; Keyes, T. Thermodynamics and Dynamics for a Model Potential Energy Landscape. J. Phys. Chem. B 2004, 108, 19786−19798. (73) Heuer, A.; Doliwa, B.; Saksaengwijit, A. Potential-Energy Landscape of a Supercooled Liquid and Its Resemblance to a Collection of Traps. Phys. Rev. E 2005, 72, 021503. (74) Pierce, M. E.; Manousakis, E. Monolayer Solid 4He Clusters on Graphite. Phys. Rev. Lett. 1999, 83, 5314−5317. (75) Calvo, F.; Doye, J. P. K.; Wales, D. J. Quantum Partition Functions from Classical Distributions. Application to Rare-Gas Clusters. J. Chem. Phys. 2001, 114, 7312−7329. (76) Briant, C. L.; Burton, J. J. Molecular Dynamics Study of the Structure and Thermodynamic Properties of Argon Microclusters. J. Chem. Phys. 1975, 63, 2045−2058. (77) Cheng, H.-P.; Berry, R. S. Surface Melting of Clusters and Implications for Bulk Matter. Phys. Rev. A 1994, 49, 7969−7980. (78) Calvo, F.; Labastie, P. Evidence for Surface Melting in Clusters Made of Double Icosahedron Units. Chem. Phys. Lett. 1996, 258, 233− 238. (79) Schiffer, J. P. Melting of Crystalline Confined Plasmas. Phys. Rev. Lett. 2002, 88, 205003. (80) Calvo, F.; Labastie, P. Melting and Phase Space Transitions in Small Ionic Clusters. J. Phys. Chem. B 1998, 102, 2051−2059. (81) Craig, I. R.; Manolopoulos, D. E. Quantum Statistics and Classical Mechanics: Real Time Correlation Functions from Ring Polymer Molecular Dynamics. J. Chem. Phys. 2004, 121, 3368−3373. (82) González, B. S.; Hernández-Rojas, J.; Bretón, J.; Gomez Llorente, J. M. Global Potential Energy Minima of (H2O)n Clusters on Graphite. J. Phys. Chem. C 2007, 111, 14862−14869.

L

dx.doi.org/10.1021/jp510799h | J. Phys. Chem. A XXXX, XXX, XXX−XXX