ARTICLE pubs.acs.org/JPCA
Noncovalent Interactions in Paired DNA Nucleobases Investigated by Terahertz Spectroscopy and Solid-State Density Functional Theory Matthew D. King,† Wayne Ouellette,‡ and Timothy M. Korter*,† †
Department of Chemistry, Syracuse University, 1-014 Center for Science and Technology, Syracuse, New York 13244-4100, United States ‡ Naval Surface Warfare Center, Indian Head Division, 4104 Evans Way, Suite 102, Indian Head, Maryland 20640, United States
bS Supporting Information ABSTRACT: Cocrystallized adenine and thymine derivatives, along with the pure monomeric crystals, were investigated by terahertz spectroscopy and solid-state density functional theory (DFT). The methylated nucleobase derivatives crystallize in planar hydrogen-bonded adenine-thymine pairs similar to the manner found in DNA. The spectra obtained for 1-methylthymine, 9-methyladenine, and the 1:1 cocrystal in the range of 10100 cm1 clearly demonstrate that absorptions in this spectral range originate from the uniquely ordered assembly and the intermolecular interactions found in each individual crystal system. The quality of spectral reproduction for the DFT simulations of each system was clearly improved by the inclusion of an empirical correction term for London-type dispersion forces to the calculations. Notably, it was found that these weak dispersion forces in the adenine-thymine cocrystal were necessary to produce a properly converged crystal structure and meaningful simulation of the terahertz vibrational spectrum.
1. INTRODUCTION There is growing interest in the prospects of using terahertz spectroscopy in the advanced study of molecular dynamics and structural elucidation of large biological systems, such as proteins and long-chain deoxyribonucleic acid (DNA) molecules.14 The application of terahertz technology and the interpretations of spectra obtained from such systems is not yet able to accommodate this level of complexity, but as the field continues to expand with rigorous investigations of increasingly larger biosystems and their theoretical analyses, it may be possible to someday achieve such high aspirations. This potential has initiated the investigations into the terahertz spectra of many small biomolecules as it is thought that in order to understand the terahertz spectra of complex systems, one must first have a complete understanding of the vibrational motions exhibited by the constituent molecules.58 By logically starting with less convoluted systems, the intermolecular interactions giving rise to the absorptions in the terahertz region may be better understood. Understanding the specific physical origins of absorptions in the terahertz region is still in its infancy. It is well-known that these absorptions arise due to the intermolecular contacts of molecule in the solid state; however, uncovering the associated molecular motions and interaction energies is not a trivial task.9,10 Computational methods such as solid-state density functional theory (DFT) have shown a great deal of promise as reliable approaches to studying these low-energy lattice vibrational motions.1114 Solid-state DFT methods with periodic boundary r 2011 American Chemical Society
conditions can be quite computationally demanding, but they provide unparalleled results in overall computational performance in the reproduction of crystalline structures and terahertz spectra. A significant setback often encountered in the application of the majority of DFT methodologies is the inherent deficiency in the proper modeling of weak noncovalent interactions. This puts limitations on the breadth of problems that can be addressed by DFT, as experimental data (i.e., unit cell parameters) are often a prerequisite for beginning such calculations. This shortcoming can often be circumvented while still producing high-quality results, but using more robust computational methodologies is desirable.1517 Depending on the application, dispersion interactions may be critical stabilizing forces necessary to accurately reproduce the crystalline structure and vibrational frequencies, and therefore DFT methods should be supplemented with correction terms accounting for these interactions. Dispersion interactions are important in biological systems, influencing the manner in which proteins maintain their intricate three-dimensional structures and how nucleotide base pairs stack in the double helical conformation of DNA. As these forces are important in active biological molecules, it seems necessary to Special Issue: A: David W. Pratt Festschrift Received: December 14, 2010 Revised: March 1, 2011 Published: March 29, 2011 9467
dx.doi.org/10.1021/jp111878h | J. Phys. Chem. A 2011, 115, 9467–9478
The Journal of Physical Chemistry A evaluate the influences of these forces when studying smaller systems composed of their monomeric subunits. In this study, the noncovalent interactions between derivatives of the nucleobases adenine and thymine were examined by terahertz spectroscopy and solid-state DFT. In DNA, the adenine and thymine molecules are covalently bound to deoxyribose components of the DNA framework. For the derivatives of adenine and thymine studied here, a methyl group is substituted in the location of the usual glycosidic attachment. These nucleobases form hydrogen bonded pairs and are arranged in the cocrystalline solid in a similar configuration as found in their natural DNA environment. The pure monomers also form ordered crystals with unique hydrogen bonding patterns. Investigation of the noncovalent interactions in these crystals can be achieved by terahertz spectroscopy to directly access the low-frequency intermolecular vibrations that are along the nucleobase binding coordinates. While some terahertz studies have been performed on DNA components,8,1820 none has investigated these particular derivatives or the mimicked DNA base pairing between cocrystallized adenine and thymine. Also important to this study is the investigation into the influence of weak intermolecular forces in the stabilization of these nucleobases in the crystalline state. The terahertz spectra of 1-methylthymine (1MT), 9-methyladenine (9MA), and their 1:1 cocrystal (AT) were measured in the spectral range of 10100 cm1. Structural and normal mode calculations were performed on each system using solid-state DFT with periodic boundary conditions. Calculations performed both with and without corrections for London-type dispersion interactions were compared for their ability to reproduce the crystalline structures and terahertz spectra of these systems. The terahertz spectra of 1MT and 9MA pure crystals were well reproduced by both computational approaches. However, it was discovered that the AT cocrystallized system required the incorporation of dispersion forces in order to obtain a meaningful calculation of vibrational frequencies. Structural variations between the systems involving weak forces acting on the methyl rotors of the 1MT and 9MA molecules necessitated the inclusion of dispersion corrections to the DFT calculations for the proper stabilization of the methyl groups in the AT system, leading to the successful simulation of the terahertz spectrum. With the proper application of the dispersion corrections to DFT methods, the quality of the terahertz spectral simulations for each of the systems was significantly improved. This demonstrates that incorporating dispersion corrections in DFT methods provides a more meaningful and valid description of the physical interactions governing the crystalline structure and the lattice vibrations of molecular crystals.
2. EXPERIMENTAL SECTION Experimental spectra were recorded using a home-built timedomain pulsed terahertz spectrometer based on an amplified Ti: sapphire femtosecond laser system. ZnTe crystals were used for generation of terahertz radiation by optical rectification21 and detection by free-space electro-optic sampling.22 A detailed description of the terahertz spectrometer has been reported elsewhere.13 All chemicals were purchased from Santa Cruz Biotechnology, Inc. and used without further purification. Monomeric 1MT and 9MA crystals were grown by slow evaporation of individual saturated aqueous solutions at room temperature. The AT cocrystals were grown by evaporation of an aqueous solution containing a 1:1 molar ratio of 1MT and 9MA at room temperature.
ARTICLE
X-ray measurements were performed on a Bruker Kappa APEXII diffactometer at low temperature (100 K) using graphite-monochromated Mo KR radiation (λMoKR = 0.71073 Å) controlled using the APEX 2.0 software package.23,24 The data were corrected for Lorentz and polarization effects and absorption using SADABS.25 The structures were solved by direct methods. All non-hydrogen atoms were refined anistropically. After all of the non-hydrogen atoms were located, the model was refined against F2, initially using isotropic and later anistropic thermal displacement parameters. Hydrogen atoms were introduced in calculated positions and refined isotropically. Neutral atom scattering coefficients and anomalous dispersion corrections were taken from the International Tables, Vol. C. All calculations were performed using SHELXTL crystallographic software packages.26,27 Atomic positional parameters, full tables of bond lengths and angles, and anisotropic temperature factors are available in the Supporting Information. Samples for terahertz measurements were mixed with poly(tetrafluoroethylene) (PTFE) powder at concentrations of approximately 2% by mass for the individual 9MA and 1MT samples and 4% for the AT cocrystal. The samples were pulverized using a stainless-steel ball mill (Dentsply Rinn 3110-3A) to minimize particle size and reduce scattering of terahertz radiation.28 Approximately 0.55 g of the sample mixtures was pressed into pellets under a pressure of 2000 psi using a hydraulic press (ICL EZ-Press 12) equipped with a 13 mm stainless steel die, giving final pellet dimensions of 13 mm 2.2 mm. Pure PTFE pellets for use as “blank” references were prepared in the same manner. The samples and blanks for measurement were held in a cryostat (Janis Research Systems) equipped with 3 mm thick polymethylpentene windows. Data were acquired at 293 and 78 K with samples under vacuum. Samples and blanks were scanned 32 times, and data were averaged for each individual data set. A 32 ps scan window consisting of 3200 data points was used to capture the terahertz waveform, which was then symmetrically zero-padded to a total of 6000 data points for the data transforms. The resulting instrument resolution was approximately 1.0 cm1. The ratio of the power spectra obtained from the Fourier-transformed data sets of the sample and blank yields the terahertz absorption spectrum. Each terahertz spectrum presented in this work is the average of four individual terahertz spectra, each representing a complete set of sample and blank measurements.
3. THEORETICAL SECTION All DFT calculations were performed using the CRYSTAL09 software package.29,30 Calculations were performed using the PBE density functional31 with the atom-centered triple-ζ 6-311G(d, p) basis set.32 The PBE functional was selected for use in this study as this functional has been shown to consistently produce high-quality simulations of solid-state structures and terahertz spectra.13,15,33 Initial atomic positions and lattice dimensions used in the DFT calculations were taken from the X-ray crystallographic data. Calculations without corrections for dispersion forces were performed within fixed unit cell geometries. Fullgeometry optimizations were performed for calculations employing corrections for dispersion forces. In these full unit cell geometry optimizations, all lattice parameters were allowed to freely optimize within the constraints of preservation of space group symmetry. 9468
dx.doi.org/10.1021/jp111878h |J. Phys. Chem. A 2011, 115, 9467–9478
The Journal of Physical Chemistry A
ARTICLE
Table 1. Parameters Used for Dispersion Corrections in DFT Calculationsa C6b
atom
Atomic Dispersion Parameters Rvdwc (Grimme40)
Table 2. X-ray Crystallographic Data Obtained at 100 K for 1-Methylthymine, 9-Methyladenine, and the 1:1 Cocrystal 1-methylthymine
41
Rvdw (scaled )
9-methyladenine
1:1 cocrystal
empirical formula
C6H8N2O2
C6H7N5
C12H15N7O2
f.w. crystal system
140.14 monoclinic
149.17 monoclinic
289.31 monoclinic
H
0.14
1.001
1.301
C
1.75
1.452
1.525
N
1.23
1.397
1.467
space group
P21/c
P21/n
P21/m
1.409
a (Å)
7.1579(8)
7.4738(6)
8.2751(3)
b (Å)
12.0979(15)
12.2720(9)
6.3547(2)
c (Å)
7.5436(9)
7.6927(6)
12.8257(4)
R (deg)
90
90
90
β (deg)
91.731(7)
112.110(6)
106.818
γ (deg) V (Å3)
90 652.94(13)
90 653.68(9)
90 645.60(4)
Z
4
4
2
μ (mm1)
0.109
0.105
0.108
T (K)
100(2)
100(2)
100(2)
λ (Å)
0.71073
0.71073
0.71073
R1a
0.0530
0.0448
0.0405
wR2b
0.1302
0.1162
0.1234
O
0.70
1.342 s6 Values (PBE Functional)
Grimme40
0.75
d
optimized 1MT
0.59
9MA
0.55
AT
0.57
a
b
Additional atomic parameters and s6 values can be found in ref 40. In J nm6 mol1. c In angstroms. d This work.
The sampling rate as a function of k points used in defining the density matrix and commensurate grid in reciprocal space were determined for total energy convergence criteria of ΔE < 108 hartree for geometry optimizations and ΔE < 1011 hartree for normal mode calculations.34,35 In order to meet these convergence criteria, shrinking factors of 4 were used for 1MT and 9MA calculations and 5 for the AT calculations. A total of 446773, 461220, and 901860 grid points were used for calculations with fixed unit cell parameters for the 1MT, 9MA, and AT systems, respectively. For calculations using full-geometry optimizations, a total of 447964, 462209, and 903996 grid points were used for the 1MT, 9MA, and AT systems, respectively. The radial and angular distributions of points were defined by a pruned (75,974) integration grid for each calculation. Truncation tolerances used for Coulomb and HF exchange integral series were set to 108, 108, 108, 108, and 1016 hartree.36 Frequencies of normal modes were calculated within the harmonic approximation by numerical differentiation of the analytical gradient of the potential energy with respect to atomic position.37 The IR intensities for normal modes were calculated from the dipole moment derivatives (dμ/dQ) determined by the Berry phase approach of calculating Born charges as polarization differences between equilibrium and distorted geometries.38 Calculations corrected for London-type dispersion forces employed the empirical energy correction term originally proposed by Grimme.39,40 The dispersion correction is of the form
R1 = Σ|Fo| |Fc|/Σ|F0|. b wR2 = Σ[w(Fo2 Fc2)2]/ Σ[w(Fo2)2]1/2; w = 1/[σ2(Fo2) þ (a*P)2 þ b*P]; P = [max(Fo2,0) þ 2Fc2]/3. a
fdmp ðRij, g Þ ¼ ð1 þ edððRij, g =Rvdw Þ 1Þ Þ1
ð2Þ
where Rvdw is the sum of the atomic van der Waals radii and d determines the steepness of the dampening function. A value of d = 20 was used in all calculations. The C6 and Rvdw values were taken from Grimme.40 The Rvdw were scaled by 1.05 for heavy atoms and by 1.30 hydrogen atoms as suggested by Civellari et al.41 However, these dispersion parameters have been shown to be too aggressive for the application in the simulation of crystal structures at nonzero temperatures and their associated terahertz spectra.42 To compensate for the overestimation of dispersive energies, the global s6 scaling factor was optimized independently for each system. The s6 value was optimized by modifying the scalar to minimize the change in unit cell dimensions with respect to the experimental dimensions supplied by the cryogenic X-ray crystallographic measurements. The s6 values and atomic parameters used for the dispersion-corrected DFT calculations are summarized in Table 1.
4. RESULTS AND DISCUSSION Edisp ¼ s6
N1
N
∑ ∑ ∑
i¼1 j¼i þ 1 g
ij
C6 fdmp ðRij, g Þ Rij6, g
ð1Þ
where the energy is the summation over all atom pairs and g lattice vectors, Cij6 is the dispersion coefficient for atom pair ij, Rij,g is the interatomic distance between atoms i and j, and S6 is the functional-dependent scaling factor. Because these calculations are performed on periodic systems, the dispersion energies are summed over all atoms i in the reference unit cell and atoms j in the unit cell at distance |g|, excluding the self-interactions of i = j contributions for g = 0. A cutoff distance of 25 Å was used to truncate the summation over lattice vectors. A dampening function is included in the dispersion correction term to avoid near-singularities for small values of R, which is given by
4.1. Experimental Section. 4.1.1. X-ray Crystallography. The crystal structures for 1MT, 9MA, and AT at 100 K were measured by X-ray crystallography. The results of these measurements are summarized in Table 2. The molecular structures and atomic labeling scheme adopted for molecular descriptions are shown in Figure 1. The 1MT system was found to be monoclinic with the space group P21/c. The unit cell contains four molecules (Z = 4), with one symmetry-unique molecule per unit cell (Figure 2a). The 9MA structure was found to be monoclinic with the space group P21/n, containing four molecules per unit cell (Figure 2b). The AT cocrystal was also monoclinic with space group P21/m. The unit cell contains four molecules, consisting of two each of the 1MT and 9MA molecules giving a Z = 2 unit cell (Figure 2c). The asymmetric unit contains one each of the 1MT and 9MA 9469
dx.doi.org/10.1021/jp111878h |J. Phys. Chem. A 2011, 115, 9467–9478
The Journal of Physical Chemistry A
ARTICLE
Figure 1. Atomic labeling scheme for (a) 1-methylthymine and (b) 9-methyladenine.
molecules. As this structure was solved having P21/m symmetry, the AT pair of the asymmetric unit is constrained to a planar geometry, including the methyl groups on both molecules being restricted to a mirrored geometry about the molecular plane. The general unit cell parameters and atomic positions of each system are in agreement with previously reported higher temperature structures.4347 The crystal structures reported here show minor constrictions in unit cell volumes as expected with the decreased thermal motion of the molecules at reduced temperatures. (CCDC’s 804320, 804321, 804322 contain supplemental crystallographic data for 1-methylthymine, 9-methyladenine, and the 1:1 cocrystal, respectively. These data can be obtained free of charge from the Cambridge Crystallographic Data Centre via http://www.ccdc.cam.ac.uk/data_request/cif.) 4.1.2. Terahertz Spectroscopy. The experimental terahertz spectra of 1MT, 9MA, and AT from 10 to 100 cm1 were obtained at 293 and 78 K (Figure 3). The vibrational frequencies for the absorptions observed in each spectrum are listed in Table 3. For the spectra of each system, broad absorption features are observed at room temperature. When cooled to liquid-nitrogen temperatures, the spectral features become narrower and better defined, uncovering additional underlying modes due to the reduction in baseline and peak overlap. This observation is due to the reduction in thermally populated vibrational states at lower temperature. The features were also observed to increase in energy at the lower temperature as a result of the constriction of unit cell volume upon cooling. The 1MT spectrum exhibited four absorptions at room temperature centered at 55.8, 69.9, 80.1, and 88.5 cm1 (Figure 3a). An underlying absorption was evident on the low-energy side of the peak located at 69.9 cm1 but was not able to be resolved at this temperature. When the sample was cooled to 78 K, the two modes constituting this absorption feature were resolvable at 68.9 and 72.5 cm1. The other peaks were shifted to higher energies, including the peak at 88.5 cm1 which was shifted beyond the experimental bandwidth. The barely visible feature at 80.1 cm1 in the 293 K spectrum became much more defined at 87.1 cm1 in the 78 K spectrum. The 9MA spectra obtained at 293 and 78 K each consisted of three absorptions (Figure 3b). The peaks located at 43.3, 81.1, and 87.9 cm1 in the 293 K spectrum were blue-shifted to 46.5, 87.4, and 95.4 cm1 in the 78 K spectrum. The observed intensities of these absorptions were greatly increased at lower temperature due to the reduction in peak width while conserving the overall integrated intensities. The cocrystallized 1MT and 9MA sample showed a markedly different terahertz spectrum than those of the pure constituent
monomers (Figure 3c). This clearly demonstrates that the terahertz spectrum is a product of the intermolecular interactions and not of the intramolecular structure of the cocrystallized molecules. Although some of the features were nearly coincident, the AT spectrum was distinctly different from either the 1MT or 9MA spectra. The terahertz spectrum obtained at 293 K consists of three features located at 43.8, 69.7, and 92.0 cm1. These modes were shifted to 47.1, 73.9, and 93.5 cm1 in the 78 K spectrum. An additional vibrational mode was evident on the low-energy side of the peak at 93.5 cm1, although the precise vibrational frequency of this mode could not be determined. As observed in the other terahertz spectra, the absorption intensities were increased at lower temperature due to the reduction in peak width. 4.2. Theoretical Section. 4.2.1. Structural Analysis. Thorough analysis of structures calculated by DFT is necessary to ensure the high-quality reproduction of crystalline structures and vibrational spectra. Calculated structures were evaluated and reported in terms of root-mean-squared deviations (RMSDs) of the calculated bond lengths and bond angles compared to those determined experimentally by X-ray crystallographic measurements. The quality of reproduction of intermolecular contacts was assessed by comparing hydrogen bond heavy atom separations, as well as non-hydrogen bonding short-contact distances as necessary, between the calculated and experimental structures. Hydrogen bonding parameters were evaluated only in terms of heavy atom separations as hydrogen positions are not accurately measured by X-ray crystallography. Calculations were performed both in fixed unit cell dimensions, using those determined by X-ray measurements, and in fully optimized unit cell dimensions in which all lattice parameters were allowed to freely optimize while preserving space group symmetries. The calculations performed without dispersion corrections are denoted as PBE, while those incorporating the dispersion corrections are denoted as PBE-D. Complete lists of heavy atom bond lengths and bond angles for the experimental and calculated structures are available as Supporting Information for this paper. The AT unit cell optimizations and normal mode calculations were performed in the P21 space group symmetry as opposed to the P21/m space group determined by the X-ray measurements. The P21/m space group constrains the 1MT and 9MA molecules of the asymmetric unit to a planar geometry. This higher symmetry does not allow for deviations from planarity and also restricts the positioning of the methyl side groups to a mirrored geometry about this molecular plane. The heavy atom positions measured by X-ray are the vibrationally averaged positions and 9470
dx.doi.org/10.1021/jp111878h |J. Phys. Chem. A 2011, 115, 9467–9478
The Journal of Physical Chemistry A
ARTICLE
Figure 3. Experimental terahertz spectra obtained at 293 K (dashed) and 78 K (solid) for (a) 1-methylthymine, (b) 9-methyladenine, and (c) the 1:1 AT cocrystal.
Figure 2. Representations of unit cell packing arrangements and hydrogen bonding for (a) 1-methylthymine, (b) 9-methyladenine, and (c) the 1:1 AT cocrystal.
are therefore correct as the experimentally observed atomic positions. In addition, the hydrogen positions are not accurately determined by X-ray crystallography and are mathematically placed in the crystal structure. Due to the low energy barriers for the rotation of these methyl groups, the minimum energy configuration likely deviates slightly from the restricted P21/m geometry. Relaxing the symmetry to P21 allows for free rotation of the methyl groups and the 9MA amine group to the true minimum energy and allows the molecules to optimize out of the exactly planar geometry without breaking symmetry. PBE calculations performed in both P21/m and P21 symmetries produced a
Table 3. Experimental Terahertz Frequencies (10100 cm1) for 1-Methylthymine, 9-Methyladenine, and the 1:1 Cocrystal 1-methylthymine
a
9-methyladenine
1:1 cocrystal
293 K
78 K
293 K
78 K
293 K
78 K
55.8a
59.5
43.3
46.5
43.8
47.1
69.9
68.9
81.1
87.4
69.7
73.9
80.1 88.5
72.5 87.1
87.9
95.4
92.0
93.5
Frequencies listed in cm1.
P21 structure that was approximately 1.11 103 kJ mol1 (4.23 107 hartree) lower in energy. While this is a small energy difference, it is greater than the convergence criteria of 108 9471
dx.doi.org/10.1021/jp111878h |J. Phys. Chem. A 2011, 115, 9467–9478
The Journal of Physical Chemistry A
ARTICLE
hartree used in the geometry optimizations. Likewise, the same calculations using PBE-D produced a P21 structure considerably lower in energy by 9.47 kJ mol1 (3.61 103 hartree). This clearly shows that optimizations restricted to P21/m symmetry would not be at a true energy minimum, and this large change in energy is in part due to the positioning of the methyl groups and associated CH 3 3 3 X interactions. The difference of nearly 9.5 kJ mol1 is a considerable change in overall energy, but when taking into account the two asymmetric units of two molecules each within the unit cell, this gives a reasonable maximum stabilization energy per methyl rotor by long-range London forces to be approximately 1.6 kJ mol1. For the proper application of dispersion corrections to the DFT methods, the global s6 scalar was first optimized for each system by modifying the value in order to minimize the changes in lattice dimensions upon full-geometry optimizations of the unit cells as compared to the 100 K X-ray data. These optimal values were determined to be 0.59, 0.55, and 0.57 for the 1MT, 9MA, and AT systems, respectively. These are considerably lower than the proposed default value of 0.75 suggested as the functional-dependent scalar for use with the PBE functional.40 However, a similarly reduced value of 0.52 was recently reported in the evaluation of the application of this empirical dispersion correction for the molecular crystal naproxen, where this approach was demonstrated to produced a high-quality simulation of the terahertz spectrum of this system.42 The percent changes in unit cell parameters for the fully optimized unit cells compared to the experimentally measured parameters are shown in Table 4. The resulting unit cell structures were quite similar to the experimental structures. The percent changes in unit cell volumes were largest for the 1MT and 9MA systems, which decreased in volume by only 0.23% and 0.24%, respectively, while the AT cocrystal experienced a negligible 0.04% decrease in volume. The PBE and PBE-D calculations produced similar RMSDs for both bond lengths and bond angles for each system. The RMSDs for bond lengths obtained from the PBE calculations were smaller overall than those from the PBE-D calculations (Table 5). For the 1MT molecules in the monomer crystal, the Table 4. Percent Changes in Unit Cell Dimensions for 1-Methylthymine, 9-Methyladenine, and the 1:1 Cocrystal for the Full Unit Cell Optimizations Using PBE-D in Comparison to the Experimental X-ray Data 1-methylthymine
9-methyladenine
1:1 cocrystal
a b
0.02 þ0.70
þ0.44 0.28
þ0.78 1.04
c
0.92
0.34
þ0.19
β
0.32
þ0.08
0.06
V
0.23
0.24
0.04
bond-length RMSDs were 0.0111 and 0.0118 Å for the PBE and PBE-D calculations, respectively. The 9MA calculations produced RMSDs of 0.0084 and 0.0085 Å by PBE and PBE-D, respectively. These values were improved in the calculations of the AT cocrystal, which had an overall bond-length RMSDs of 0.0075 and 0.0082 for the PBE and PBE-D calculations, respectively, considering both the 1MT and 9MA molecules. No overall trend was realized for the RMSDs for bond angles between the two methods, although the PBE-D calculations did produce lower RMSDs for all but the 1MT system. As observed with the bond lengths, the RMSDs for bond angles were decreased in the AT system for both the 1MT and 9MA molecules. The inclusive bond-angle RMSDs were 0.415 and 0.380° for the PBE and PBE-D calculations, respectively. The RMSD values obtained from all PBE and PBE-D calculations demonstrate the high-quality structural reproductions obtainable by DFT methods both with and without corrections for dispersion interactions. Previous studies dedicated to the study of DFT methods in the reproduction of the terahertz spectra of molecular crystals have shown that this level of structural reproduction is sufficient to accurately simulate terahertz spectra.5,11,33 In the P21 optimizations of the AT cocrystal, the molecules slightly deviated from the planarity forced by the measured P21/ m symmetry. The 1MT and 9MA molecules were moved out of plane with respect to one another by only 0.84 and 0.91° in the PBE and PBE-D calculations, respectively. The methyl and amine groups optimized similarly between the two methods. For the 1MT molecule, the C5 methyl group was rotated approximately 2.35° and the C6 methyl group 1.33° in the PBE calculation, while in the PBE-D optimized structure, these methyl groups were rotated 2.11° and 1.33°, respectively. The C6 methyl group of the 9MA molecule was rotated 1.56° in the PBE optimized structure and 1.74° in the PBE-D structure. The 9MA amine group was rotated about 0.25° in both the PBE and PBE-D calculations. While these structural modifications are small, they do influence the overall energy minimization of the AT system, allowing for more precise calculations of vibrational frequencies using these optimized structures. In addition to the assessment of intramolecular structure, it is also important to evaluate the reproduction of intermolecular contacts as absorptions in the terahertz region originate as a result of these interactions. The symmetry-unique hydrogen bonds of each system are listed in Table 6 in terms of heavy atom distances. Each 1MT molecule contains one unique intermolecular hydrogen bond with two total hydrogen bonds being made with a neighboring 1MT molecule forming a dimer (Figure 2a). These dimers are in turn held together in the solid by weaker dispersion interactions. The hydrogen-bond heavy atom separation was measured at 2.8174 Å. This distance was slightly underestimated in both the PBE and PBE-D calculations by 0.0407 and 0.0251 Å, respectively. The deviation was larger in
Table 5. RMSDs for Calculated Bond Lengths (Å) and Bond Angles (deg) for 1-Methylthymine and 9-Methyladenine in Pure Crystals and in the 1:1 AT Cocrystal monomers 1-MT
AT cocrystal 9-MA
1-MT
9-MA
PBE
PBE-D
PBE
PBE-D
PBE
bond lengths
0.0111
0.0118
0.0084
0.0085
0.0082
0.0089
0.0070
0.0075
0.0075
0.0082
bond angles
0.470
0.504
0.381
0.357
0.510
0.455
0.247
0.237
0.415
0.380
9472
PBE-D
PBE
overall PBE-D
PBE
PBE-D
dx.doi.org/10.1021/jp111878h |J. Phys. Chem. A 2011, 115, 9467–9478
The Journal of Physical Chemistry A
ARTICLE
Table 6. Hydrogen Bond Heavy Atom Distances (Å) for 1-Methylthymine, 9-Methyladenine, and the 1:1 Cocrystala hydrogen bond
exptlb
PBE
PBE-D
2.7767
2.7923
3.0307 2.9256
3.0113 2.9327
1-Methylthymine N1T(H) 3 3 3 O2T
2.8174 9-Methyladenine
N1A 3 3 3 (H)N5A N4A 3 3 3 (H)N5A
3.0500 2.9504
N1A 3 3 3 (H)N1T N5A(H) 3 3 3 O1T
2.9024
2.8865
2.8805
2.8200
2.7656
2.7758
N5A(H) 3 3 3 O2T
2.8680
2.8341
2.8789
A/T Cocrystal
a
Hydrogen positions are indicated in parentheses. b Experimental.
the PBE calculation possibly due to the overcompensation of these strong electrostatic interactions in the absence of dispersion forces which bind neighboring dimers in the crystalline solid. In the 9MA solid, the molecules form staggered hydrogenbonded chains that are stacked together by dispersion interactions (Figure 2b). Each molecule is engaged in four intermolecular hydrogen bonds, two of which are unique by symmetry. These bonds were measured at 3.0500 and 2.9504 Å, and were underestimated by an average of 0.0221 and 0.0282 Å for the PBE and PBE-D calculations, respectively. These hydrogen bond geometries are less affected by the presence or absence of dispersion interactions as there is a greater number per molecule than in the 1MT system, and they are arranged in nearly opposing directions so the overcompensation of the relative hydrogen bonding strength in the absence of dispersion is negated. Another contributing factor may be that the longer NH 3 3 3 N contacts of the 9MA system are weaker than the NH 3 3 3 O contacts made by the 1MT dimers. The hydrogen bonding scheme of the AT cocrystal involves staggered chains of 1MT molecules bound to 9MT molecules (Figure 2c). Each molecule forms three hydrogen bonds with two molecules of the other species. The N1A 3 3 3 HN1T and N5AH 3 3 3 O1T hydrogen bonds are in the same geometry as those found in adeninethymine base pairing in DNA. The third N5AH 3 3 3 O2T hydrogen bond links the base pairs in the bound chain. The N1A 3 3 3 HN1T and N5AH 3 3 3 O1T hydrogen bonds were measured at distances of 2.9024 and 2.8200 Å, respectively, and the N5AH 3 3 3 O2T at 2.8680 Å. These hydrogen bonds were all underestimated by the theory with the exception of the N5AH 3 3 3 O2T bond, which was overestimated by 0.0109 Å by the PBE-D calculation. The average deviations were 0.0347 and 0.0257 Å in the PBE and PBE-D calculations, respectively. To demonstrate the requirement of dispersion interactions to calculate accurate unit cell dimensions when performing fullgeometry optimizations, the unit cells were optimized without dispersion corrections. Very large increases in unit cell volumes were observed for all three systems. The AT unit cell underwent an expansion in volume of 14.9%, with similarly large increases for 1MT and 9MA of 12.4% and 13.3%, respectively. The resulting optimized structures are no longer physically meaningful representations of the unit cells and should therefore not be used in subsequent calculations of vibrational frequencies. Performing full unit cell optimizations without dispersion corrections are
useful, however, in determining the degree of basis set superposition error and estimating the amount of intermolecular interactions not accounted for in typical DFT calculations. Allowing the unit cells to relax along all dimensions will always lead to a lower energy optimized structure as the initial unit cell parameters induce strain on the unit cell contents due to the deficiencies in the representations of forces which is true of any computational approach. Examining the differences in energies between the calculated structures obtained from fixed and fully optimized unit cell parameters provides an estimation of these deficiencies. The decreases in unit cell energies associated with the volume changes observed for the full optimizations of the systems studied here ranged from 25 to 35 kJ mol1. With each system having four molecules per unit cell, this is approximately 69 kJ mol1 in molecular terms, which is a reasonable quantity for dispersive interaction energies for molecules in the solid state. It was demonstrated that corrections for dispersion interactions are required to reproduce the crystalline structure of the 1MT, 9MA, and AT systems when optimizing all unit cell parameters. The underestimation of dispersion interactions in DFT methods has been well documented, and thus the inclusion of dispersion corrections creates a more realistic representation of the physical forces governing molecular solid-state configurations when these corrections are appropriately applied. However, by removing the physical constraints of the fixed unit cell parameters, the molecules are permitted greater degrees of freedom, often leading to somewhat lower quality intramolecular structural reproduction as evident by the RMSD values. Nevertheless, the quality of the structures calculated using both the PBE and PBE-D methods are sufficient for use in subsequent normal mode calculations. 4.2.2. Simulated Terahertz Spectra and Mode Assignments. The quality of reproduction of the terahertz spectra by solid-state DFT methods can be assessed visually by comparison of the simulated 1MT, 9MA, and AT spectra with the experimental spectra (10100 cm1) obtained at 78 K (Figure 4). The correlation was quantitatively evaluated by calculating RMSDs between the frequencies of experimental peaks and the assignable calculated modes. Lorentzian line shapes were convolved into the calculated frequencies using values for the full width halfmaximum (fwhm) determined empirically through least-squares fitting of the experimental spectra. These fwhm values obtained from the 1MT, 9MA, and AT spectra were 4.1, 1.9, and 4.8 cm1, respectively. The calculated frequencies and intensities of the IRactive modes are listed in Table 7. IR intensities are reported in units of km mol1 but are normalized to the experimental intensities in Figure 4 for ease of visual comparison. The locations of low-intensity modes that may be present in the spectral simulations are not highlighted in the representations of the simulations, as these do not contribute significantly to the spectral absorption intensity; however, the frequencies and intensities of these modes can be found in Table 7. Calculated intensities listed as 0.00 km mol1 are of IR-active modes but with intensities below the level of precision presented here. Descriptions of the vibrational characters of modes that are assigned to experimental features are provided in Table 8. 4.2.2.1. PBE Simulations with Fixed Lattice Dimensions. The spectral simulations utilizing the PBE functional without dispersion corrections are shown in panels ac Figure 4. These calculations were performed using fixed unit cell dimensions obtained from X-ray measurements (Table 2). The simulated 9473
dx.doi.org/10.1021/jp111878h |J. Phys. Chem. A 2011, 115, 9467–9478
The Journal of Physical Chemistry A
ARTICLE
Figure 4. Simulated terahertz spectra calculated using PBE in fixed unit cell dimensions for (a) 1-methylthymine, (b) 9-methyladenine, and (c) the 1:1 AT cocrystal, and using PBE-D in fully optimized unit cell dimensions for (d) 1-methyladenine, (e) 9-methyladenine, and (f) the 1:1 AT cocrystal.
1MT spectrum is clearly assignable to the experimental spectrum (Figure 4a). The vibrational mode calculated at 58.8 cm1 was assigned to the lowest-energy experimental absorption centered at 59.5 cm1. The slight asymmetry of this experimental feature can be attributed to the underlying mode calculated at 54.4 cm1. The modes calculated at 72.0 and 76.1 cm1 were tentatively assigned to the absorptions at 68.9 and 72.5 cm1. Basing the mode assignments on only a single simulation, there is some ambiguity as to whether these two frequencies are equally overestimated in addition to the flawed intensities or if the mode calculated at 76.1 cm1 is alone overestimated and actually corresponds to the experimental absorption at 68.9 cm1. However, the tentative assignments were confirmed in the PBE-D calculation, which is discussed in more detail in the following section. The final mode in the experimental spectrum at 87.1 cm1 was logically assigned the remaining calculated mode at 91.2 cm1. The RMSD for frequencies was 3.2 cm1 taking into account the four assignable spectral features. The vibrational modes calculated for the 9MA spectrum by PBE could be tentatively assigned to the observed experimental
aborptions (Figure 4b). The only certain assignment is that of the lowest-energy absorption centered at 46.5 cm1 to the mode calculated at 42.6 cm1. The two larger absorptions could be assigned to the calculated modes at 90.1 and 93.7 cm1, although no definitive statement could be made without the aid of the additional PBE-D calculation. The relative intensities of the modes at 42.6 and 93.7 cm1 are greatly underestimated or the mode at 90.1 cm1 is overestimated by the theory for reasons that are unknown. Because the simulated spectrum is normalized to the latter, it gives the appearance as though the other two modes are underestimated. Despite the incorrect relative spectral intensities, the RMSD value for frequency location was quite good at 2.9 cm1. The frequency calculation for the AT cocrystal using the PBE functional was of poor quality, yielding imaginary frequencies involving motions of the methyl rotors. The imaginary values obtained for these methyl rotor displacements indicate that the optimized structure was not at an absolute minimum energy conformation, which stem from the ill-represented methyl positions in the absence of dispersion forces, as suggested in the 9474
dx.doi.org/10.1021/jp111878h |J. Phys. Chem. A 2011, 115, 9467–9478
The Journal of Physical Chemistry A
ARTICLE
Table 7. Calculated Vibrational Frequencies and Intensities for the IR-Active Modes of 1-Methylthymine, 9-Methyladenine, and the 1:1 Cocrystal below 100 cm1 1-methylthymine PBE freqa
9-methyladenine PBE-D
PBE
1:1 cocrystal PBE-D
PBE
PBE-D
intb
freq
int
freq
int
freq
int
freq
int
freq
int
39.6
0.03
39.3
0.00
42.6
0.07
44.3
0.19
16.7
1.03
14.3
0.00
54.4
2.02
51.1
2.90
74.6
0.02
70.2
0.02
28.9
0.98
22.7
0.04
58.8 72.0
8.35 15.96
58.5 68.8
8.12 13.78
88.7 90.1
0.00 7.24
86.5 87.3
0.13 6.01
43.1 44.6
0.02 0.00
45.6 46.8
0.02 0.02
76.1
11.88
72.6
11.91
93.7
4.56
98.3
3.40
91.2
5.85
90.1
4.28
45.5
0.00
47.6
0.00
50.5
1.48
50.3
0.01
57.4
0.01
51.5
1.65
68.8
0.24
56.7
0.05
74.0
7.71
70.1
0.05
74.0
7.98
76.7
7.35
78.2 86.4
0.23 25.07
79.3 84.9
0.10 1.53
86.8
0.08
90.8
33.74
97.3
10.36
98.3
0.04
RMSDc 3.2 a
1.6
2.9
2.1
—
7.3
Frequency. b Intensity. c RMSDs for calculated frequencies.
structural analysis. Without the incorporation of proper corrections for dispersion forces in the DFT methods, the structure of the AT cocrystal cannot be calculated. Because of the incorrect structure, the vibrational energies of all calculated modes were miscalculated, evident in the simulated terahertz spectrum (Figure 4c). It is important to note that this PBE calculation cannot be considered as a real spectral simulation since imaginary frequencies were obtained from the normal mode calculation. Although some calculated modes could be correlated with experimental absorptions, such as the two degenerate modes calculated at 74.0 cm1, it was unable to be determined whether these represent the actual vibrational character of the experimental absorption or if the locations of the coincident modes were purely by chance. Additional low-energy vibrational modes involving methyl rotations were calculated at 16.7 and 28.9 cm1 in a region of no experimental absorptions. Given the results of the normal mode calculation using the PBE functional without corrections for dispersion energy, it was concluded that no positive assignments could be made with confidence and that this spectral simulation must be disregarded. 4.2.2.2. PBE-D Simulations with Dispersion Corrections. The approach originally proposed by Grimme40 to supplement dispersion corrections to DFT methods was incorporated in the PBE-D simulations of 1MT, 9MA, and the AT cocrystal. To adjust the magnitude of the overall dispersion energy correction, the s6 parameter of eq 1 was modified independently for each system in the manner previously demonstrated for molecular crystals.42 This method was shown to increase the quality of the spectral simulations while using a more physically accurate model of the interaction energies influencing the terahertz spectra of molecular crystals. The terahertz spectral simulations performed using the PBE-D functional are shown in Figure 4df. The 1MT spectrum is readily assignable using the very high quality PBE-D simulation, as it was using the PBE functional
(Figure 4d). The PBE-D simulation verifies the mode assignments made in the PBE calculation, supporting the assignments of the two highest-intensity experimental features with the and ν1MT calculated at 68.8 and 72.6 cm1 (72.0 modes ν1MT 3 4 1 and 76.1 cm in the PBE calculation). The relative intensities of these modes are improved, although they remain inverse of the experimental intensities. The locations of ν1MT and ν1MT at 58.5 2 5 1 and 90.1 cm allowed for the likewise assignment of the experimental absorptions at 59.5 and 87.1 cm1, respectively. The vibrational mode ν1MT calculated at 51.1 cm1 is consider1 ably underestimated in frequency as this mode is likely an underlying mode of the absorption at 59.5 cm1. The RMSD of the calculated frequencies was improved in the PBE-D simulation to 1.6 cm1. This value, however, does not take into because the experimental account the underestimation of ν1MT 1 location of this mode could not be accurately determined. The vibrational characters of the modes lying in this spectral region are primarily of external character (Table 8). The assignment of the experimental absorptions in the 9MA spectrum became more apparent with the PBE-D simulation (Figure 4e). While the discrepancies in relative peak intensities remain, the mode locations are more accurately represented. The at 87.3 cm1 is nearly coincident with frequency of mode ν9MA 2 the experimental feature centered at 87.4 cm1. Very slight deviations between experimental and calculated frequencies exist and ν9MA calculated at 44.3 and 98.3 cm1. in the modes ν9MA 1 3 9MA was Mode ν1 was underestimated by 2.2 cm1, while ν9MA 3 overestimated by 2.9 cm1. The resulting RMSD for peak location was 2.1 cm1, which was an improvement of 0.8 cm1 over the PBE simulation. The vibrational motions observed in this spectral region of interest were mostly external in nature (Table 8). The most notable improvement resulting from the inclusion of dispersion corrections to the DFT calculations was observed in 9475
dx.doi.org/10.1021/jp111878h |J. Phys. Chem. A 2011, 115, 9467–9478
The Journal of Physical Chemistry A
ARTICLE
Table 8. Descriptions of Primary Mode Characters for Calculated IR-Active Vibrational Modes (1—100 cm1) Labeled in Figure 4 Based PBE-D Calculations 1-Methylthymine symmb
mode
freqa
ν1MT 1
51.1
Bu
ER about c axis
ν1MT 2 ν1MT 3
58.5 68.8
Au Bu
ET along a and c axes ER about a axis
ν1MT 4
72.6
Au
ER about c axis
ν1MT 5
90.1
Bu
ER about a and c axes
mode
freq
ν9MA 1 ν9MA 2 ν9MA 3
44.3
Au
ET along a axis
87.3
Bu
ER about b axis
98.3
Bu
ER about a and c axes
9-Methyladenine symm
mode descriptionc
mode description
1:1 Cocrystal mode description mode
freq
symm
νAT 1
51.5
B
ET along a axis; ER about b axis
ER about b axis
νAT 2
76.7
B
ER about b axis
ER about b axis
νAT 3
84.9
A
ER about a axis;
ER about a axis;
IM (methyl rotations)
IM (NH2 inversion)
νAT 4
90.8
A
1-methylthymine
ER about c axis
9-methyladenine
ER about c axis; IM (ring flex)
a Frequency, cm1. b Mode symmetry. c Key: ET = external translation; ER = external rotation; IM = internal motion.
the PBE-D simulation of the terahertz spectrum of the AT cocrystal (Figure 4f). While the calculation using the PBE functional without dispersion corrections produced a dismal spectral simulation for this system in that all attempts yielded imaginary frequencies, the PBE-D simulation generated a highquality terahertz spectrum by which all experimental absorptions could be unambiguously identified. A large number of IR-active vibrational modes were calculated, although only four modes are responsible for nearly all the observed absorption intensity in this spectral range (Table 7). The vast number of superfluous modes was due to the relaxation of symmetry constraints making modes that would typically be IR-inactive under P21/m symmetry, active in P21 symmetry, but with little or no spectral intensity. Of the 1 was assigned to experimental primary modes, νAT 1 at 51.5 cm 1 AT absorption at 47.1 cm , ν2 at 76.7 cm1 was assigned to the peak at 73.9 cm1, and νAT 4 was assigned to the highest-intensity experimental feature at 93.5 cm1. The underlying mode νAT 3 clearly contributes to the asymmetry observed on the low-energy side of the absorption at 93.5 cm1. The RMSD for frequencies was 7.6 cm1 considering the three resolvable absorptions. This value is larger than those obtained for the 1MT and 9MA systems; however, the experimental spectrum could be easily assigned by the calculated vibrational modes. The vibrational characters of these modes consist primarily of external translations and rotations, with the two higher-energy modes exhibiting some intramolecular motions (Table 8).
The success of the PBE-D calculation in the reproduction of the AT terahertz spectrum can be attributed to the stabilization of the methyl rotors through dispersion interactions. The short contacts made by the methyl groups of the 1MT and 9MA molecules in the AT cocrystal are stabilized by weak dispersion interactions that are not accounted for in the PBE calculations. The C5T methyl groups of the 1MT molecules form short contacts with the C5A hydrogen atoms of the 9MA molecules, and the C6T methyl groups of 1MT are located adjacent to the C6A methyl groups of 9MA. These methyl groups are also stabilized intramolecularly. While the equilibrium geometries of the PBE and PBE-D optimized structures are nearly identical, the behaviors of these methyl rotors during the vibrational motions are greatly affected by dispersive interactions. The artificially low energy barriers for the methyl rotations due to the absence of dispersion forces in the PBE calculations leads to the erroneous potential energy surfaces for the calculated vibrational motions. The inclusion of dispersion interactions stabilizes these methyl rotors and prevents the realization of false lower energy configurations during the dynamic motions of the lattice vibrations. These same incidences involving misrepresented methyl rotations were not observed in either the pure 1MT and 9MA systems as the methyl groups make stronger electrostatic interactions with neighboring oxygen and nitrogen atoms in these crystalline solids. The additional structural stabilization provided by the dispersion interactions is demonstrated by the differences in energy between the P21/m and P21 optimizations of the AT crystal structure. Without corrections for dispersion, the difference in energy between the P21/m and P21 optimized structures was on the order of 103 kJ mol1. When the dispersion corrections were included in the PBE-D optimizations, this difference was nearly 10 kJ mol1. The stabilization of the methyl rotors cannot account for the entire difference between the P21/m and P21 structures; however, the largest structural difference between these two crystal structures was the rotation of the 1MT and 9MA methyl groups out of their respective molecular planes. This suggests that this energy difference must in part be due to the positioning of the methyl side groups. The results obtained from the simulations of the 1MT, 9MA, and AT THz spectra provides convincing evidence for the importance of dispersion forces in the stabilization of molecular structures in the solid state even with the presence of strong hydrogen bonding interactions. Without the application of the dispersion corrections, accurate calculations of normal modes of vibration would not be possible for systems such as the AT cocrystal, where dispersion interactions are critical in the descriptions of the potential energy surfaces for the calculated lattice vibrations. Even in systems such as 1MT and 9MA, where spectral simulations were of adequate quality without inclusion of dispersion corrections, providing a better representation of the true physical forces governing solid-state interactions by incorporating these dispersion corrections yields more accurate calculation of vibration frequencies. This is evident by the lower RMSD values obtained for the PBE-D calculations of both the 1MT and 9MA systems. The results of this study demonstrate the elevated capacity to simulated crystalline structures and the terahertz spectra of molecular crystals through the proper application of modified dispersion corrections to conventional DFT methods. It is also important to note that accurate descriptions of all noncovalent interactions, including dispersion forces, are necessary for the calculation of structures and vibrational 9476
dx.doi.org/10.1021/jp111878h |J. Phys. Chem. A 2011, 115, 9467–9478
The Journal of Physical Chemistry A frequencies of cocrystallized nucleobase pairs. It is therefore reasonably concluded that expanded systems comprised of DNA components will also have this same computational requirement. Uncovering the vibrational nature of the terahertz absorptions of these simplified base pair analogues is important in understanding the more complex motions giving rise to absorptions in the terahertz spectra of larger DNA sequences. The exact energies of the concerted motions of the cocrystallized 1MT and 9MA molecules will differ from the adeninethymine base pair found in native DNA because of the differences in hydrogen bonding and weak noncovalent interactions. However, it has been undoubtedly shown that the near-rigid body vibrations of paired nucleobases lie in the frequency range readily accessible by terahertz spectroscopy. With the motions of the cocrystallized base pair being primarily rotational in character and modulated by the hydrogen bonds between the 1MT and 9MA molecules, it is expected that the same types of vibrational motions are likely observed at comparable frequencies between base pairs in natural DNA. This study serves as a benchmark in deciphering the intricacies of noncovalent interactions in DNA base pairing by terahertz spectroscopy and solid-state DFT. Further analysis of DNA base pair mimics and increasingly larger nucleotide fragments by these experimental and computational techniques could eventually lead to a comprehensive understanding of the sequence-specific terahertz spectra of DNA, which would be of great value to many areas of science, such as nondestructive biosensing and biomedical terahertz imaging.
5. CONCLUSIONS The terahertz spectra of 1-methylthymine, 9-methyladenine, and the 1:1 cocrystal were measured from 10 to 100 cm1. Solidstate DFT simulations of the crystal structures and the terahertz spectra were performed both with and without corrections for London-type dispersion interactions. The 1MT and 9MA spectra were well reproduced using the PBE and PBE-D functionals. However, the AT spectrum was unable to be adequately simulated in the PBE calculation excluding the dispersion corrections since the optimized structure yielded imaginary modes. When the correction term was incorporated in the DFT methods, a high-quality terahertz spectral simulation was obtained. The inability to calculate the lattice vibrational frequencies of the terahertz spectrum without corrections for dispersion interactions was attributed to the lack of stabilization of the methyl rotors of the 1MT and 9MA molecules in this particular crystalline system. This study demonstrates the importance of including proper dispersion energy corrections in DFT calculations of structures and vibrational frequencies of crystalline solids, where dispersion forces play a significant role in accurately describing the potential energy surfaces governing vibrational motions in the terahertz region. ’ ASSOCIATED CONTENT
bS
Supporting Information. Complete lists of experimental and calculated heavy-atom bond lengths and bond angles for these compounds. This material is available free of charge via the Internet at http://pubs.acs.org.
’ AUTHOR INFORMATION Corresponding Author
*E-mail:
[email protected].
ARTICLE
’ ACKNOWLEDGMENT This research was funded by a grant from the National Science Foundation CAREER Program (CHE-0847405). The authors also thank Syracuse University for continued support. ’ REFERENCES (1) Knab, J. R.; Chen, J.-Y.; He, Y.; Markelz, A. G. Proc. IEEE 2007, 95, 1605. (2) Globus, T. R.; Woolard, D. L.; Khromova, T.; Crowe, T. W.; Bykhovskaia, M.; Gelmont, B. L.; Hesler, J.; Samuels, A. C. J. Biol. Phys. 2003, 29, 89. (3) Markelz, A. G.; Roitberg, A.; Heilweil, E. J. Chem. Phys. Lett. 2000, 320, 42. (4) Kim, S. J.; Born, B.; Havenith, M.; Gruebele, M. Angew. Chem., Int. Ed. 2008, 47, 6486. (5) King, M. D.; Hakey, P. M.; Korter, T. M. J. Phys. Chem. A 2010, 114, 2945. (6) Korter, T. M.; Balu, R.; Campbell, M. B.; Beard, M. C.; Gregurick, S. K.; Heilweil, E. J. Chem. Phys. Lett. 2006, 418, 65. (7) Kutteruf, M. R.; Brown, C. M.; Iwaki, L. K.; Campbell, M. B.; Korter, T. M.; Heilweil, E. J. Chem. Phys. Lett. 2003, 375, 337. (8) Shen, Y. C.; Upadhya, P. C.; Linfield, E. H.; Davies, A. G. Appl. Phys. Lett. 2003, 82, 2350. (9) Beard, M. C.; Turner, G. M.; Schmuttenmaer, C. A. J. Phys. Chem. B 2002, 106, 7146. (10) Taday, P. F. Philos. Trans. R. Soc., A 2004, 362, 351. (11) Allis, D. G.; Fedor, A. M.; Korter, T. M.; Bjarnason, J. E.; Brown, E. R. Chem. Phys. Lett. 2007, 440, 203. (12) Allis, D. G.; Zeitler, J. A.; Taday, P. F.; Korter, T. M. Chem. Phys. Lett. 2008, 463, 84. (13) Hakey, P. M.; Allis, D. G.; Hudson, M. R.; Ouellette, W.; Korter, T. M. ChemPhysChem 2009, 10, 2434. (14) King, M. D.; Korter, T. M. J. Phys. Chem. A 2010, 114, 7127. (15) Allis, D. G.; Korter, T. M. ChemPhysChem 2006, 7, 2398. (16) Hakey, P. M.; Allis, D. G.; Ouellette, W.; Korter, T. M. J. Phys. Chem. A 2009, 113, 5119. (17) King, M. D.; Buchanan, W. D.; Korter, T. M. J. Phys. Chem. A 2010, 114, 9570. (18) Fischer, B. M.; Walther, M.; Jepsen, P. U. Phys. Med. Biol. 2002, 47, 3807. (19) Wilk, R.; Stewing, F.; Rutz, F.; Kleine-Ostmann, T.; Koch, M. Int. J. Nanotechnol. 2005, 2, 303. (20) Woolard, B.; Zhao, P. Proc. SPIE 2007, 6549, 65490H/1. (21) Rice, A.; Jin, Y.; Ma, X. F.; Zhang, X. C.; Bliss, D.; Larkin, J.; Alexander, M. Appl. Phys. Lett. 1994, 64, 1324. (22) Wu, Q.; Litz, M.; Zhang, X. C. Appl. Phys. Lett. 1996, 68, 2924. (23) APEX2, Version 2009.9; Bruker AXS, Inc.: Madison, WI, 2009. (24) SAINT, Version 7.68A; Bruker AXS, Inc.: Madison, WI, 2009. (25) SADABS: Program for Empirical Absorption Corrections, Version 2008/1; Sheldrick, G. M., Ed.; Bruker, AXS, Inc., Madison, WI, 2008. (26) XPREP, Version 2008/2; Sheldrick, G. M., Ed.; Bruker AXS, Inc.: Madison, WI, 2008. (27) Sheldrick, G. M. Acta Crystallogr., Sect. A 2008, A64, 112. (28) Shen, Y. C.; Taday, P. F.; Pepper, M. Appl. Phys. Lett. 2008, 92, 051103/1. (29) Dovesi, R.; Orlando, R.; Civalleri, B.; Roetti, C.; Saunders, V. R.; Zicovich-Wilson, C. M. Z. Kristallogr. 2005, 220, 571. (30) Dovesi, R.; Saunders, V. R.; Roetti, C.; Orlando, R.; Zicovich-Wilson, C. M.; Pascale, F.; Civalleri, B.; Doll, K.; Harrison, N. M.; Bush, I. J.; D’Arco, P.; Llunell, M. Crystal09 User’s Manual; University of Torino: Torino, 2009. (31) Perdew, J. P.; Burke, K.; Ernzerhof, M. Phys. Rev. Lett. 1996, 77, 3865. (32) Krishnan, R.; Binkley, J. S.; Seeger, R.; Pople, J. A. J. Chem. Phys. 1980, 72, 650. 9477
dx.doi.org/10.1021/jp111878h |J. Phys. Chem. A 2011, 115, 9467–9478
The Journal of Physical Chemistry A
ARTICLE
(33) Hakey, P. M.; Allis, D. G.; Hudson, M. R.; Korter, T. M. IEEE Sens. J. 2010, 10, 478. (34) Gilat, G. J. J. Comput. Phys. 1972, 10, 432. (35) Monkhorst, H. J.; Pack, J. D. Phys. Rev. 1976, 13, 5188. (36) Dovesi, R.; Pisani, C.; Roetti, C.; Saunders, V. R. Phys. Rev. B: Condens. Matter 1983, 28, 5781. (37) Pascale, F.; Zicovich-Wilson, C. M.; Gejo, F. L.; Civalleri, B.; Orlando, R.; Dovesi, R. J. Comput. Chem. 2004, 25, 888. (38) Dall’Olio, S.; Dovesi, R.; Resta, R. Phys. Rev. B: Condens. Matter 1997, 56, 10105. (39) Grimme, S. J. Comput. Chem. 2004, 25, 1463. (40) Grimme, S. J. Comput. Chem. 2006, 27, 1787. (41) Civalleri, B.; Zicovich-Wilson, C. M.; Valenzano, L.; Ugliengo, P. CrystEngComm 2008, 10, 1693. (42) King, M. D.; Buchanan, W. D.; Korter, T. M. Phys. Chem. Chem. Phys. 2011, 13, 4250. (43) Hoogsteen, K. Acta Crystallogr. 1959, 12, 822. (44) Hoogsteen, K. Acta Crystallogr. 1963, 16, 28. (45) Kistenmacher, T. J.; Rossi, M. Acta Crystallogr., Sect. B 1977, B33, 253. (46) McMullan, R. K.; Benci, P.; Craven, B. M. Acta Crystallogr., Sect. B 1980, B36, 1424. (47) Stewart, R. F.; Jensen, L. H. J. Chem. Phys. 1964, 40, 2071.
9478
dx.doi.org/10.1021/jp111878h |J. Phys. Chem. A 2011, 115, 9467–9478
