Structural Features of Small Benzene Clusters (C6H6) - American

Article pubs.acs.org/JPCA

Structural Features of Small Benzene Clusters (C6H6)n (n ≤ 30) As Investigated with the All-Atom OPLS Potential Hiroshi Takeuchi* Division of Chemistry, Graduate School of Science, Hokkaido University, 060-0810 Sapporo, Japan S Supporting Information *

ABSTRACT: The structures of the simplest aromatic clusters, benzene clusters (C6H6)n, are not well elucidated. In the present study, benzene clusters (C6H6)n (n ≤ 30) were investigated with the all-atom optimized parameters for liquid simulation (OPLS) potential. The global minima and low-lying minima of the benzene clusters were searched with the heuristic method combined with geometrical perturbations. The structural features and growth sequence of the clusters were examined by carrying out local structure analyses and structural similarity evaluation with rotational constants. Because of the anisotropic interaction between the benzene molecules, the local structures consisting of 13 molecules are considerably deviated from regular icosahedron, and the geometries of some of the clusters are inconsistent with the shapes constructed by the interior molecules. The distribution of the angle between the lines normal to two neighboring benzene rings is anisotropic in the clusters, whereas that in the liquid benzene is nearly isotropic. The geometries and energies of the low-lying configurations and the saddle points between them suggest that most of the configurations previously detected in supersonic expansions take different orientations for one to four neighboring molecules.

1. INTRODUCTION Since interactions including π electrons are important in biomolecules and molecular crystals,1−3 a lot of investigations on the benzene dimer have been performed (see refs 1−6 and references cited therein). Lee et al.1 examined aromatic−aromatic interactions in the benzene dimer and substituted benzene− benzene dimers with coupled cluster theory calculations. The results show that the T-shaped configuration of the benzene dimer is slightly more stable than the parallel displaced one. The dominant interaction in the aromatic−aromatic interaction is due to the dispersive force but the dispersion energy tends to be canceled out by the exchange repulsion energy. The energies of configurations mainly depend on the electrostatic energies, and thereby, the relative stability of the configurations is closely related to the polarity of their surroundings. In apolar environments, the parallel displaced form is preferred over the T-shaped one. This indicates that the relative stability of the parallel displaced configuration of two benzene rings in the hydrophobic core of a protein decreases with the increasing polarity of the surroundings.1 The energy difference between the T-shaped and parallel displaced configurations of the dimer is small (∼0.4 kJ mol−1),1,2 and the existence of substituents on the benzene ring leads to preference of the parallel displaced configuration over the T-shaped one.1 This generally causes the preference of the parallel displaced structures of substituted benzene rings in organic crystals. According to the experimental investigations,4−6 the benzene dimer takes a floppy T-shaped structure in agreement with the theoretical investigations.1−3 Benzene clusters with more than © XXXX American Chemical Society

two molecules (C6H6)n have been investigated to reveal the structures governed by the aromatic−aromatic interactions. The benzene trimer and tetramer were examined by mass-selective resonantly enhanced two-photon ionization (R2PI) and ultraviolet−ultraviolet hole burning spectroscopy.7,8 The experimental data are consistent with a cyclic trimer structure with C3 symmetry7 and a distorted tetrahedral structure with S4 symmetry.8 For benzene clusters with n = 7−19, Easter and coworkers9−13 studied these structures formed in free jet expansions by R2PI spectroscopy. The numbers of configurations of the clusters with 13−19 molecules were estimated to be one to three referring to the transitions observed in the spectra of isotopically substituted benzene clusters.11 Moreover, the configurations were interpreted with aid of the global-minimum geometries of the Lennard-Jones atomic clusters. However, orientational degrees of freedom of molecules lacking in the study must be taken into account to interpret the structures of the benzene clusters. Many theoretical investigations14−30 on benzene clusters have been reported. Ab initio calculations are limited to relatively small clusters (five or less molecules)28−30 because of huge computational cost. Recently, Gadre and co-workers17 extended ab initio calculations to benzene clusters with up to 10 molecules; a fragment-based-strategy was employed to study the growth sequence of the clusters at the MP2 level. For the 13 molecule cluster, Easter14,15 carried out the Monte Carlo simulations with Received: June 17, 2012 Revised: September 12, 2012

A

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

The Journal of Physical Chemistry A

Article

In the equation, k and l represent atoms in the molecules i and j, respectively, and rkl means the interatomic distance. The values of the potential coefficients used in this work are listed in Table 1.

several model potentials and found that the existence of two distinct configurations explains the experimental results. Chakrabarti et al.16 examined the potential energy surface of (C6H6)13 using the model potential developed from symmetry-adapted perturbation theory based on density functional theory. The obtained results support the conclusion reported by Easter.14,15 The present author developed an efficient method for geometry optimization of molecular clusters and applied it to benzene clusters to examine its efficiency.24 In the study, the potential function of Williams and Starr (WS)31 was employed since the benzene cluster (C6H6)n (n ≤ 15) expressed by the WS potential has been regarded as a benchmark problem20,21 for optimization methods. It was possible to find new global minima for (C6H6)11, (C6H6)14, and (C6H6)15, and the minima for n = 16−30 were first reported. The global minima proposed in ref 24 were confirmed by Llanio-Trujillo et al.25 with a different optimization method. Recently, Fu and Tian26 carried out molecular dynamics (MD) simulations for the liquid benzene with eight potentials consisting of Lennard-Jones and Coulomb terms. By referring to the results observed by high-resolution neutron diffraction for the liquid,32 they recommended the all-atom optimized parameters for liquid simulation (OPLS) potential proposed by Jorgensen and Severance27 as the best model. According to the recommendation, the present author calculated the lowest-energy geometry of the benzene dimer with the all-atom OPLS (OPLS-AA) model. The geometry is shown in Figure 1. It is consistent with the tilted T geometry reported in

Table 1. Coefficients in the OPLS-AA Potential atom pair

B (kJ mol−1 Å12)

C (kJ mol−1 Å6)

Q (kJ mol−1 Å)

C···C C···H H···H

4.69343 × 106 3.08336 × 105 2.02562 × 104

2.34488 × 103 4.86287 × 102 1.00847 × 102

18.37422 −18.37422 18.37422

The rigid structure with D6h symmetry is assumed for the benzene molecule; the bond lengths r(C−C) and r(C−H) are 1.40 and 1.08 Å, respectively.27 To show anisotropy of the intermolecular potential, potential energies of the benzene dimer were calculated keeping molecular orientations in the T-shaped and parallel (face-to-face) forms and varying the intermolecular distance defined by the centers of mass of the two molecules. The configurations and the calculated potential curves are shown in Figure 2. The T-shaped

Figure 2. Intermolecular potentials of the T-shaped and parallel configurations of the dimer calculated with the OPLS-AA model. The relative orientations of the molecules are fixed, and the intermolecular distance r is varied.

configuration shows a deeper valley than the parallel one. The minimum energies of the parallel and T-shaped configurations are −7.09 and −9.01 kJ mol−1, respectively. These values are consistent with the values reported in the literature,26 −1.69 and −2.15 kcal mol−1. The two configurations in Figure 2 do not correspond to local minima on the potential energy surface. In the studies by Fu and Tian26 and Headen et al.,32 perpendicular configurations are divided into the T-shaped and Y-shaped ones where one C−H and two C−H bonds of one molecule are directed toward the ring of the other molecule, respectively. According to the classification, the most stable configuration of the dimer (Figure 1) is the Y-shaped one. 2.2. Geometry Optimization. At present, the global optimization of atomic/molecular clusters is one of the most difficult problems. Recently, the present author encountered a difficulty relating to molecular orientations in ethane clusters33 and modified the optimization algorithm. The global optimization is not performed routinely, and many articles have been published to improve efficiency of the global optimization.25,33−38 Optimal geometries of the WS benzene clusters were successfully searched with the heuristic method combined with geometrical perturbations (HMGP).24 Accordingly, the same method was applied to the OPLS-AA clusters. Starting geometries were generated randomly and refined with geometricalperturbation operators followed with a quasi-Newton method

Figure 1. Stereographic views of the global minima of the benzene dimer described by the Williams−Starr (WS) and OPLS-AA potentials.

the literature27 but is considerably different from the WS dimer geometry.20,23 In the present study, therefore, the geometries of the benzene clusters (C6H6)n with n ≤ 30 were reinvestigated with the OPLS-AA model. The main purpose of the present study is to compare the results obtained in the present study with those found in the previous investigations11,24,26,32 to report structural features of the benzene clusters.

2. CALCULATIONS 2.1. All-Atom OPLS-AA Potential. The potential energy of (CH4)n takes the pairwise-additive form n−1

Vn =

n

∑ ∑

V ( i , j) (1)

i=1 j=i+1

Here, V(i, j) denotes the potential energy between molecules i and j and is expressed as V (i , j ) =

⎡B kl

∑∑⎢ k

l

⎣ rkl

12

−

Ckl rkl

6

+

Q kl ⎤ ⎥ rkl ⎦

(2) B



Article

(the L-BFGS method39). The number of the starting geometries is 20 for n = 2, 100 for n = 3−15, 200 for n = 16−20, 500 for n = 21−24, 1000 for n = 25−29, and 3000 for n = 30, respectively. For each of the clusters with n = 2−5, the global minimum was easily located after local optimizations39 of starting geometries because several optimizations yielded the same lowest energy. However, it was difficult to obtain the same lowest energy for a larger cluster after local optimizations of starting geometries. Hence, geometrical-perturbation operators (the interior, surface, and orientation operators)24 with subsequent local optimizations were used for searching lower-energy configurations. The interior operator perturbs a cluster configuration by moving some molecules to the neighborhood of the center of mass of a cluster. The surface operator partially modifies a cluster configuration by moving them to the most stable positions on the surface of a cluster. The orientation operator randomizes orientations of molecules in a cluster. The global minimum of each of the clusters with n ≥ 6 was located at least 4 times by means of HMGP. Calculation was executed with dual core 3 GHz Intel Xeon 5160 processors; some calculations were independently performed, and one core was used for one job (no parallel mode was employed). The lowest energies of the clusters are listed in Table 2. The potential energy of the dimer is consistent with the

the local structures with Nmax = 12 is also shown in the same figure. In the experimental study on the liquid benzene,32 the number of the nearest neighbor molecules around a molecule is approximately 12. This is in good agreement with the above result and the results obtained by the MD simulations.26 In the literature,26 the populations of the parallel and perpendicular configurations of two nearest neighbor molecules are reported. Comparison of the configurations in the liquid with those in the clusters may clarify features of the benzene clusters. Hence, molecular pairs with intermolecular distances smaller than 6 Å were classified into parallel and perpendicular configurations. The angle θ is defined by the angle between the two vectors normal to the benzene rings. Parallel configurations have the condition of 0° ≤ θ ≤ 40°, whereas perpendicular configurations satisfy 50° ≤ θ ≤ 90°. The results of the classification are shown in Figure 5. 2.4. Structural Similarity Analyzed with Rotational Constants. Geometrical features are also obtained from gm gm rotational constants, Agm n , Bn , and Cn , of the global-minimum geometry of (C6H6)n. The asymmetry parameter κ (= (2B − A − C)/(A − C))40 was calculated from the rotational constants to identify shapes of the clusters. The results are shown in Figure 6. The cluster with n = 13 is prolate (B = C), and the clusters with n = 2, 5, 8, 14, 15, 18−22, 25, and 28−30 are nearly prolate. Oblate structures (A = B) are found for n = 3, 4, and 26, and the n = 6, 7, 9, 12, and 23 clusters are nearly oblate. The clusters with n = 10, 11, 16, 17, 24, and 27 take nearly spherical shapes. The rotational constants are used to evaluate the similarity between the global-minimum geometries of (C6H6)n and (C6H6)n+1. The following analysis was performed for the geometries of the clusters. A molecule was removed from the (n + 1) molecule cluster, and the rotational constants were calculated from the resulting geometry. The obtained constants, gm gm gm An+1→n , Bn+1→n , and Cn+1→n , were compared with the rotational gm gm constants of (C6H6)n, Agm n , Bn , and Cn , using the equation

Table 2. Lowest Potential Energies (in kJ mol−1) of Benzene Clusters (C6H6)n Obtained with the OPLS-AA Model n

Vn

n

Vn

n

Vn

n

Vn

2 3 4 5 6 7 8 9

−9.693 −29.066 −51.984 −74.206 −100.554 −126.298 −153.621 −183.052

10 11 12 13 14 15 16

−211.953 −241.068 −273.667 −312.932 −341.339 −373.898 −405.224

17 18 19 20 21 22 23

−435.644 −470.118 −507.571 −540.559 −572.053 −607.625 −645.808

24 25 26 27 28 29 30

−678.308 −712.947 −751.957 −784.659 −817.795 −854.537 −888.906

⎡⎛ gm ⎞2 ⎛ Bngm − Bngm ⎞2 A n − A ngm +1→n +1→n ⎢ dn , n + 1 = ⎜ ⎟ +⎜ ⎟ ⎢⎣⎝ A ngm Bngm ⎠ ⎝ ⎠

value reported in ref 27 (−2.31 kcal mol−1). The lowest energy of (C6H6)13 in ref 15, −312.873 kJ mol−1, is slightly higher than that in Table 2 because of differences in the potential coefficients adopted in the calculations. The global-minimum configurations are shown in Figure 3. The Cartesian coordinates of the globalminimum configurations of the OPLS-AA clusters are reported in the Supporting Information. The geometrical data of the other low-lying configurations are available from the author upon request. Average computational times required for one hit of the global minimum of each of the n = 10, 20, and 30 clusters were approximately 21 min, 5 h, and 222 h on a single core, respectively. 2.3. Local Structure Analyses. To clarify structural features of the clusters, local structures were analyzed using relative positions of the centers of mass of molecules (CMMs). For each of the CMMs in a cluster, distances of other CMMs from it were calculated, and the number of the distances smaller than a cutoff distance was obtained. The cutoff distance was arbitrarily set at a value larger than the intermolecular distance of the equilibrium structure of the dimer; since the equilibrium distance is 4.94 Å, a value of 6 Å was used for the cutoff distance. From the numbers for all CMMs in a cluster, the maximum number Nmax was determined. The obtained results are plotted in Figure 4. The value of Nmax is 12 for n ≥ 13. The number of

⎛ C gm − Cngm ⎞2 ⎤ +1→n ⎥ +⎜ n ⎟ Cngm ⎝ ⎠ ⎥⎦

1/2

(3)

If the distance dn,n+1 is zero, the geometry of (C6H6)n+1 is constructed from that of (C6H6)n and an additional benzene molecule. This indicates that the structural change from (C6H6)n to (C6H6)n+1 is continuous. All the ways of deleting a molecule from (n + 1) molecules were tried, and the smallest value of dn,n+1 was adopted as a measure of geometrical similarity (see Figure 7). The parameter dn,n+1 is an indicator sensitive to the structural difference between two clusters. For example, when the molecule on the extreme right of the n = 14 cluster is removed, the resulting geometry is similar to the geometry of the n = 13 cluster (see Figure 3), and then, the value of d13, 14 is 0.01. The corresponding value for the n = 14 cluster where another molecule is removed is 0.10 to 0.24. To discuss the effect of the intermolecular potentials on the geometries, the global-minimum configurations of the WS clusters with n ≥ 424 were compared with the configurations of the OPLS-AA clusters. It should be noted that the WS and OPLS-AA clusters with the same size take different equilibrium intermolecular distances because of the potential difference. C



Article

Figure 3. Stereographic views of the global-minimum geometries of the benzene clusters for n = 3−30.

Here, rWS denote the vectors directed to the ⃗ MM and rWS,corrected ⃗ MM C C uncorrected and corrected CMM in the WS cluster, and WS rOPLS−AA CMM...CMM,ave and rCMM...CMM,ave are the average values of all the intermolecular distances in the OPLS-AA and WS clusters,

Hence, positions of CMMs in the WS cluster were corrected for the potential difference by using the equation WS,corrected OPLS−AA WS WS rCMM = (rCMM ⃗ ⃗ ···CMM,ave/ rCMM···CMM,ave) × rCMM

(4) D



Article

Figure 4. Results of local structure analysis for the global-minimum geometries of the benzene clusters: (1) the maximum number Nmax of nearest neighboring molecules surrounding a molecule; (2) the number of local structures with Nmax = 12.

Figure 7. Geometrical distances between the OPLS-AA global-minimum geometries of (C6H6)n and (C6H6)n+1. These values are calculated according to eq 3.

Figure 8. Geometrical distances between the OPLS-AA and Williams− Starr global-minimum geometries with the same cluster size. These values are calculated according to eq 5.

Figure 5. Populations X of the perpendicular and parallel configurations of the molecular pairs with intermolecular distances smaller than 6 Å in the global-minimum geometries of the benzene clusters. See the text for the definitions of the perpendicular and parallel configurations.

2.5. Potential Barriers between Low-Lying Configurations. The experimental study11 estimated the numbers of the stable configurations of the clusters for n = 13−19. To discuss the numbers in this article, it is necessary to examine the potential barriers between the lowest-energy configurations. These were calculated from the geometries of the two configurations (end points) using the nudged elastic band (NEB) method proposed by Henkelman and Jónsson.41 The linear interpolation is used to generate the value of the geometrical parameter i of nth intermediate configuration ⎛ n⎞ n g int , i(n) = ⎜1 − ⎟gend1, i + gend2, i ⎝ ⎠ m m

Here, m is a constant (18 or 36), and n represents integers between 1 and m − 1. The geometrical parameters i of the two end points (Cartesian coordinates of CMMs and Euler angles of the molecules) are denoted by gend1,i and gend2,i, respectively. It is crucial to carefully prepare the geometrical data of the two end points. If the data are improper, the interpolation generates unreasonable geometries for the intermediates.

Figure 6. Asymmetry parameter κ calculated from the global-minimum geometries of the benzene clusters. The values of the parameter are calculated from the rotational constants (see text).

respectively. After the correction, the similarity between the two geometries was calculated with the equation

3. DISCUSSION 3.1. Relative Stability of Benzene Clusters. To calculate the relative stability of the clusters, the potential energies (Table 2) were substituted into the equation Sn = Vn − 1 + Vn + 1 − 2Vn (7)

dOPLS−AA,WS 2 ⎡⎛ OPLS−AA ⎛ BOPLS−AA − B WS ⎞2 − A WS ⎞ A = ⎢⎜ n OPLS−AA n ⎟ + ⎜ n OPLS−AA n ⎟ ⎢⎝ An Bn ⎠ ⎝ ⎠ ⎣

⎛ C OPLS−AA − C WS ⎞2 ⎤ + ⎜ n OPLS−AA n ⎟ ⎥ Cn ⎝ ⎠ ⎥⎦

(6)

1/2

Since Sn represents the curvature of Vn, a positive value of Sn indicates that the n molecule cluster is more stable than the clusters with the sizes of n ± 1. The values of Sn obtained for the OPLS-AA and WS clusters are plotted in Figure 9. The OPLS-AA

(5)

The results are shown in Figure 8. E



Article

3.3. Structures of All-Atom OPLS Clusters. The OPLS-AA model shows that the dimer evolves to triangular trimer and tetrahedral tetramer (see Figure 3). These configurations take C3 and S4 symmetry, respectively, in accordance with the experimental results.7,8 However, the MP2/6-31G calculations30 locate a trigonal C3h structure and a tetrahedral C3 structure. This suggests that the OPLS-AA results are more reliable than the MP2 calculations. As described by Sinnokrot and Sherrill,28 calculations beyond the MP2 level are necessary to obtain accurate results for intermolecular interaction energies. As mentioned in the introduction, Gadre and co-workers17 investigated benzene clusters with n ≤ 10 by the ab initio calculations at the MP2 level. From a lot of geometries generated for the clusters with n = 3−8, some selected geometries were optimized with the 6-31G* or 6-31++G** basis sets. For the n = 9 and 10 clusters, the geometry optimization was performed with the 6-31+G* basis set by referring to the geometry of the precedent cluster (C6H6)n−1. The lowest-energy geometry of the tetramer takes the cyclic form consisting of perpendicular dimer geometries. The other optimized geometries of (C6H6)4 include no tetrahedral shape; the CMMs in the OPLS-AA geometry construct the tetrahedron as shown in Figure 3. The trigonal bipyramidal geometry is shown for the CMMs in the OPLS-AA geometry of the n = 5 cluster, but it is not taken into account in the optimized geometries of the pentamer in ref 17. In ref 17, the tetrahedral geometry is missed as a structural basis of the benzene clusters. As shown in Figure 3, the pentamer is built from tetramer by adding a molecule on the bottom facet of the tetrahedron and is also generated by subtracting the right-side molecule from the hexamer. Hence, the structural growth sequence of these clusters is continuous. The n = 6, 7, and 8 clusters lack continuity of geometries. The cluster with n = 9 is decomposed into (C6H6)8 and a benzene molecule, and similar decomposition is observed for the n = 9 and 10 clusters. The growth sequence of the n = 10−13 clusters is again discontinuous. The clusters with n = 13 and 14 are similar to each other if the right-side molecule in the n = 14 cluster is ignored. The above discussion is consistent with the results that small values of dn,n+1 (Figure 7) are obtained for n = 4, 5, 8, 9, and 13. The values of dn,n+1 suggest that the growth sequence of the clusters with n ≥ 14 is continuous with a few exceptions. This is confirmed by seeing the stereographic views of these configurations in Figure 3. The geometry of the cluster with n = 13 is prolate (see Figure 6), in disagreement with that of the LJ cluster with 13 atoms (spherical shape). The anisotropic intermolecular interaction must influence the geometry of the cluster as mentioned before. The prolate geometry of the n = 13 cluster develops to a different prolate one by adding 6 molecules. The intermediate clusters (n = 14−18) take prolate or nearly spherical shape. A prolate geometry is observed for each of the clusters with n = 19−22 where two interior molecules exist. The interior molecules form a prolate shape, and thereby, it is reasonable that outer molecules form a prolate shape. The clusters with n = 23−25 have three interior molecules as shown in Figure 3. Since a triangular-prism-like shape is formed by the interior molecules, outer molecules are considered to form an oblate (pancake) shape. An oblate geometry is observed for (CH4)23 as expected, but the shape of (C6H6)25 is nearly prolate. For the clusters with 26 ≤ n ≤ 28, the number of interior molecules is 4. The interior molecules construct oblate tetrahedral shapes similar to the structure of the tetramer. This may make the shapes of the clusters oblate. However, the oblate

Figure 9. Relative stability Sn (eq 7) of the global-minimum geometries expressed by the OPLS-AA (closed circles) and Williams−Starr (open circles) potentials.

clusters with n = 13, 19, 23, 26, and 29 (magic numbers) are relatively stable in agreement with the magic numbers of the LJ atomic clusters whose potential energies are taken from ref 42. However, the WS clusters show no remarkable feature at the sizes of 26 and 29. The results obtained with the OPLS-AA potential is more reliable than those with the WS potential since the dimer geometry of the former model is much closer to the experimental geometry4−6 than the geometry of the latter model (see Figure 1). The number of the local structures with Nmax = 12 increases stepwise at the sizes of 13, 19, 23, and 26 (Figure 4). The n = 29 cluster takes a local structure with Nmax = 11 in addition to the local structures. The construction of these local structures leads to enhanced stability of the clusters. As well-known,43 the LJ atomic clusters take icosahedra as local structures; each of them consist of 13 nearest neighboring atoms. The local structures found in the benzene clusters are significantly deviated from regular icosahedra. For example, in the 13 molecule cluster distances between the central molecule and surrounding molecules range from 4.95 to 5.41 Å (these distances are equal to each other for a regular icosahedron). A regular icosahedron is built from 20 slightly distorted tetrahedra.43 The benzene tetramer takes a considerably distorted tetrahedron because of the anisotropy of the intermolecular potential; r12 = r13 = r24 = r34 = 4.97 Å and r23 = r14 = 5.64 Å. This makes the local structures of 13 nearest neighboring molecules in the benzene clusters irregular. Hence, the anisotropy causes the distorted structural motif in the benzene clusters. In the previous study on the ethylene clusters expressed by the Morse potential,44 they show no icosahedral motif because of the anisotropy of the intermolecular potential. Anisotropy of intermolecular potential plays an important role in determining structural features of molecular clusters, and its effect can be monitored by examining the geometry of the tetramer. 3.2. Comparison of Williams−Starr Geometries with All-Atom OPLS Ones. Figure 8 shows the distances between the global-minimum WS and OPLS-AA clusters. For the cluster with the distances of approximately 0.01, the similarity between the OPLS-AA and WS geometries could be confirmed by seeing stereographic views of these configurations. The large distances (n = 8, 10, 21, 24, 25, and 27−30) indicate that the globalminimum WS geometries are different from the global-minimum OPLS-AA ones. Large geometrical differences observed for n = 24, 25, and 27−30 suggest that the error introduced by the WS potential would be enlarged by the increasing size. F



Article

shape is observed for only the n = 26 cluster. Hence, the shapes of the clusters do not always depend on those formed by the interior molecules. For the methane clusters,45 the geometries of the clusters are closely related to the shapes constructed by interior molecules. This is not valid for the benzene clusters as found for n = 24, 25, 27, and 28. Hence, the anisotropy of the intermolecular interaction in the benzene clusters must affect these geometries. As discussed before, the icosahedron consisting of 13 molecules is considered as a structural motif, and its geometry is distorted because of the anisotropy of the intermolecular potential. The stepwise construction of the icosahedron stabilizes the clusters with n = 13, 19, 23, 26, and 29. These numbers are equal to the magic ones of the LJ clusters. The anisotropy of the intermolecular potential between benzene molecules is not strong so that the above magic numbers appear. The construction of the icosahedron is a building-up principle of the benzene clusters and results in the continuity found for the growth sequence of the benzene clusters. The anisotropy affects the geometries of two nearest-neighboring molecules as described later. 3.4. Comparison of the Clusters with the Liquid. The experimental results on the liquid32 show that the number of the nearest neighbor molecules around a molecule is approximately 12 in agreement with the results of the MD simulations.26 The present study indicates that the maximum number of molecules surrounding a molecule in the clusters is 12 (Figure 4). The crystallographic geometry of solid benzene9 shows that the number of nearest neighboring molecules is 12. Hence, the maximum number of the nearest neighbor molecules is independent of the states of molecular aggregations. Figure 5 shows that, in the n = 2, 3, and 4 clusters, two nearest neighboring molecules take perpendicular geometries. The population of the perpendicular configuration in the clusters (>75%) is larger than that in the liquid (57%) obtained by the MD simulation of liquid benzene26 with the OPLS-AA potential. The population of the parallel configuration in the clusters is smaller than the corresponding value of the liquid 24% except for the n = 10 cluster (25%). In both the clusters and the liquid, two nearest neighboring molecules prefer the perpendicular configuration to the parallel one. To further examine geometrical features, the number of the configurations of the molecular pairs is plotted as a function of θ in Figure 10. The numbers of the molecular pairs with θ = 59°

shows that these are parallel displaced ones. The relative orientations of two molecules were investigated for all the configurations, and it was found that no face-to-face configuration is observed in the benzene clusters. The neutron diffraction study32 shows no face-to-face configuration in the liquid and reports that the most probable angle of θ of the parallel configuration is around 20°. The value is slightly different from the value of 28° obtained for the parallel configuration in the present study. It is difficult to explain the origin of the difference since the results are based on different backgrounds, the experiment on the liquid and theoretical calculation on the clusters. The configurational distribution experimentally obtained for the liquid32 is close to be isotropic, and thereby, the maximal distribution is observed at θ = 90° (see Figure 4d in ref 26 and Figure 10a in ref 32). However, the configurational distribution obtained for the clusters is considerably anisotropic as shown in Figure 10. The maximum is observed for θ ≅ 60°, and this value is close to the value of the dimer (67°). This suggests that, in the small benzene clusters under the investigation, molecular pairs tend to take configurations similar to the dimer. In aggregations of benzene molecules, a configuration of two nearest neighboring molecules is perturbed by aromatic− aromatic interactions between them and surrounding molecules. Consequently, the distribution of θ for the clusters must be different from that for the dimer. In the dimer, the global minimum (Y-shaped), T-shaped, parallel displaced, and face-toface configurations take energies of −9.69, −9.01, −8.79,26 and −7.09 kJ mol−1, respectively. The last three configurations do not correspond to energy minima on the potential energy surface. The perturbation would enable the T-shaped and parallel displaced configurations to have considerable populations. It is considered that the perturbation is not so strong as the highestenergy face-to-face configuration appears. In the liquid, the perturbation is more significant, and thereby, the distribution is nearly isotropic. 3.5. Prediction of Structures of Benzene Clusters with n = 13−19. Easter and co-workers11 measured the R2PI spectra of the isotopically substituted benzene clusters (C6H6)(C6D6)n−1. The transitions due to the C6H6 molecule at the cluster interior site were observed, and the numbers of the configurations were determined to be one to three from the transitions. The present author located a lot of local energy minima as well as global minima on the potential energy surfaces of the clusters; selected configurations of the n = 13−15, 17, and 18 clusters used later are shown in Figure 11. The most stable three configurations of the 13 molecule cluster in Figure 11 are consistent with those reported by Easter.15 The orientations of three lower molecules in the red rectangle in the second stable configuration (minimum 2) are different from those in the most stable configuration (minimum 1). According to the NEB method,41 the saddle point between the two minima takes the energy higher than the global minimum by 6.8 kJ mol−1; the geometry calculated for the saddle point (see Figure 11) indicates that the minimum 2 is converted into the minimum 1 by the concerted movements of the four molecules in the red rectangle. These results indicate that the two energy-minimum configurations are separately detected in the experiment. The third stable configuration (minimum 3) is similar to the minimum 2; a slight difference is detected for the orientation of the left-side molecule in the green rectangle. The geometry for the saddle point between the minima 2 and 3 is the intermediate between the geometries for these minima. The barrier height between the minima 2 and 3 is calculated to be 0.2 kJ mol−1;

Figure 10. θ dependence of the numbers (N) of the molecular pairs with intermolecular distances smaller than 6 Å in the global minimum configurations of the benzene clusters. The angle θ is calculated from the two vectors normal to the benzene rings.

and 62° are the largest, and a small peak at θ = 28° is observed. Careful examination of the configurations with the peak at θ = 28° G



Article

Figure 11. Stereographic views of the lowest-energy configurations and the configurations corresponding to some saddle points of the clusters with n = 13−15, 17, and 18. The global minimum geometries are also shown. The energies of the minima α relative to the lowest energy are denoted by ΔEα. The main differences between the two minima are displayed with red or green rectangles. The energies of the saddle points between the minima α and β relative to the lowest energy are denoted by ΔSαβ. The main geometrical differences between the minima and those between the saddle points and the related minima are indicated with the rectangles. For the blue rectangle and arrows in the minimum 1 of the n = 14 cluster, see text.

the value suggests that the two configurations can be easily interconverted. Easter15 reported that one of the two configurations detected in the experiment is the lowest-energy one and that the other is the minimum 2 or 3 without considering the transition states. The results of the present study, however, suggest that the configurations of the minima 2 and 3 cannot be separately identified because of the low barrier height between them. As the second stable configuration, the experiment must detect the configuration averaged over the large-amplitude motion including interconversion of the two forms. According to the interpretation of the experimental spectra observed for the 14 molecule cluster,11 three distinct clusters are detected, and each of the configurations takes a different location and/or orientation of the 14th molecule. Figure 11 shows that geometrical differences between the minima 1 and 2 are found for the position of the molecule in the small red rectangle and the orientations of the three molecules in the large rectangle. The relative energy of the saddle point for the transition from the minimum 1 to 2 (1-to-2 transition) is 9.4 kJ mol−1. The molecule in the right-side red rectangle in the minimum 1 moves to the different position according to the arrow shown in the figure, and the molecule in the blue rectangle simultaneously moves down according to the another arrow. Hence, the saddle-point geometry is considerably different from the geometry of the

minimum 1; the structural feature of the saddle-point geometry is eventually found for the four molecules in the red rectangle. Geometrical differences between the minima 1 and 3 are mainly ascribed to the orientations of the three lower molecules in the green rectangle. The relative energy of the saddle point for the 1-to-3 transition is 5.8 kJ mol−1 and its geometry shows the movement of the three molecules in the rectangle. According to the results for the two saddle points, the three energy-minimum configurations can be separately detected in the free jet expansion. This indicates that the three or four molecules are moved when the configurations are interconverted. This prediction is in disagreement with the speculation that three configurations take a different location and/or orientation of the 14th molecule.11 In the literature, the LJ geometries were assumed, while the present author carried out the structural analysis with the more realistic model. Consequently, the present prediction is considered to be more reliable than the previous one. Easter et al.11 also explained the configurations of other benzene clusters with n = 15 and 19 by referring to the geometry of the LJ clusters. The configurations proposed in ref 11 are not discussed hereafter. The spectra of the 15 molecule cluster are explained by the existence of two configurations.11 The configurations for the minima 1 and 2 are similar to each other except for the orientation of the upper molecule in the red rectangle. The energy of the saddle point between them is 0.5 and 0.1 kJ mol−1 higher than the H



Article

clusters are inconsistent with the core structures formed by the interior molecules. The results obtained for the liquid benzene are compared with the present data. For the molecular pairs with intermolecular distances smaller than 6 Å, the distribution of the configurations as a function of θ is anisotropic in the clusters, whereas that in the liquid is nearly isotropic. The configuration of two nearest neighboring molecules is perturbed by aromatic−aromatic interactions between them and surrounding molecules. The perturbation would enable the T-shaped and parallel displaced configurations to have considerable populations in the clusters. The configurations of (C6H6)n with n = 13−19 were experimentally investigated by Easter and co-workers.11 The present results suggest that the stable configurations detected in the experiment take different orientations for one to four nearestneighboring molecules. Drastic geometrical changes are generally accompanied with large energy differences, and thereby, the above geometrical changes would occur in the clusters generated in free jet expansions.

energies of them, respectively. Hence, the configuration 2 is easily convertible into the configuration 1 by the orientational motion of the benzene molecule in the red rectangle. The configuration 3 differs from the configuration 1 by the orientation of the top molecule in the green rectangle. The minimum energy path of the 3-to-1 transition shows a saddle point with the relative energy of 3.5 kJ mol−1. Accordingly, the cluster would exist in the mixture of the configurations 1 and 3. For the n = 16 cluster, one configuration predominantly exists.11 The lowest-energy configuration shown in Figure 3 is considered to be exclusively present. For the clusters with n = 17, 18, and 19, a pair of configurations are detected.11 The configuration 2 of the n = 17 cluster has different orientations for the two lower molecules in the rectangle compared with the configuration 1. The relative energy of the saddle point between the minima 1 and 2 is 2.4 kJ mol−1, and thereby, these configurations would be separately detected in the spectra. Figure 11 shows that the two molecules in the rectangle move independently of other molecules when the minimum 1 is changed to the minimum 2. The barrier height is smaller than the heights required for the movements of more molecules (ΔS12 for n = 13 and 14). The structures of the minima 1 and 2 of (C6H6)18 are similar to each other; the difference between them is found for the orientation of the right-side molecule in the red rectangle. The second stable configuration is easily converted into the most stable configuration since the barrier for the 2-to-1 transition is smaller than 0.1 kJ mol−1. The configuration 3 has the four molecules whose positions and/or orientations are considerably different from those in the configuration 1. These molecules are shown in the green rectangle for the minimum 3 and in the two green rectangles for the minimum 1. The saddle point between the two minima shows the movements of the four molecules with the relative energy of 12.6 kJ mol−1. Hence, the configuration 3 would coexist with the configuration 1. The 19 molecule cluster takes two interior molecules as shown in Figure 3. Since the two molecules are inequivalent, these molecules must be separately detected in the spectra. This conclusion was obtained for the WS cluster with n = 19.24 The experimental observation11 was examined with the energy minima and the barrier heights between them. The orientational change of one molecule generally gives no large energy increase, but the changes of two or more neighboring molecules result in the large barrier heights separating two nearest local minima as mentioned for the clusters with n = 13, 14, 17, and 18. The translational movement of a molecule is found for the n = 14 and 18 clusters (see the 1-to-2 transition for n = 14 and the 1-to-3 transition for n = 18). In these cases, the orientational changes of three molecules near the moved molecule are accompanied with the translational movement. Orientational degrees of freedom of a few molecules must be taken into account for estimating stable configurations of molecular clusters.

■

ASSOCIATED CONTENT

S Supporting Information *

Cartesian coordinates of the global-minimum geometries of (C6H6)n. This material is available free of charge via the Internet at http://pubs.acs.org.

■

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. Notes

The authors declare no competing financial interest.

■

REFERENCES

(1) Lee, E. C.; Kim, D.; Jurečka, P.; Tarakeshwar, P.; Hobza, P.; Kim, K. S. J. Phys. Chem. A 2007, 111, 3445−3457. (2) Kim, K. S.; Karthikeyan, S.; Singh, N. J. J. Chem. Theory Comput. 2011, 7, 3471−3477. (3) Gamieldien, M. R.; Strümpfer, J.; Naidoo, K. J. J. Phys. Chem. B 2012, 116, 324−331. (4) Kusaka, R.; Inokuchi, Y.; Ebata, T. J. Chem. Phys. 2012, 136, 044304. (5) Arunan, E.; Gutowsky, H. S. J. Chem. Phys. 1993, 98, 4294−4296. (6) Henson, B. F.; Hartland, G. V.; Venturo, V. A.; Felker, P. M. J. Chem. Phys. 1992, 97, 2189−2208. (7) Iimori, T.; Aoki, Y.; Ohshima, Y. J. Chem. Phys. 2002, 117, 3675− 3686. (8) Iimori., T.; Ohshima, Y. J. Chem. Phys. 2002, 117, 3656−3674. (9) Easter, D. C.; Whetten, R. L.; Wessel, J. E. J. Chem. Phys. 1991, 94, 3347−3354. (10) Easter, D. C.; Li, X.; Whetten, R. L. J. Chem. Phys. 1991, 95, 6362− 6370. (11) Easter, D. C.; Khoury, J. T.; Whetten, R. L. J. Chem. Phys. 1992, 97, 1675−1682. (12) Easter, D. C.; Baronavski, A. P.; Hawley, M. J. Chem. Phys. 1993, 99, 4942−4951. (13) Easter, D. C.; Mellott, J.; Weiss, T. J. Chem. Phys. 1998, 109, 8365−8373. (14) Easter, D. C. J. Phys. Chem. A 2003, 107, 2148−2159. (15) Easter, D. C. J. Phys. Chem. A 2003, 107, 7733−7742. (16) Chakrabarti, D.; Totton, T. S.; Kraft, M.; Wales, D. J. Phys. Chem. Chem. Phys. 2011, 13, 21362−21366. (17) (a) Mahadevi, A. S.; Rahalkar, A. P.; Gadre, S. R.; Sastry, G. N. J. Chem. Phys. 2010, 133, 164308. (b) Yeole, S. D.; Gadre, S. R. J. Chem. Phys. 2011, 134, 084111. (18) Dulles, F. J.; Bartell, L. S. J. Phys. Chem. 1995, 99, 17100−17106.

4. CONCLUSIONS The magic numbers of the clusters are found to be 13, 19, 23, 26, and 29, in good agreement with the results of the LJ atomic clusters. This is explained by considering the local structure constructed by 13 or 12 nearest-neighboring molecules as a building unit. The structural growth sequence of the globalminimum geometries is examined with the rotational constants. The anisotropy of the interaction between benzene molecules affects the geometries of benzene clusters; the local structures take distorted icosahedra, and the geometries of some of the I



Article

(19) Bartell, L. S.; Dulles, F. J. J. Phys. Chem. 1995, 99, 17107−17112. (20) Pullan, W. J. J. Chem. Inf. Comput. Sci. 1997, 37, 1189−1193. (21) White, R. P.; Niesse, J. A.; Mayne, H. R. J. Chem. Phys. 1998, 108, 2208−2218. (22) Williams, D. E. Acta Crystallogr., Sect. A: Found. Crystallogr. 1980, 36, 715−723. (23) van de Waal, B. W. Chem. Phys. Lett. 1986, 123, 69−72. (24) Takeuchi, H. J. Chem. Inf. Model. 2007, 47, 104−109. (25) Llanio-Trujillo, L. J.; Marques, M.; Pereira, F. B. J. Phys. Chem. A 2011, 115, 2130−2138. (26) Fu, C.-F.; Tian, S. X. J. Chem. Theory Comput. 2011, 7, 2240− 2252. (27) Jorgensen, W. L.; Severance, D. L. J. Am. Chem. Soc. 1990, 112, 4768−4774. (28) Sinnokrot, M. O.; Sherrill, C. D. J. Phys. Chem. A 2004, 108, 10200−10207. (29) DiStasio, R. A., Jr.; von Hilden, G.; Steele, R. P.; Head-Gordon, M. Chem. Phys. Lett. 2007, 437, 277−283. (30) Gonzalez, C.; Lim, E. C. J. Phys. Chem. A 2001, 105, 1904−1908. (31) Williams, D. E.; Starr, T. L. Comput. Chem. 1977, 1, 173−177. (32) Headen, T. F.; Howard, C. A.; Skipper, N. T.; Wilkinson, M. A.; Boron, D. T.; Super, A. K. J. Am. Chem. Soc. 2010, 132, 5735−5742. (33) Takeuchi, H. J. Comput. Chem. 2011, 32, 1345−1352. (34) Hartke, B. WIREs Comput. Mol. Sci. 2011, 1, 879−887. (35) Tao, Y.; Ruchu, X.; Wenqi, H. J. Chem. Inf. Model. 2011, 51, 572− 577. (36) Kazachenko, S.; Thakkar, A. J. Chem. Phys. Lett. 2009, 476, 120− 124. (37) Schönborn, S.; Goedecker, S.; Roy, S.; Oganov, R. J. Chem. Phys. 2009, 130, 144108. (38) Wu, X.; Cai, W.; Shao, X. Chem. Phys. 2009, 363, 72−77. (39) Liu, D. C.; Nocedal, J. Math. Prog. 1989, 45, 503−528. (40) Wilson, E. B., Jr.; Decius, J. C.; Cross, P. C. Molecular Vibrations; McGraw-Hill: New York, 1955. (41) Henkelman, G.; Jónsson, H. J. Chem. Phys. 2000, 113, 9978− 9985. (42) Wales, D. J.; Doye, J. K.; Dullweber, A.; Hodges, M. P.; Naumkin, F.; Calvo, F.; Hernandez-Rojas, J.; Middleton, T. F. The Cambridge Cluster Database. http://www-wales.ch.cam.ac.uk/CCD.html. (43) Wales, D. J. Energy Landscapes; Cambridge University Press: Cambridge, U.K., 2003. (44) Takeuchi, H. Comput. Theor. Chem. 2011, 970, 48−53. (45) Takeuchi, H. Comput. Theor. Chem. 2012, 986, 48−56.

J


Structural Features of Small Benzene Clusters (C6H6) - American

Recommend Documents