Microsolvation within the Systematic Molecular Fragmentation by

Dec 21, 2016 - The significance of these fragment sizes can be noted in Table 2, which shows the equivalent timings for an optimization step, run on 1...
1 downloads 12 Views 3MB Size
Article pubs.acs.org/JPCA

Microsolvation within the Systematic Molecular Fragmentation by Annihilation Approach Published as part of The Journal of Physical Chemistry virtual special issue “Mark S. Gordon Festschrift”. Rika Kobayashi,*,†,‡ Roger Amos,‡ and Michael A. Collins§ †

International Centre for Quantum and Molecular Structure, College of Sciences, Shanghai University, Shanghai 200444, China Australian National University, Leonard Huxley Bldg 56, Mills Road, Canberra, ACT 2601, Australia § Research School of Chemistry, Australian National University, Canberra, ACT 2601, Australia ‡

ABSTRACT: We have applied the systematic molecular fragmentation by annihilation (SMFA) fragmentation technique to glycine and DNA base pairs in water clusters, systems for which explicit solvation is believed to be important. The SMFA method was found to be capable of describing the structures, especially in handling the complexity of hydrogen bonding, with energies produced being comparable with those from full molecule results. Thus, the ability to break down large calculations into a manageable time without loss of accuracy shows promise for application to real biological systems for which these effects are relevant.

1. INTRODUCTION Solvation is known to play an important role in the interactions of many chemical systems, particularly biological ones. For example, Rueda et al. have shown that water is crucial to the stability of the DNA double helix1 and Nguyen et al. showed that highly structured water molecules mediate between ligands and the DNA minor groove.2 However, biological systems are often of a size too large to be handled by high-accuracy quantum chemical methods, thus giving rise to the popularity of fragmentation approaches, divide-and-conquer strategies to reduce the cost of the calculation. The “traditional” way of including solvent effects into ab initio quantum chemical calculations is by polarizable continuum solvent models (for a review, see Mennucci3) in which the solvent is modeled as a continuum reaction field with fixed dielectric. However, these reaction field methods are unable to pick up subtle solvation effects, arising from the effect of hydrogen bonding, thus requiring the inclusion of some measure of explicit solvation. It should be noted that the structural integrity of the aforementioned double helix is maintained by hydrogen bonds between its base pairs.4 Various fragmentation methods have been devised and applied to the properties of a variety of systems, including water clusters.5−8 An earlier paper by Pruitt et al.9 compared the fragment molecular orbital (FMO)10 and systematic molecular fragmentation (SMFA)11 methods to the binding energies of water clusters. They demonstrated that both methods could be viable and computationally efficient for water clusters despite their complex hydrogen-bonded networks. Collins further studied water clusters, particularly with respect to embedding charges, in a refinement of the SMFA method.8 However, little work has © 2016 American Chemical Society

been done applying fragmentation methods to solvation of ions and neutral molecules in water. Here we investigate the application of SMFA to microsolvation for various examples that demonstrate the effect of water molecules on the structure of model biological systems: glycine and DNA base pairs, to determine whether SMFA is precise enough to carry out geometry optimizations and energetics to capture these effects.

2. METHOD 2.1. Systematic Molecular Fragmentation by Annihilation. The basic rationale behind SMFA is that segments of a large molecule that are far apart, in terms of bonded connectivity, make almost independent contributions to the total electronic energy of the molecule. This idea is implemented in the following procedure. The atoms in a molecule (and solvent) system are collected into chemical functional groups. Atoms connected by double bonds are contained in the same group. Single bonds are defined in terms of covalent radii. Hydrogen bonds are defined between a hydrogen atom and an atom (A) from groups 5, 6, or 7 in the periodic table if the atom−atom distance RHA obeys RHA ≤ VdW(H) + VdW(A) − 0.32

(1)

where VdW denotes the van der Waals radius of an atom (in Ångstrom). For example, a hydrogen bond exists between two water molecules if an intermolecular H···O distance does not exceed 2.4 Å. A molecule and solvent is then represented by a Received: October 31, 2016 Revised: December 20, 2016 Published: December 21, 2016 334

DOI: 10.1021/acs.jpca.6b10919 J. Phys. Chem. A 2017, 121, 334−341

Article

The Journal of Physical Chemistry A

fragments are defined by the user on the basis of “chemical intuition”. This diversity leads, inter alia, to variations in accuracy and computational cost. A useful summary of the characteristics of several approaches is provided in Table 1 of the recent review.14 The table shows that the connectivitybased approach of SMFA leads to smaller fragment sizes than other methods (and higher computational efficiency), while maintaining similarly high accuracy compared to the generalized energy-based fragmentation (GEBF) method.8,15 Even higher accuracy can be obtained by methods such as molecular tailoring,16 but at the cost of much larger fragment sizes and computational cost. The total ab initio computational cost of fragmentation methods is only linearly proportional to the size of the system, as the number of fragments is linearly proportional to the size of the system. However, because all fragment calculations are independent of one another, the actual calculation time is only proportional to the size of the largest fragment, given sufficient processing units. Hence, SMFA provides a very suitable approach to a study of the solvation of large molecules. 2.2. Systems of Interest. The systems chosen were those that have already had their microsolvation studied thoroughly using conventional methods. This allowed comparison of our results with the conventional results and suggested refinements of the fragmentation approaches. First, we revisited zwitterionic glycine, followed by, as representatives of DNA, the guaninecytosine and adenine-thymine base pairs, for which, as mentioned above, hydrogen bonding is an important feature. Nonionized and Zwitterionic Glycine. Glycine is an amino acid that carries out many important functions in the human body, such as the biosynthesis of nucleic and other acids, breaking down of fat, production of heme, and as an inhibitory neurotransmitter in the central nervous system. In the gas phase the neutral form is most stable but in aqueous solution it is predominantly found as a zwitterionic tautomer. Gordon and others17−23 have investigated how many water molecules are necessary to stabilize the zwitterion tautomer. In the gas phase the zwitterion is a high-energy local minimum. Aikens and Gordon18 studied many structures using Hartree−Fock and MP2 level of theory with 6-311+G(d,p) basis. They found that, using MP2 and both explicit solvent molecules and a continuum

set of functional groups that are connected by single bonds or by hydrogen bonds. The system is then fragmented in a sequence of steps that are governed by a parameter denoted as Level, which takes integer values from 1 upward. Starting with some arbitrary group (call it group A): (i) Remove group A from the molecule, storing the resulting fragments. (ii) Leave A untouched, but remove all groups separated from A by more than a sequence of Level bonds (single bonds or hydrogen bonds), storing the resulting fragments. (iii) Remove A and all the groups in (ii) from the molecule, storing the resulting fragments. (iv) The molecule is represented by the sum of the fragments produced by steps (i) and (ii), minus the fragments produced by step (iii). (v) Take each fragment in (iv), repeat steps (i) to (iii), and continue with every new fragment produced until there is no fragment that has groups separated by more than Level bonds. For a general molecule and solvent system, M, the fragmentation into a total of Nfrag fragments at some value of Level is written as Nfrag

M→

∑ fn Fn n=1

(2)

where the number of fragments is denoted by Fn and f n is an integer coefficient (usually ±1). The number of fragments is proportional to the number of functional groups in the system. Where fragments are terminated by broken single bonds, hydrogen atom “caps” are added to restore the valence, as previously described.12 In the fragmentation procedure, hydrogen bonds are treated equivalently to covalent single bonds, except that hydrogen caps are not required when hydrogen bonds are broken. The electronic energy, E, of the original molecule and solvent system is given by the right-hand side of these fragment equations as Nfrag

E=

∑ fn E(Fn) n=1

(3)

Equation 3 neglects the energy of interaction of functional groups that are separated by more than Level bonds. If such groups are close in space, their interaction energy can be evaluated by using a “nonbonded” interaction scheme that employs a molecular fragmentation at Level = 1.13 More distant interactions can be treated by perturbation theory to account for electrostatic interactions, dispersion interactions, and induction effects.13 Alternatively, long-range electrostatic and induction interactions can be accurately accounted for using embedded charges in all ab initio calculations of fragment energies.8 By systematically increasing the value of the parameter Level, one obtains a systematic sequence of approximations to the total electronic energy of the molecule. Similar approximations to any molecular property (including energy gradients and Hessians) can be obtained from equations similar to eq 3. It is clear from this brief description of SMFA that fragmentation is based primarily on the bonded connectivity of the functional groups in a molecule, with a later account for “through space” interactions. This and several other current approaches to molecular fragmentation have been recently reviewed.14 Some other approaches also base the fragmentation mainly on bonded connectivity, some produce fragments that are mainly defined by collecting atoms that are close together in space. Some methods are not automated, in the sense that

Table 1. Number of Fragments and Range of Fragment Sizes from the Level = 1 SMFA Fragmentation of Neutral Glycine and Its Zwitterion N7-a Z7-a 8N1-a 8N8-a 8Z-a

no. of atoms

no. of fragments

fragment size

31 31 34 34 34

22 18 34 32 23

3−11 3−16 3−8 3−11 3−16

Table 2. Walltime (in hh:mm:ss) for the Whole Molecule Calculations Compared with Timings for the Largest of the Fragments N7-a Z7-a 8N1-a 8N8-a 8Z-a 335

walltime for whole molecule

walltime for largest fragment

08:26:05 09:22:52 12:13:52 13:03:51 13:55:21

00:05:34 00:31:28 00:02:13 00:05:41 00:31:05 DOI: 10.1021/acs.jpca.6b10919 J. Phys. Chem. A 2017, 121, 334−341

Article

The Journal of Physical Chemistry A

Figure 1. Optimized geometries from the SMFA method for the neutral glycine and zwitterionic structures: (i) N7-a, (ii) Z7-a, (iii) 8N1-a, (iv) 8N8-a, and (v) 8Z-a.

Adenine-Thymine Base Pair. Jacquemin et al.32 examined the adenine-thymine (AT) base pair to see whether there was a water catalytic effect in the double proton transfer that could lead to spontaneous mutation into the rarer tautomeric forms. This mechanism appeared to be strongly affected by the solvation environment and, with the inclusion of two water molecules in the reaction, they could find a stable double proton transfer product. 2.3. Computational Details. Coordinates for the systems mentioned above were fed into an automated geometry optimizer using the GAMESS-US33 program, version 18 AUG 2016 (R1), to obtain energies and gradients. The structures were analyzed and are displayed here with MacMolPlt.34 The molecules were fragmented under the SMFA procedure described above at Level = 2, and gradient and Hessian calculations were

solvent, seven or eight water molecules were required for zwitterionic glycine to become the global minimum. Bachrach19 found with the PBE024 density functional that seven water molecules were needed for the systems to become isoenergetic. Guanine-Cytosine Base Pair. The microsolvation of the guanine-cytosine (GC) base-pair system has been studied extensively by many groups.25−31 Kumar et al.25 investigated the change to hydration structure and electron affinity on solvation to understand radiation-induced damage, primarily through ionization, to the structure of DNA. They found that six to 11 water molecules were enough to have significant effects on the geometry of the base pair, in agreement with X-ray crystallographic data.26 In a later paper, Kumar and Sevila27 also looked at the possibility of mutations involving proton transfer between the bases without finding anything definite. 336

DOI: 10.1021/acs.jpca.6b10919 J. Phys. Chem. A 2017, 121, 334−341

Article

The Journal of Physical Chemistry A

impractical for a geometry optimization. Fragmentation at Level = 2 without including hydrogen bonds broke the system down into basically glycine and the individual water molecules, equivalent to a simple pairwise interaction model, which was deemed insufficient for getting the correct structures. Previous experience8 indicated that the inclusion of hydrogen bonding was important for obtaining the correct structures and so the fragmentation was carried out at Level = 1. This level of fragmentation produced a greater number of fragments but of much more reasonable fragment size, as shown in Table 1. The significance of these fragment sizes can be noted in Table 2, which shows the equivalent timings for an optimization step, run on 16 processors of an x86_64 cluster. Bearing in mind that as the fragments are independent and can be run simultaneously on a multicore machine, the time to solution for the whole calculation can be reduced to the time taken for the largest fragment. Note that in larger biological systems the time savings will be even more dramatic. These systems are probably too small to show the full benefits of the fragmentation method in that the fragments, especially when hydrogen bonds are considered, can almost be as large as the full system but were chosen as there were previous studies that could be compared with. For a real biological system of interest, which is our ultimate aim for this procedure, the fragments will be much smaller than the full system. The optimized geometries of the nonionized form, shown in Figure 1 i,iii,iv resemble those found as low-energy structures by Aikens and Gordon,18 though they are not identical as we have used different basis sets and have performed optimizations at the MP2 level. Aikens and Gordon optimized the structures at the SCF level and added MP2 corrections without reoptimizing. Our structures are also close to the full molecule geometry optimizations on the complete system. The zwitterion structures (Figure 1 ii,v) also resemble those from earlier papers, and likewise match the structures obtained by a full optimization on the whole cluster, except that the latter structures are slightly more compact with shorter hydrogen bonds. The zwitterions are local minima on the surface and retained their form throughout the geometry optimization. The images in the figure show that the waters are forming several hydrogen bonds with the amide group, including bridges that could facilitate a proton transfer back to the oxygen atoms, but no such transfer occurs spontaneously. Though SMFA at this level (Level = 1) was able to reproduce the explicitly solvated structures for these systems, the effect on the energetics was more subtle, and so further refinement of the method was necessary. To make the N and Z fragmentation calculations similar, the CN, CC, and CO bonds in glycine were defined to be double bonds. All the heavy atoms were assigned charges (usually zero). This was found to be very useful in the water cluster paper by the Collins and Gordon groups.9 This also means that no induction correction is required, so the group polarizabilities need not be calculated. No long-range dispersion corrections were required in any of the calculation (as the waters and glycine were too close in all structures). Larger clusters will no doubt require calculation of long-range dispersion interactions. Within SMFA, this is done using perturbation theory to calculate the dispersion coefficients between functional groups, using the imaginary frequency polarizability.13 Given results for water clusters as large as 64 molecules,9 it appears that long-range electrostatic interactions, including induction, are best treated using embedded charges,

Table 3. Total Energies, in hartree, of Various Structures of Nonionized and Zwitterionic Glycine with 7 or 8 Water Molecules Calculated for the Full System and for Level = 1 and Level = 2 Fragmentation Schemes and Relative Energies in kJ/mol of the Nonionized Form Compared to the Zwitteriona 8N1 8N8 8Z 8Z-8N8 N7 Z7 Z7-N7

full

level = 1

level = 2

−894.138684 −894.143132 −894.142267 2.3 −817.844095 −817.845989 −5.0

−894.136615 −894.141474 −894.140816 1.7 −817.842825 −817.845454 −6.9

−894.139239 −894.143241 −894.142983 2.8 −817.844145 −817.845972 −4.8

a

The notation for the isomers corresponds to Figure 1 and earlier papers.

Table 4. Number of Fragments and Largest Fragment Size from the Level = 2 SMFA Fragmentation of the DNA Base Pairs gu-cy.5water gu-cy.10water gu-cy.12water gu-cy.17water ad-thy.5water ad-thy.10water ad-thy.12water ad-thy.17water

no. of atoms

no. of fragments

largest fragment size

44 59 65 80 45 60 66 81

12 25 27 41 14 32 32 55

29 26 29 29 27 26 26 27

carried out on the fragments. The gradients and Hessians of the fragments are combined into values for the full molecule, and the geometry is updated via an augmented Hessian method.35 The glycine calculations were carried out at the MP2 level with 6-311++G(d,p) basis. The starting geometries were chosen to resemble those structures with seven or eight water molecules that were the lowest energy forms of nonionized glycine and the zwitterion, denoted N7-a, Z7-a, 8N1-a, 8N8-a, and 8Z-a, from previous papers.18,19 The base pairs were computed at the ωB97XD/6-311+ +G(2d,p) level starting with the optimized geometries of Kumar et al.25 and Jacquemin et al.32 but adding more solvent, eventually up to 17 initially randomly placed water molecules whose positions were optimized. At this number of water molecules, an almost complete ring around the base pair is formed. As these systems were larger, we used density functional theory for this part of the study. The ωB97XD36 functional was chosen as it includes a dispersion term, which is needed for the description of hydrogen bonds.37 Only selected structures have been chosen, as unlike previous papers, the intention here is not to search for a global minimum on the surface but to demonstrate that the SMFA fragment method is capable of accurate calculations on systems including microsolvation.

3. RESULTS AND DISCUSSION Fragmentation of all the systems was initially carried out at Level = 2, including the algorithm taking hydrogen bonding into account. However, for glycine, being a relatively small molecule, using Level = 2 produced several fragments that were almost equivalent to the full system and so glycine was deemed 337

DOI: 10.1021/acs.jpca.6b10919 J. Phys. Chem. A 2017, 121, 334−341

Article

The Journal of Physical Chemistry A

Figure 2. Optimized geometries from the SMFA method for the guanine-cytosine base pair with (a) 5, (b) 10, (c) 12, and (d) 17 explicit waters.

as herein. When calculated for the full cluster, the energies of the zwitterion forms were almost isoenergetic with the nonionized form; for the eight-water system the zwitterion was 0.0009 au (2.3 kJ/mol) higher, whereas for the seven-water structure it was 5.0 kJ/mol lower. Energy differences as small as this are quite challenging for a fragment approach, but the results are in agreement with the full molecule calculation, though depending on the level of the fragmentation. The eightwater zwitterion structure is 1.7 kJ/mol higher than the nonionized (8N8) form when Level = 1 fragmentation is used, and 2.8 kJ/mol higher with Level = 2. The seven-water system is 6.9 kJ/mol lower than the nonionized form with a Level = 1 fragmentation and 4.8 kJ/mol lower at Level = 2. The full energies of the structures are given in atomic units in Table 3. A larger case-study, and more representative of biological molecules of ultimate interest, is of the DNA base pairs guanine-cytosine (G-C) and adenine-thymine (A-T). There are many local minima on the potential energy surfaces of these

systems, and it was not the intention to do an exhaustive search for the global minimum, which would in any case not be a very meaningful objective as in practice these systems would be described as a thermal average over several local minima. Rather, the intention was to demonstrate that the SMFA approach has sufficient precision, for systems with weak hydrogen bonds, to obtain structures corresponding to the geometry optimization of the full G-C-water or A-T-water complexes. The fragmentation process again, at Level = 2, produced fragments (as shown in Table 4) that were basically the individual base pairs and water molecules, due to the necessity of keeping the rings intact. The optimized structures for the G-C and A-T complexes with 5, 10, 12, and 17 explicit waters, obtained via SMFA fragmentation at Level = 2, are given in Figures 2 and 3, respectively. These are larger numbers of water molecules than considered in previous studies. The structures obtained for both complexes are physically reasonable. For example, they all have direct hydrogen bonds 338

DOI: 10.1021/acs.jpca.6b10919 J. Phys. Chem. A 2017, 121, 334−341

Article

The Journal of Physical Chemistry A

Figure 3. Optimized geometries from the SMFA method for the adenine-thymine base pair with (a) 5, (b) 10, (c) 12, and (d) 17 explicit waters.

described. The hydrogen bonds are more variable. One first point is that the definition of an H-bond, as mentioned earlier, is that it is less than 2.40 Å. The images in the figures were produced using macMolPlt,34 which uses a different definition, a length of 1.95 Å. There are thus fewer H-bonds in the images than were included in the optimization. If the bond lengths are all explicitly calculated, then it is found that every bond specified as an H-bond in the fragments is also an H-bond in the full structure. There are strong H-bonds between the basepair dimers, and also where water molecules directly bond to one of the target molecules. These bonds are all reproduced by the fragment calculation, differing from the full molecule results by, on average, 0.02 Å. The full molecule optimization is producing slightly more compact structures. There are further H-bonds where water molecules are bound to other water molecules. These molecules are generally in about the right place; i.e., they form the H-bond networks seen in the figures, but their orientations can differ. Finally, there are a small

from water molecules to the unsaturated oxygen and nitrogen atoms in the G-C dimer, with the other water molecules forming a hydrogen-bonded network. One interesting feature, not obvious in the “face-on” views in Figures 2 and 3 is that the G-C structures remain essentially planar, but the A-T structures are bent, with the water molecules clustered on one side of the dimer. This is more obvious in the “side-on” views in Figure 4. This feature has been noticed previously.38 It is possible to quantify the statement that the geometries obtained from the fragment optimization are close to those obtained from a geometry optimization of the full molecule. If one takes the bond length of all covalent bonds in these systems, i.e., the bond within each molecule, including the water molecules, then the fragment optimized values have an RMSD (root-mean-square deviation) of 0.006 Å relative to the full optimization, and a MAD (median absolute deviation) of 0.004 Å. This includes the O−H bonds that are the donors in the hydrogen bonds, so these elongations are correctly 339

DOI: 10.1021/acs.jpca.6b10919 J. Phys. Chem. A 2017, 121, 334−341

Article

The Journal of Physical Chemistry A

Figure 4. Side view of the optimized geometries from the SMFA method for the (a) guanine-cytosine and (b) adenine-thymine base pairs with 17 explicit waters.

Notes

number of cases where the positions of water molecules are not well-defined; for example, in adenine-thymine plus 12 waters there are two water molecules that are not really bound at all. All the other waters in this example are bound.

The authors declare no competing financial interest.



ACKNOWLEDGMENTS This research was mostly carried out at Shanghai University under a High-End Foreign Expert Grant for RK. The computations were undertaken on the NCI National Facility in Canberra, Australia, which is supported by the Australian Commonwealth Government.

4. CONCLUSION Microsolvation can be a significant factor in modeling biological systems, affecting both structure and energetics, due to the importance of hydrogen bonding. Traditional solvation models, such as polarizable continuum methods cannot take into account such bonding, necessitating inclusion of some explicit water molecules. Indeed, even within the SMFA approach, added refinement was required to handle hydrogen bonding correctly. We could thus reproduce the dramatically changing energetics of glycine on solvation and structural changes to the adenine-thymine base-pair complex. These calculations indicate that the SMFA approach is able to handle the complexity of hydrogen-bonded systems, certainly with respect to obtaining explicitly solvated structures and the more subtle and sensitive effects on their energetics. Our investigations on these systems, especially with the benefit of reduction of calculation cost, give us confidence to start looking at real biological systems for which an explicitly solvated picture is believed to be necessary. This work is currently underway.





REFERENCES

(1) Rueda, M.; Kalko, S.; Luque, F.; Orozco, M. The structure and dynamics of DNA in the gas phase. J. Am. Chem. Soc. 2003, 125, 8007−8014. (2) Nguyen, B.; Neidle, S.; Wilson, W. D. A role for water molecules in DNA-ligand minor groove recognition. Acc. Chem. Res. 2009, 42, 11−21. (3) Mennucci, B. Polarizable continuum model. WIREs Comput. Mol. Sci. 2012, 2, 386−404. (4) Watson, J. D.; Crick, F. H. C. Molecular structure of nucleic acids: a structure for deoxyribose nucleic acid. Nature 1953, 171, 737− 738. (5) Dahlke, E. E.; Truhlar, D. G. Electrostatically embedded manybody expansion for large systems, with applications to water clusters. J. Chem. Theory Comput. 2007, 3, 46−53. (6) Furtado, J. P.; Rahalkar, A. P.; Shanker, S.; Bandyopadhyay, P.; Gadre, S. R. Facilitating minima search for large water clusters at MP2 level via molecular tailoring. J. Phys. Chem. Lett. 2012, 3, 2253−2258. (7) Yang, Z.; Hua, S. G.; Hua, W. J.; Li, S. H. Structures of neutral and protonated water clusters confined in predesigned hosts: a quantum mechanical/molecular mechanical study. J. Phys. Chem. B 2011, 115, 8249−8256.

AUTHOR INFORMATION

Corresponding Author

*R. Kobayashi. E-mail: [email protected]. Phone: +612 6125 5986. ORCID

Rika Kobayashi: 0000-0002-0672-833X 340

DOI: 10.1021/acs.jpca.6b10919 J. Phys. Chem. A 2017, 121, 334−341

Article

The Journal of Physical Chemistry A (8) Collins, M. A. Molecular forces, geometries and frequencies by systematic molecular fragmentation including embedded charges. J. Chem. Phys. 2014, 141, 094108. (9) Pruitt, S. R.; Addicoat, M. A.; Collins, M. A.; Gordon, M. S. The fragment molecular orbital and systematic molecular fragmentation methods applied to water clusters. Phys. Chem. Chem. Phys. 2012, 14, 7752−7764. (10) Kitaura, K.; Ikeo, E.; Asada, T.; Nakano, T.; Uebayasi, M. Fragment molecular orbital method: an approximate computational method for large molecules. Chem. Phys. Lett. 1999, 313, 701−706. (11) Collins, M. A. Systematic fragmentation of large molecules by annihilation. Phys. Chem. Chem. Phys. 2012, 14, 7744−7751. (12) Netzloff, H. M.; Collins, M. A. Ab initio energies of nonconducting crystals by systematic fragmentation. J. Chem. Phys. 2007, 127, 134113. (13) Addicoat, M. A.; Collins, M. A. Accurate treatment of nonbonded interactions within systematic molecular fragmentation. J. Chem. Phys. 2009, 131, 104103. (14) Collins, M. A.; Bettens, R. P. A. Energy-based molecular fragmentation methods. Chem. Rev. 2015, 115, 5607−5642. (15) Hua, W.; Fang, T.; Li, W.; Yu, J.-G.; Li, S. Geometry optimizations and vibrational spectra of large molecules from a generalized energy-based fragmentation approach. J. Phys. Chem. A 2008, 112, 10864−10872. (16) Sahu, N.; Gadre, S. R. Molecular tailoring approach: a route for ab initio treatment of large clusters. Acc. Chem. Res. 2014, 47, 2739− 2747. (17) Jensen, J. H.; Gordon, M. S. J. Am. Chem. Soc. 1995, 117, 8159− 8170. (18) Aikens, C. M.; Gordon, M. S. Incremental solvation of nonionized and zwitterionic glycine. J. Am. Chem. Soc. 2006, 128, 12835−12850. (19) Bachrach, S. M. Microsolvation of glycine: A DFT study. J. Phys. Chem. A 2008, 112, 3722−3730. (20) Choi, C. H.; Re, S.; Feig, M.; Sugita, Y. Quantum mechanical/ effective fragment potential molecular dynamics (QM/EFP-MD) study on intra-molecular proton transfer of glycine in water. Chem. Phys. Lett. 2012, 539, 218−221. (21) Yogeswari, B.; Kanakaraju, R.; Abiram, A.; Kolandaivel, P. Molecular dynamics and quantum chemical studies on incremental solvation of glycine. Comput. Theor. Chem. 2011, 967, 81−92. (22) Jacquemin, D.; Michaux, C.; Perpete, E. A.; Frison, G. Comparison of microhydration methods: protonated glycine as a working example. J. Phys. Chem. B 2011, 115, 3604−3613. (23) Alonso, J. L.; Pena, I.; Sanz, M. E.; Vaquero, V.; Mata, S.; Cabezas, C.; Lopez, J. C. Observation of dihydrated glycine. Chem. Commun. 2013, 49, 3443−3445. (24) Perdew, J. P.; Ernzerhof, M.; Burke, K. Rationale for mixing exact exchange with density functional approximations. J. Chem. Phys. 1996, 105, 9982−9985. (25) Kumar, A.; Sevilla, M. D.; Suhai, S. Microhydration of the guanine-cytosine (GC) base pair in the neutral and anionic radical states: a density functional study. J. Phys. Chem. B 2008, 112, 5189− 5198. (26) Rosenberg, J. M.; Seeman, N. C.; Day, R. O.; Rich, A. RNA double-helical fragments at atomic resolution. II. The crystal structure of sodium guanylyl-3′,5′-cytidine nonahydrate. J. Mol. Biol. 1976, 104, 145−167. (27) Kumar, A.; Sevilla, M. D. Influence of hydration on proton transfer in the guanine−cytosine radical cation (G•+−C) base pair: a density functional theory study. J. Phys. Chem. B 2009, 113, 11359− 11361. (28) Moroni, F.; Famulari, A.; Raimondi, M. Ab initio study of the crystallographic solvation pattern of the cytosine-guanine base pair in DNA. J. Phys. Chem. A 2001, 105, 1169−1174. (29) Urashima, S.-H.; Asami, H.; Obha, M.; Saigusa, H. Microhydration of the guanine-guanine and guanine-cytosine Base Pairs. J. Phys. Chem. A 2010, 114, 11231−11237.

(30) Komura, I.; Ishikawa, Y.; Tsukamoto, T.; Natsume, T.; Kurita, N. Density-functional calculations of hydrated structures and electronic properties for G−C and A−T base pairs. J. Mol. Struct.: THEOCHEM 2008, 862, 122−129. (31) Kabelác,̌ M.; Hobza, P. Hydration and stability of nucleic acid bases and base pairs. Phys. Chem. Chem. Phys. 2007, 9, 903−917. (32) Cerón-Carrasco, J. P.; Requena, A.; Michaux, C.; Perpète, E. A.; Jacquemin, D. Effects of hydration on the proton transfer mechanism in the adenine−thymine base pair. J. Phys. Chem. A 2009, 113, 7892− 7898. (33) Schmidt, M. W.; Baldridge, K. K.; Boatz, J. A.; Elbert, S. T.; Gordon, M. S.; Jensen, J. H.; Koseki, S.; Matsunaga, N.; Nguyen, K. A.; Su, S.; et al. General atomic and molecular electronic structure system. J. Comput. Chem. 1993, 14, 1347−1363. (34) Bode, B. M.; Gordon, M. S. Macmolplt: a graphical user interface for GAMESS. J. Mol. Graphics Modell. 1998, 16, 133−138. (35) Lengsfield, B. H., III General 2nd Order MCSCF Theory - A density-matrix directed algorithm. J. Chem. Phys. 1980, 73, 382−390. (36) Chai, J. D.; Head-Gordon, M. Long-range corrected hybrid density functionals with damped atom-atom dispersion corrections. Phys. Chem. Chem. Phys. 2008, 10, 6615−6620. (37) Thanthiriwatte, K. S.; Hoehenstein, E. G.; Burns, L. A.; Sherrill, C. D. Assessment of the performance of DFT and DFT-D methods for describing distance dependence of hydrogen-bonded interactions. J. Chem. Theory Comput. 2011, 7, 88−96. (38) Kabelác,̌ M.; Hobza, P. Hydration and stability of nucleic acid bases and base pairs. Phys. Chem. Chem. Phys. 2007, 9, 903−917.

341

DOI: 10.1021/acs.jpca.6b10919 J. Phys. Chem. A 2017, 121, 334−341