HCN, Formamidic Acid, and Formamide in Aqueous Solution: A Free

Mar 25, 2016 - ... in aqueous solution–water acts as both solvent and a participant in the ..... So far we have found relatively good agreement with...
0 downloads 0 Views 3MB Size
Article pubs.acs.org/JPCB

HCN, Formamidic Acid, and Formamide in Aqueous Solution: A FreeEnergy Map Jeremy Kua* and Kyra L. Thrush Department of Chemistry and Biochemistry, University of San Diego, 5998 Alcala Park, San Diego, California 92110, United States S Supporting Information *

ABSTRACT: What chemical species might be found if water or ammonia reacts with HCN in aqueous solution under neutral conditions? Is it energetically favorable for formamidic acid, the first hydration product of HCN, to tautomerize into formamide under standard conditions? Do these molecules form stable oligomers in solution? To answer these questions, we constructed a Gibbs free-energy map of the molecules that might be present to evaluate their relative thermodynamic and kinetic stability. Our protocol utilizes density functional theory calculations, Poisson−Boltzmann implicit solvent, and thermodynamic corrections. We find that for C1 species, formamide is indeed the thermodynamic sink, although the initial barrier to hydration is ∼30 kcal/ mol. Molecules with one carbon and three heteroatoms are less stable. We also find that for HCN trimerization, although the planar sp2 six-membered ring is more stable compared to its monomers, the reverse is true for the nonplanar sp3 six-membered rings formed by trimerization of formamidic acid or formamide.



INTRODUCTION Hydrogen cyanide (HCN) and formamide are two small molecules implicated in the chemistry of the origin of life. Both molecules have been observed1−3 in the interstellar medium (ISM) and are likely to have been present on the early earth as building blocks for more complex structures. HCN was recently observed among the volatiles formed in the protoplanetary disk around a young star.4 The nucleobase adenine is formally a pentamer of HCN. Adenine was first synthesized directly by Oro and Kimball starting from relatively concentrated aqueous solutions of HCN (1−15 M).5 HCN is also a building block in Urey−Miller experiments and participates in the formation of amino acids via the Strecker reaction.6 More recently, HCN was implicated as the key building block in reaction networks that potentially lead to a range of amino acids, nucleotides, and lipids.7 Formamide is the simplest amide. Its tautomer, formamidic acid, is formally a hydrate of HCN (see Figure 1). Formamide

meteoritic material led to the formation of amino acids, sugars, nucleobases, and even nucleosides.11 The goal of this study is to generate a Gibbs free-energy map in a system starting with the small molecules HCN and NH3 in aqueous solution−water acts as both solvent and a participant in the network of reactions. The present study complements our recent work generating a free-energy map for the cooligomerization reactions of CH2O and NH3 in aqueous solution.12 Both HCN and CH2O are formed in situ in the classic Urey−Miller experiment, and we are interested in building a free-energy map of small molecules interacting in solution that can self-oligomerize or co-oligomerize. The present study also lays the groundwork for current and future work involving a more complex mixture containing CH2O, HCN, and NH3 in aqueous solution. As mentioned in previous work, “computational chemistry is particularly useful in studying complex mixtures with large and varied product distributions because we can tease out and identify the energetic contributions of the many different species present in a reaction mixture”. We had previously developed a relatively fast protocol to map the free-energy landscape (both thermodynamics and kinetics) for the preliminary steps of self-oligomerization of several watersoluble aldehydes in solution”.13−15 We validated our protocol by comparing our computed results with NMR measurements for the self-oligomerization of a 1 M solution of glycolaldehyde.16 As mentioned in previous work,12 “our calculated equilibrium concentrations of the dominant species in solution

Figure 1. Formamidic acid and formamide.

was also recently detected in a protostar.8 Two extensive reviews by Saladino et al. outline how formamide is particularly well-suited as a building block to all the natural nucleobases found in extant DNA and RNA.9,10 In recent experiments, irradiation of formamide in the presence of powdered © XXXX American Chemical Society

Special Issue: J. Andrew McCammon Festschrift Received: February 18, 2016 Revised: March 23, 2016

A

DOI: 10.1021/acs.jpcb.6b01690 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B Table 1. Energies of Stable Speciesa Eelec (a.u.) H2O −76.44744 NH3 −56.57604 1 (HCN) −93.45202 2 −169.86409 3 −150.04902 4 −169.92569 5 −169.94622 6 −189.81948 7 −206.62585 8 −226.50607 9 −246.38805 10 −266.26154 six-membered rings 11 −280.43198 12 −356.87817 13 −433.33621 14 −509.78877 15 −509.79219 acyclic formamide dimer and trimer 16 −339.87939 17 −509.81508 a

Esolv

Hcorr

Gcorr

-0.5TScorr

G298

ΔGr

−8.11 −6.29 −2.47 −5.46 −12.81 −7.68 −10.28 −6.88 −16.57 −16.73 −13.26 −13.95

15.74 23.91 12.52 30.70 39.23 31.62 31.39 23.82 65.94 58.09 50.42 42.27

2.30 10.20 −1.81 12.68 20.90 13.79 13.03 6.14 44.58 37.12 30.11 21.91

−6.72 −6.86 −7.16 −9.01 −9.16 −8.92 −9.18 −8.84 −10.68 −10.48 −10.15 −10.18

−47970.59 −35491.25 −58639.15 −106575.11 −94139.95 −106614.98 −106630.95 −119105.45 −129621.01 −142103.86 −154583.86 −167063.53

0.0 0.0 0.0 34.63 −9.64 −5.73 −21.21 −16.27 0.46 −2.96 −3.53 −3.77

−7.51 −18.81 −22.71 −29.19 −23.65

44.00 62.56 81.44 100.16 99.59

24.74 39.43 56.97 73.73 72.60

−9.63 −11.57 −12.24 −13.22 −13.49

−175946.90 −223912.30 −271876.13 −319839.60 −319837.05

−29.44 −24.25 −17.49 −10.37 −7.82

−18.16 −25.21

65.47 99.25

41.71 70.59

−11.88 −14.33

−213242.15 −319854.14

−22.66 −24.91

All energies are in kcal/mol except where indicated (Eelec is in a.u.).

radius of 1.40 Å. As in previous work, the solvation energy was calculated at the optimized gas-phase geometry because in most cases there is practically no change between the gas-phase and implicit solvent optimized geometries. The electronic energy of the optimized gas-phase structures and the solvation energy are designated Eelec and Esolv, respectively, in Table 1. Even though the solvation energy contribution is to some extent a freeenergy correction, it certainly does not account for all of the free energy. The analytical Hessian was calculated for each optimized structure, and the gas-phase energy is corrected for zero-point vibrations. Negative eigenvalues in transition-state calculations were not included in the zero-point energy (ZPE). The temperature-dependent enthalpy correction term is straightforward to calculate from statistical mechanics where we assume that translational and rotational corrections are a constant times kT, that low-frequency vibrational modes will generally cancel out when calculating enthalpy differences, and that the vibrational frequencies do not change appreciably in solution. The combined ZPE and enthalpy corrections to 298 K are designated Hcorr and the corresponding gas-phase Gibbs freeenergy correction to 298 K is designated Gcorr in Table 1. The corresponding free-energy corrections in solution are much less reliable.24−26 Changes in free-energy terms for translation and rotation are poorly defined in solution, particularly as the size of the molecule increases. Additional corrections to the free energy for concentration differentials among species (to obtain the chemical potential) can be significant, especially if the solubility varies among the different species in solution. Furthermore, since the reactions being studied are in solution, the free energy being accounted for comes from two different sources: thermal corrections and implicit solvation. Neither of these parameters is easily separable, nor do they constitute all the required parts of the free energy under our approximations of the system. To estimate the free energy, we followed the approach of Lau and Deubel27 who assigned the solvation entropy of each

(monomers and dimers) agreed very well with experiment. We also did relatively well in estimating rate constants and activation barriers.” Our present study highlights a weakness in our previous protocol as we have ventured to study a wider range of molecules and functional groups. In particular, two tautomerization reactions seem to yield a “negative” barrier, requiring us to consider establishing a modified protocol. In the Results and Discussion section, we compare the original and modified protocols. There are two reaction networks that form the scope of this study. The first involves all reasonable C1 species that might be formed starting from HCN in the presence of NH3 and H2O as potential coreactants. The second explores both the direct and stepwise trimerization of HCN, formamide, and formamidic acid into six-membered rings. The present work also extends the free-energy map of inter-related chemical species, allowing us to directly compare the solution-phase free energies across multiple studies.12−16 In these studies, pH 7 and 25 °C are the baseline standard conditions we have chosen for the majority of our previous work. Current and future work involves extending our maps to different pH and temperature regimes.



COMPUTATIONAL METHODS The protocol described below is similar to our previous studies, and much of the text in this section is reproduced from earlier work for clarity and ease of the reader.12−16 All calculations were carried out using Jaguar 6.017 at the B3LYP18−21 flavor of density functional theory (DFT) with the 6-311G** basis set. A comparison of our chosen level of theory, basis set, and implicit solvent scheme, with other methods can be found in our previous work.14 To maximize the probability of finding the global minimum, multiple conformers of each species are calculated. The Poisson−Boltzmann (PB) continuum approximation22,23 was used to describe the effect of water acting as an implicit solvent with a dielectric constant of 80.4 and a probe B

DOI: 10.1021/acs.jpcb.6b01690 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B Table 2. Transition-State Energiesa Eelec (a.u.) ts1a −302.91486 ts1b −322.79003 ts3a −359.52212 ts3b −379.41258 ts3c −302.96904 ts4a −379.38602 ts4b −399.27193 ts4c −322.85383 ts5a −379.41665 ts5b −399.29461 ts6a −399.30458 ts6b −419.16533 ts6c −342.74320 six-membered rings ts1t −280.29281 ts4t −509.71026 ts5t −509.77560 ring dehydration ts12b −509.77651 ts13b −586.23635 ts14b −662.69557 stepwise trimerization ts5d −492.78295 ts16d −662.72584 ts17d −662.70112 a

Esolv

Hcorr

Gcorr

-0.5TScorr

G298

ΔGr

−22.37 −15.90 −22.87 −25.47 −15.85 −26.31 −22.97 −12.89 −24.38 −16.88 −20.46 −20.56 −14.80

67.69 59.21 97.45 88.81 68.63 86.75 78.61 60.72 88.77 79.37 81.39 71.31 53.98

42.84 34.71 69.85 61.06 44.81 59.60 51.49 37.41 63.04 53.33 56.43 45.41 30.73

−12.42 −12.25 −13.80 −13.87 −11.91 −13.58 −13.56 −11.66 −12.86 −13.02 −12.48 −12.95 −11.63

−190049.09 −202522.78 −225542.80 −238035.57 −190075.11 −238021.51 −250504.89 −202557.70 −238036.07 −250511.73 −250519.01 −262992.47 −215047.10

22.40 28.14 19.85 6.51 −3.62 20.57 16.62 −6.78 6.01 9.78 2.50 8.47 −16.75

−6.27 −25.98 −22.77

39.14 96.77 97.47

15.98 68.88 69.63

−11.58 −13.95 −13.92

−175865.14 −319791.23 −319828.30

52.32 38.00 0.93

−21.65 −25.18 −37.13

90.31 109.45 130.84

61.33 78.37 96.83

−14.49 −15.54 −17.01

−319835.48 −367800.21 −415771.13

−6.25 −0.39 −0.73

−22.26 −28.77 −34.23

98.60 131.46 130.67

67.62 96.72 99.07

−15.49 −17.37 −15.80

−309165.18 −415781.51 −415770.67

−4.52 −11.10 −0.27

All energies are in kcal/mol except where indicated (Eelec is in a.u.).

species as half (or 0.5) of its gas-phase entropy. This was based on work by Wertz and Abraham proposing that upon dissolving in water, molecules lose a constant fraction of their entropy, typically close to 0.5. Wertz28 had originally proposed a factor of 0.46 on the basis of a small suite of molecules, but as the diversity of molecules was expanded, Abraham29 found a wider range between 0.4 and 0.6. Following our established protocol, we began this study by using a factor of 0.5 for the entropic correction, designated −0.5TScorr in Table 1 and calculated as 0.5(Gcorr − Hcorr). We are not the only ones using this approach; there are other recent computational studies using a similar factor of 0.5 in unrelated systems.30−32 The free energy of each species, designated G298 in Table 1, is the sum of Eelec, Esolv, Hcorr, and −0.5TScorr. Although we calculated multiple conformers, only the most stable conformer for each unique molecular species (both minima and transition states) is reported in our free-energy map. ΔG values are calculated from the difference in G298 between the reactants and products, and therefore, they include the zero-point energy, enthalpic, and entropic corrections to 298 K, for a reaction in aqueous solution. In Table 1, the rightmost column (ΔGr) is the relative free energy of each species with respect to HCN, H2O, and NH3 as the reference states. These three molecules are always assigned values of ΔGr = 0.0 kcal/mol. This allows us to quickly and easily visualize a map of the energy landscape, for the myriad reactions that can take place. Although water is a reactant in hydration reactions, concentration corrections are not included in this landscape, the advantages and disadvantages of which are discussed in our previous work.16 For transition-state calculations, additional water molecules were explicitly added to the system to find the lowest energy barrier. We tried different numbers of water molecules and

report the free energies of the optimal structures with the lowest barriers. All calculated transition states have one large negative eigenvalue corresponding to the reaction coordinate involving bond breaking/forming and accompanying proton transfer. The corresponding energy components of the transition states are listed in Table 2. In our previous study examining the equilibrium distribution of glycolaldehyde and its oligomers in solution,16 we found that this protocol for calculating free energies and estimating equilibrium concentrations were in good agreement with NMR experimental results. The activation barriers compared to experiment are reasonable but not as close, and our protocol typically overestimated transition-state barriers by 2−3 kcal/ mol. On the other hand, our protocol performs well in calculating the relative free energies of stable species, typically within 0.5 kcal/mol of experimental results, or an uncertainty of within a factor of 2.3 in terms of equilibrium constant ratios. The barriers were calculated in reference to the separated reactants rather than a preassociated complex because in previous cases we found differences of 0−2 kcal/mol between the two methods,12 and we were willing to tolerate this error in favor of a cleaner protocol that simplified the task of determining relative free energies in a network of chemical species and reactions between them. This study alerted us to potential outliers (tautomerization reactions) requiring us to consider a slightly modified protocol involving rescaling the entropy factor of transition states. Details are provided in the Results and Discussion section.



RESULTS AND DISCUSSION The presentation of our results and subsequent analysis will be divided into three sections. (1) First, we will present the energy map for all relevant C1 species starting from the reaction of C

DOI: 10.1021/acs.jpcb.6b01690 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B

(see Figure 3, left). The reaction barrier is G298(ts1b) − [G298(1) + 3 G298(H2O)] = +28.14 kcal/mol. Since the reactants are all reference molecules, ΔGr(ts1b) = +28.14 kcal/ mol.

HCN with H2O or NH3. Free energies of all species will be calculated using the 0.5 entropy factor as described in the Computational Methods section. (2) Next, we will examine the problem that arises with the tautomerization barrier between formamidic acid and formamide. We will propose a potential modified protocol for scaling the relative barriers to avoid the artifact of a “negative” barrier. (3) Finally, we present the energy map for the trimerization of HCN, formamidic acid, and formamide into six-membered rings. Both the direct trimolecular collision and stepwise pathways will be considered. Free-Energy Map of HCN Reacting with NH3 and H2O. Figure 2 shows the network of all relevant C1 species starting

Figure 3. Transition states for hydration of HCN (left) and formamide (right).

Next, consider the hydration of formamide, 5 + H2O → 9 with transition-state ts5b. The free-energy change for this specific reaction is G298(9) − [G298(5) + G298(H2O)] = +17.68 kcal/mol. However, since the starting point in this specific reaction is formamide with ΔGr(5) = −21.21 kcal/mol, the relative free energy of 9 with respect to the reference molecules is ΔGr(9) = G298(9) − [G298(HCN) + 2 G298(H2O)] = −3.53 kcal/mol as shown in Table 1. The free energy for this specific reaction can also be found by comparing ΔGr(5) and ΔGr(9), i.e., their difference is +17.68 kcal/mol. Two additional catalytic water molecules are also involved in the optimal transition state (see Figure 3, right). The barrier for this reaction, G298(ts5b) − [G298(5) + 3 G298(H2O)] = +30.99 kcal/mol, but relative to the reference states, ΔGr(ts3b) = +6.51 kcal/mol as indicated in Table 2. The optimum transition state with the lowest energy barrier for HCN hydration involves two additional water molecules (in addition to reactants) to facilitate proton transfer, as shown in Figure 3. Because four bonds are being broken and four bonds are being formed, we refer to this as an 8-center transition state. The CN triple bond has only lengthened marginally to 1.19 Å, the newly forming C−O bond is 1.75 Å, while N−H and O−H bonds directly involved in the transition state are between 1.12 and 1.38 Å. With the exception of CN, compared to equilibrium distances, all the bonds being broken and formed are, as expected, about 20−40% longer. The HCN bond angle is 140.2°. The optimum transition state for formamide hydration (also in Figure 3) shows similar features. In previous work on aldehyde formation, similar 8-center optimal transition states were found. As mentioned in previous work,12 “the calculated energy barriers for the corresponding 4-, 6-, and 10center transition states, involving zero, one, and three catalytic water molecules have higher barriers”.15 We were not able to find much in terms of experimental values for barriers or equilibrium distributions for this network of reactions. The best we could do was a data extrapolation from a study by Miller and co-workers measuring the hydrolytic stability of HCN and formamide.33 Because the reaction is slow at room temperature and neutral pH, the data points available across a wide pH range are only at higher temperatures. At lower temperatures, there are some data points either under strongly acidic or basic conditions. Details of our extrapolation method are provided in Supporting Information. In short, we estimated the experimental barriers for HCN and formamide

Figure 2. C1 species in the reaction network of HCN, H2O and NH3.

with HCN. Stable species with energy minima are labeled 1 through 10. Single-headed arrows represent the reaction coordinate for addition of NH3 and H2O. Double-headed arrows represent tautomerization reactions. Transition states are labeled with the prefix ts followed by a numeral representing the “reactant” species. The suffix letters represent (a) addition/removal of NH3, (b) addition/removal of H2O, and (c) tautomerization. The energies of the stable species and transition states are in Table 1 and 2 respectively. The three small molecule “reactants” (HCN, H2O, NH3) are the reference states and assigned a relative free energy of zero. The relative energies are listed in the rightmost column by ΔGr in Tables 1 and 2. Here are two examples of how these values are calculated. The hydration reaction of HCN to formamidic acid is 1 + H2O → 4 with transition-state ts1b. The free-energy change for this reaction (referring to the values in Table 1) is G298(4) − [G298(1) + G298(H2O)] = −5.73 kcal/mol. Because the reactants are all reference molecules and assigned values of zero, ΔGr(4) = −5.73 kcal/mol. The optimal transition state for this reaction has two additional catalytic water molecules D

DOI: 10.1021/acs.jpcb.6b01690 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B hydration at pH 7 and 25 °C to be 27.9 and 27.0 kcal/mol, respectively. Our calculated values as shown above are 28.1 and 31.0 kcal/mol, respectively. While our calculated HCN hydration barrier is very close to the estimated experimental value, this is probably fortuitous. The larger difference in the formamide hydration is more typical of what we have seen in our previous work; that is, our calculated values are overestimates (typically 2−3 kcal/mol, but in this case it is 4 kcal/mol). There is however much uncertainty in our estimated extrapolation, so one should not put too much stock on excellent or poor agreement in these two cases. Coincidentally, our value is very similar to a recent computational study using a larger cluster with additional explicit water molecules in a continuum solvent that reports a barrier of 31.4 kcal/mol.34 These values are also close to the 29−32 kcal/mol predicted for peptide hydrolysis via No Barrier Theory.35 (Computational studies in gas-phase systems report much higher barriers and are not relevant to the present study.) Another experimental study at higher temperature, which delves into the acid- and base-catalyzed mechanisms of formamide hydrolysis, estimates that the water reaction (at pH 7) extrapolated to 25 °C is in the ∼30 kcal/mol range,36 but they do not provide an exact value. One other possibility that may contribute to discrepancies is that we have chosen to restrict our energy map to consider only neutral species with no “formal” charges. Hence, we may have overestimated the formamide hydration barrier if it proceeds more favorably through a zwitterionic intermediate.34 Considering the overall free-energy landscape of the C1 network, we see that formamide is the thermodynamic sink with ΔGr(5) = −21.21 kcal/mol. The second-most stable molecule is formic acid with ΔGr(6) = −16.27 kcal/mol. 5 and 6 are the only molecules in this network containing carbonyl groups. Formation of the imines 3 and 4 by addition of NH3 and H2O, respectively, to HCN are both exergonic. On the other hand, formation of formaldehyde oxime 2, containing an N−O bond, is highly endergonic (+34.63 kcal/mol). Hence, we did not consider further isomerization to the even less stable nitrosomethane, nor did we consider any other hydration reactions involving N−O bond formation. We also chose not to consider additions of NH3 involving N−N bond formation in the present study for similar reasons. Nucleophilic addition to molecules containing CO or CN (i.e., 3−6) leads to species with three heteroatoms attached to carbon (7−10). All these reactions are endergonic; ΔGr for these molecules is in the +0.5 to −3.8 kcal/mol range (i.e., similar in relative stability to HCN). Thermodynamic stability of these species increases with number of OH groups. Thermodynamically, the most favorable route from HCN to formamide is 1 → 3 → 8 → 5. This is also the most favorable route kinetically, traversing the lowest barriers among the different pathways (i.e., going through the transition-states ts1a → ts3b → ts5a). Although the addition of NH3 to HCN is more favorable than hydration both kinetically and thermodynamically (ΔGr(ts1a) < ΔGr(ts1b) and ΔGr(3) < ΔGr(4)), the kinetics for addition to an imine favors hydration over ammonia addition. Starting from the imines 3 or 4, the hydration barriers are lower than for addition of ammonia, i.e., ΔGr(ts3b) < ΔGr(ts3a) and ΔGr(ts4b) < ΔGr(ts4a). However, starting from the carbonyls 5 or 6, ammonia addition is more favorable than hydration, i.e., ΔGr(ts5a) < ΔGr(ts5b) and ΔGr(ts6a) < ΔGr(ts6b). In four out of five pairs of reactions comparing the addition of NH3 versus H2O, the difference within each pair is

4−6 kcal/mol. The one exception is the anomalously low barrier for ΔGr(ts3b), which is 13 kcal/mol lower than ΔGr(ts3a). After examining the structures and trying many different conformers, we have no good explanation for this anomaly and simply report the results as they stand. Transitionstate structures for ts3a, ts3b, ts4a, ts4b are shown in Figure 4.

Figure 4. Transition states for NH3 and H2O hydration to 3 and 4.

Problems arise in the tautomerization reactions. Considering the tautomerization barrier between formamidic acid and formamide, we find that the relative energy of the transition state, ΔGr(ts4c) = −6.78 kcal/mol lies slightly below the energy of formamidic acid, ΔGr(4) = −5.73 kcal/mol. The optimized transition state does not indicate any anomalies. (The structure of ts4c is found in Figure 5, bottom row, middle.) In our original protocol, the barriers were calculated in reference to “infinitely”-separated reactants instead of a preassociated complex, because there was little difference in the results after the entropic factor of 0.5 was applied. We see a similar negative barrier in the self-tautomerization of formic acid. ΔGr(ts6c) = −16.75 kcal/mol, which is marginally lower than ΔGr(6) = −16.27 kcal/mol. We do not see a negative barrier in the self-tautomerization of formamidine: ΔGr(ts3c) = −3.62 kcal/mol, which is higher than ΔGr(2) = −9.64 kcal/ mol. The problem is analyzed in the next section, along with a possible solution using a modified protocol. Rescaling the Barriers in the C 1 Network by Perturbing the Entropy Factor. In our original protocol, the number of molecules in the reactants, products, and transition states vary, but our calculations still yielded good results (at least when compared with available experimental data). A perennial problem in computational chemistry is choosing the appropriate approximations for the desired system of interest, so that error cancellation yields useful results that might sufficiently represent the system. The largest approximation in our protocol was using half the gas phase entropy to estimate the free energy of each species. This allowed us to forego optimizing separate reactant complexes and product complexes for each reaction (an approach otherwise used to more accurately determine barriers in gas phase reactions while eliminating the need for counterpoise corrections, which may E

DOI: 10.1021/acs.jpcb.6b01690 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B

Figure 5. Reactant and product complexes of three representative reactions.

Table 3. Energy Comparisons using Reactant and Product Complexes in the Three Representative Reactionsa complexes 1C → 3C 1C → ts1a 3Cii → 8C 3Cii → ts3b 4C → 5C 4C → ts4c a

separated molecules

ΔH

ΔG (full ΔS)

ΔG (0.5 ΔS)

−8.94 +19.40 +6.87 +16.33 −13.60 +2.50

−5.21 +25.19 +8.15 +19.85 −14.21 +4.88

−11.71 +26.18 +5.82 +20.94 −16.02 +8.26

1 1 3 3 4 4

+ NH3 → 3 + NH3 + 2 H2O → ts1a + H2O → 8 + 3 H2O → ts3b →5 + 2 H2O → ts4c

ΔG (0.5 ΔS)

ΔG (rescaled)

−9.16 +22.40 +6.68 +16.15 −15.48 −1.05

−9.16 +21.29 +6.68 +17.52 −15.48 +2.73

All energies are in kcal/mol.

not lead to the appropriate error cancellation in our protocol when solvation is taken into account). In our protocol we are able to directly compare the relative free energies of various chemical species in solution with different molecular formulas, which would be otherwise challenging if separate reactant/ product complexes were used. Our scheme worked well for a range of reactions involving oligomerization of aldehydes and co-oligomerization of aldehydes with ammonia in aqueous solution.12−16 However, clearly there is a problem for at least two of the tautomerization reactions in the present study. Our modified protocol makes use of reactant and product complexes to scale the entropy factor for the transition-state free energies (but not the minima), and that is their only contribution to the energy map. Our approach is explained

using three representative examples for the three types of reactions in this system: (1) the addition of NH3 to HCN (1 → 3 via ts1a), (2) the subsequent hydration of 3 (3 → 8 via ts3b), and the tautomerization of formamidic acid and formamide (4 → 5 via ts4c). These three reactions were chosen because of their inherent interest in this map. The first two reactions are along the most favorable route for the conversion of HCN to formamide, while the tautomerization of formamide is the most relevant involving two different species (the other two examples are self-tautomerizations). A summary of the results can be found in Table 3. The raw data for reactant and product complexes can be found in Supporting Information. How might key energy values change if complexes are used? We first consider the addition of NH3 to F

DOI: 10.1021/acs.jpcb.6b01690 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B

exergonic ranging from −13.60 to −16.02 kcal/mol. More importantly all the barriers are now positive. ΔH has the lowest barrier of +2.50 kcal/mol, and the two ΔG values are +4.88 and +8.26 kcal/mol, respectively. So we expect the barrier to be small (in the 2−4 kcal/mol range). Is it possible to correct for the negative barrier in our original protocol but avoid the problem of different complexes representing the same species such as 3C and 3Cii? Can we retain our approach of representing all relative free energies on a single map? So far we have found relatively good agreement with experiment (where available) for the ΔGr values for nontransition-state species across several small systems. In our transition states, we often overestimate the barrier by 2−3 kcal/ mol in the original protocol. Taking all this into account, our suggested modified protocol retains the ΔGr values for nontransition-state species (i.e., 1 to 10) using the 0.5 entropy factor. However, for transition states, we modify the entropy factor by taking into account the reactant and product complexes. This may not be optimal, but it seems to work for now, and we will have to revisit this approach as we explore more varied systems. The way we choose new entropy factors for the transition states is by requiring that the free-energy change between a reactant and its complex (and also between a product and its complex) be exactly zero. This allows us to keep 1 to 10 on the same free-energy scale. Table 4 summarizes the data for our

HCN. If we do not use complexes, the net reaction is given by 1 + NH3 → 3, and the forward barrier can be calculated by considering the reaction 1 + NH3 + 2 H2O → ts1a (two catalytic water molecules are found in the optimal transition state). The free-energy changes are −9.16 and +22.40 kcal/mol, respectively (the first two entries in the Separated Molecules ΔG (0.5 ΔS) column). This is our original protocol. Note that none of the free-energy changes in Table 3 are labeled ΔGr with the subscript r because these represent free-energy changes for specific reactions and are not calculated relative to only the reference molecules. If instead we used reactant and product complexes that were identical in atomic composition to the transition state (1C and 3C, see top row of Figure 5), we can now write the net reactions as 1C → 3C and the forward barrier as 1C → ts1a (the first two entries in the leftmost column of Table 3). The complexes 1C and 3C were both generated starting from the optimized ts1a and marching downhill to the local minima on either side of the reaction coordinate. In many computational studies, enthalpies rather than free energies are reported because the latter often fare more poorly when compared to experimentally measured barriers. We examine ΔH, ΔG (full ΔS), and ΔG (0.5 ΔS) for both reactions. For ΔH (which includes the zero point energy, enthalpic correction to 298 K, and the solvation energy), there happens to be a good match for the overall energy change in the reaction, −8.94 compared to our −9.16 kcal/mol. We also see that the forward enthalpic barrier of +19.40 kcal/mol is ∼3 kcal/mol lower than our calculated +22.40 kcal/mol. As mentioned in the Computational Methods section, we find our protocol often overestimates the barrier by 2−3 kcal/mol compared to experimental results. If we compared free energies in the complexes, we now see significant variation, regardless of whether we use the full entropy change or half the entropy change (0.5 ΔS). The barriers are significantly higher at +25.19 and +26.18 kcal/mol but probably spuriously so. Consider the next reaction, the addition of water to 3 to form 8, with data listed in the middle two rows of Table 3. In our original protocol, the reaction energy and barrier are +6.68 and +16.15 kcal/mol, respectively. If we use reactant and product complexes constructed in the same way starting from ts3b (Figure 5, middle row) we see that the reaction energies vary from +5.82 to +8.15 kcal/mol, with ΔH = +6.87 very close to our original protocol of +6.68 kcal/mol. The enthalpic barrier using the complex is also very similar to the value from our original protocol, while the free-energy barriers are higher. However, now we have introduced an additional complexity. Notice that the reactant complex in this reaction is designated 3Cii to distinguish it from 3C used in the previous reaction. This is because 3Cii is designed to have the same composition as its associated transition-state ts3b. On the other hand, 3C (paired with ts1a) has a different composition than 3Cii. (3C has one less H2O compared to 3Cii.) Thus, a problem arises because it is unclear where the free energy of 3 should be in the map, given that it can participate as reactant or product in several different reactions with different complex compositions. This will be true not just for 3 but for the majority of species in this system. We will return to this point momentarily, as it impacts our choice of a modified protocol. In the tautomerization of 4 to 5 via ts4c, our original protocol gives a negative forward barrier of −1.05 kcal/mol. (The reaction is exergonic −15.48 kcal/mol.) When the complexes 4C and 5C are constructed, the reactions are

Table 4. Entropy Factors for Scaling in the Three Representative Reactions factors for C for zero free-energy change 1→3 via ts1a 3→8 via ts3b 4→5 via ts4c

1C → 1 + NH3 + 2 H2O 3C → 3 + 2 H2O 3Cii → 3 + 3 H2O 8C → 8 + 2 H2O 4C → 4 + 2 H2O 5C → 5 + 2 H2O

0.504 0.586 0.438 0.463 0.296 0.380

averaged factor for ts 0.545 0.451 0.338

three representative reactions. For example, in the first reaction, the reactant complex 1C is composed of 1, NH3, and two H2O. The free energies of the individual separated molecules are calculated as before using 0.5 as the entropy factor. However, the entropy factor for 1C is adjusted so that the change in free energy for the reaction 1C → 1 + NH3 + 2 H2O is zero. For this reaction, this factor is 0.504. We perform the same calculation for the products. The individual separated molecules 3 and two H2O use the standard factor of 0.5. When the entropy factor of 3C is adjusted to 0.586, the change in free energy for 3C → 3 + 2 H2O is zero. The entropy factor for the transition state (ts1a) is chosen to be the simple average of the entropy factors of the reactant and product complexes, i.e., 0.5(0.504 + 0.586) = 0.545 as shown in the rightmost column of Table 4. This is an arbitrary choice, although perhaps similar in arbitrariness as the original choice of 0.5 as the entropy factor in our original protocol. Should the “adjustment” for the transition-state free energy be exactly halfway between the reactant and product complexes? One might think that for highly exergonic reactions, the transition state often lies closer to the reactants along the reaction coordinate, and therefore, the factor used should be closer to the reactants. Conversely, in endergonic reactions, a transition state that lies closer to the products might be better scaled G

DOI: 10.1021/acs.jpcb.6b01690 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B

ts3c and ts6c, these values are 0.362 and 0.285, respectively. Applying these factors pushes the barriers up, eliminating all spurious “negative” barriers in our free-energy map. When an entropy factor is close to 0.5, it means that there is hardly any difference in free energy between the complex and its separated reactants. The three isomerization reactions show deviations significantly far from 0.5; that is, there is significant energetic stabilization of both reactant and product complexes relative to their separated molecules. Thus, the actual barrier should be higher due to a more significant reduction in entropy when the additional water molecules are added to the transition-state calculation relative to its separate molecules. In current and future work, as we extend our energy map to include a larger diversity of molecular structures and reaction types, we will have to pay close attention to this issue to determine if the modified protocol in the present work should be extended across the board. This may also require going back to previous work to examine if rescaling the barriers will be needed. Trimerization of HCN, Formamidic Acid, and Formamide. Because an overarching interest of our research group is to generate maps for oligomerization reactions, in this final section, we consider the trimerization of HCN, formamidic acid, and formamide into six-membered rings. Because all the reactions examined in this section only involve addition or elimination reactions and no tautomerizations are considered, all calculated free energies use the original protocol with a 0.5 entropy factor. The concerted trimerization reactions are shown in Figure 6, and the energies of the six-membered rings are found in Table

using a factor closer to the products. At present, it is unclear how such a sliding scale should be defined, so we have chosen the halfway point to keep the modified protocol as simple as possible. Using the same approach for the other two representative reactions (3 → 8 via ts3b and 4 → 5 via ts4c) leads to scaled entropy factors of 0.489 and 0.338 for ts3b and ts4c, respectively. The rightmost column of Table 3 shows the rescaled ΔG values using our modified protocol. For net reactions with no transition states involved, the free-energy changes remain the same as the original protocol since molecules 1 through 10 (and also NH3 and H2O) retain the original entropy factor of 0.5. However, for barriers involving transition states, the application of the new scaling factor to the ts species leads to rescaled ΔG values. For the first two reactions, the changes are minimal (+21.29 instead of +22.40 and +17.42 instead of +16.15 kcal/mol, respectively) because the rescaled factors (0.545 and 0.451, see Table 4) are both relatively close to 0.5. The major difference is in the tautomerization reaction. When the entropy factor of 0.338 is applied to ts4c, the new rescaled ΔG is +2.73 instead of −1.05 kcal/mol (i.e., the barrier is no longer negative). Interestingly, and possibly coincidentally, it is very close to the ΔH value of +2.50 kcal/mol when complexes are used. Other computational studies report barriers in the range of 1−5 kcal/mol for the tautomerization of formamidic acid to formamide.37−39 The utility of constructing a single map is to locate the relative energy of all species with respect to the reference molecules. To calculate G298, we still use the same Eelec, Esolv, Hcorr values as found in Table 2. In the old protocol, G298 = Eelec + Esolv + Hcorr + (−0.5TScorr), where −0.5TScorr = 0.5(Gcorr − Hcorr). In the modified protocol, G298 = Eelec + Esolv + Hcorr + (−TScorr), where −TScorr = n × (Gcorr − Hcorr) and n is the new rescaled entropy factor. Table 5 shows rescaled values of ΔGr Table 5. Rescaled Transition-State Relative Free Energiesa

ts1a ts1b ts3a ts3b ts3c ts4a ts4b ts4c ts5a ts5b ts6a ts6b ts6c a

entropy factor

rescaled -TScorr

rescaled G298

rescaled ΔGr

unscaled ΔGr

0.545 0.476 0.481 0.451 0.362 0.451 0.477 0.338 0.498 0.458 0.507 0.407 0.285

−13.54 −11.65 −13.26 −13.70 −8.63 −12.24 −12.93 −7.88 −12.82 −11.92 −12.65 −10.54 −6.61

−190050.20 −202522.18 −225542.26 −238034.20 −190071.83 −238020.18 −250504.25 −202553.92 −238036.02 −250510.63 −250519.18 −262990.06 −215042.08

21.29 28.74 20.39 7.88 −0.34 21.90 17.25 −3.00 6.06 10.88 2.33 10.88 −11.73

22.40 28.14 19.85 6.51 −3.62 20.57 16.62 −6.78 6.01 9.78 2.50 8.47 −16.75

Figure 6. Concerted trimerization of HCN, formamidic acid, and formamide.

All energies are in kcal/mol.

1. The trimerization of HCN into s-triazine (1 → 11) is significantly exergonic. ΔG of this reaction is ΔGr(11) − 3 ΔGr(1) = −29.44−3(0.00) = −29.44 kcal/mol. The entropic loss is more than compensated by delocalization and resonance stabilization in the electronic energy. (The solvation energy change is close to zero: ΔEsolv = −7.51 + 3(2.47) = −0.10 kcal/ mol.) The trimerization of formamidic acid into the triazinane (4 → 14) is slightly endergonic. ΔG of this reaction is ΔGr(14) − 3 ΔGr(4) = −10.37−3(−5.37) = +5.74 kcal/mol. The most stable conformer of 14 has all three pendant hydroxyls on one side of the ring as shown in Figure 7 (left). The three hydrogen

using the rescaled entropy factors for each transition state. For comparison, the unscaled (n = 0.5) factors are shown in the rightmost column. Most of the rescaled barriers are close to the unscaled barriers. This is because most of the entropy factors range between 0.45 and 0.55. We found this to be the case when we developed our original protocol, leading us to the simplifying approximation of a 0.5 factor.12 The tautomerization transition-state factors, however, significantly deviate from 0.5. For ts4c, this value is 0.338. For the self-tautomerizations, H

DOI: 10.1021/acs.jpcb.6b01690 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B

a stepwise trimerization of formamide into 14. Since formamide is the thermodynamic sink, and the barrier to tautomerization starting from formamidic acid is very low, formamidic acid would be transient compared to formamide. The most likely stepwise reaction is nucleophilic attack on the carbonyl group via the amine on a neighboring formamide molecule as shown in Figure 8. Each step is endergonic although decreasingly so. Dimerization of formamide is uphill 19.8 kcal/mol, owing to the thermodynamic stability of the monomers in solution. Addition of another monomer to form the trimer, via similar nucleophilic attack, is uphill 19.0 kcal/mol. The ring closure reaction is uphill 14.5 kcal/mol. The transition states are labeled with the suffix d for the stepwise pathway. A representative 8-center transition state, ts16d, is shown in Figure 8. The endergonicity of the stepwise reactions leads to relatively high barriers (37.9 and 32.8 kcal/mol), although these are significantly lower than the direct trimerization. The final ring closure has a modest barrier of 24.6 kcal/mol. One could imagine a parallel stepwise pathway from formamidic acid monomers to the oxane 15. We were unsuccessful in mapping out the energetics of this pathway because even in the first dimerization step, formamidic acid isomerizes into formamide during transition-state optimizations. This suggests that oxane 15 is unlikely to be formed (and is not observed, or at least not mentioned from a literature search). One could force the stepwise reaction toward 15 by applying additional constraints to the transition-state optimization, or by removing the catalytic water molecules, but this would result in much higher barriers. Hence we chose not to pursue this pathway in the present work.

Figure 7. Most stable ring trimer conformers of formamidic acid and formamide.

bonds between the pendant hydroxyl groups stabilizes this conformer. Formamide, the thermodynamic sink, is likely not to favorably trimerize into the oxane (5 → 15) as this reaction is highly endergonic. ΔG = ΔGr(15) − 3 ΔGr(5) = −7.82− 3(−21.21) = +55.81 kcal/mol. The most stable conformer has two pendant amine groups in the equatorial position and one in the axial (Figure 7, right). Although the triazinane 14 is formed from the less stable formamidic acid monomer, it is still relatively more stable than trioxane 15, formed from the more stable formamide monomer. This is likely due to ring strain. In previous work, we saw that the ring closure to form trioxane is endergonic,15 while the ring closure to form triazinane is exergonic.12 Transition state energies are shown in Table 2. The concerted trimerization transition states are labeled with the suffix t. The barrier for concerted HCN trimerization (via ts1t) is +52.32 kcal/mol. We did not consider nonconcerted pathways because we are examining this pathway under neutral conditions excluding ionic species, and with water being the only possible reactant molecule. For formamidic acid, the trimerization barrier (via ts4t) is +38.00−3(−15.37) = +55.19 kcal/mol, not too different than HCN. For formamide, the trimerization barrier (via ts5t) is higher at +0.93−3(−21.21) = +64.56 kcal/mol due to the stability of the monomers. Triazine 11 is the dehydrated product of triazinane 14. The pathway for stepwise dehydration (14 → 13 → 12 → 11) is stepwise exergonic. The dehydration barriers are similar to those observed in the C1 network discussed earlier, and the transition states are labeled with the suffix b because these are hydration/dehydration reactions. Since the concerted pathway requires traversing a large energy barrier (>50 kcal/mol), we investigated the possibility of



CONCLUSION One of the initial goals in our study was to generate an energy map starting with the molecules HCN, NH3, and H2O in aqueous solution. This would allow us to extend our previous free-energy maps of aqueous organics with carbon having oxidation numbers of zero (e.g., CH2O and its oligomers) to species with carbon formally having oxidation numbers of +2 (e.g., HCN, formamidic acid, formamide, formic acid). As the map is extended, we will be in a better position to understand how complex mixtures evolve in molecular composition as reaction conditions change, and as more molecules are added to the mix.

Figure 8. Stepwise trimerization of formamide. Energy values are in kcal/mol. I

DOI: 10.1021/acs.jpcb.6b01690 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B

(8) Kahane, C.; Ceccarelli, C.; Faure, A.; Caux, E. Detection of Formamide, the Simplest but Crucial Amide, in a Solar-Type Protostar. Astrophys. J., Lett. 2013, 763, L38. (9) Saladino, R.; Crestini, C.; Ciciriello, F.; Costanzo, G.; Di Mauro, E. Formamide Chemistry and the Origin of Informational Polymers. Chem. Biodiversity 2007, 4, 694−720. (10) Saladino, R.; Crestini, C.; Pino, S.; Costanzo, G.; Di Mauro, E. Formamide and the Origin of Life. Phys. Life Rev. 2012, 9, 84−104. (11) Saladino, R.; Carota, E.; Botta, G.; Kapralov, M.; Timoshenko, G. N.; Rozanov, A. Y.; Krasavin, E.; Di Mauro, E. Meteorite-Catalyzed Syntheses of Nucleosides and of Other Prebiotic Compounds from Formamide under Proton Irradiation. Proc. Natl. Acad. Sci. U. S. A. 2015, 112, E2746−2755. (12) Kua, J.; Rodriguez, A. A.; Marucci, L. A.; Galloway, M. M.; De Haan, D. O. Free Energy Map for the Co-Oligomerization of Formaldehyde and Ammonia. J. Phys. Chem. A 2015, 119, 2122−2131. (13) Kua, J.; Hanley, S. W.; De Haan, D. O. Thermodynamics and Kinetics of Glyoxal Dimer Formation: A Computational Study. J. Phys. Chem. A 2008, 112, 66−72. (14) Krizner, H. E.; De Haan, D. O.; Kua, J. Thermodynamics and Kinetics of Methylglyoxal Dimer Formation: A Computational Study. J. Phys. Chem. A 2009, 113, 6994−7001. (15) Kua, J.; Avila, J. E.; Lee, C. G.; Smith, W. D. Mapping the Kinetic and Thermodynamic Landscape of Formaldehyde Oligomerization under Neutral Conditions. J. Phys. Chem. A 2013, 117, 12658−12667. (16) Kua, J.; Galloway, M. M.; Millage, K. D.; Avila, J. E.; De Haan, D. O. Glycolaldehyde Monomer and Oligomer Equilibria in Aqueous Solution: Comparing Computational Chemistry and NMR Data. J. Phys. Chem. A 2013, 117, 2997−3008. (17) Jaguar, version 6.0; Schrodinger, LLC: Portland, OR, 2005. (18) Vosko, S. H.; Wilk, L.; Nusair, M. Accurate Spin-Dependent Electron Liquid Correlation Energies for Local Spin Density Calculations: A Critical Analysis. Can. J. Phys. 1980, 58, 1200−1211. (19) Becke, A. D. Density-Functional Exchange-Energy Approximation with Correct Asymptotic Behavior. Phys. Rev. A: At., Mol., Opt. Phys. 1988, 38, 3098−3100. (20) Lee, C.; Yang, W.; Parr, R. G. Development of the Colle-Salvetti Correlation-Energy Formula into a Functional of the Electron Density. Phys. Rev. B: Condens. Matter Mater. Phys. 1988, 37, 785−789. (21) Becke, A. D. Density-Functional Thermochemistry. III. The Role of Exact Exchange. J. Chem. Phys. 1993, 98, 5648−5652. (22) Tannor, D. J.; Marten, B.; Murphy, R.; Friesner, R. A.; Sitkoff, D.; Nicholls, A.; Honig, B.; Ringnalda, M.; Goddard, W. A., III Accurate First Principles Calculation of Molecular Charge Distributions and Solvation Energies from Ab Initio Quantum Mechanics and Continuum Dielectric Theory. J. Am. Chem. Soc. 1994, 116, 11875− 11882. (23) Marten, B.; Kim, K.; Cortis, C.; Friesner, R. A.; Murphy, R. B.; Ringnalda, M. N.; Sitkoff, D.; Honig, B. A New Model for Calculation of Solvation Free Energies: Correction of Self-Consistent Reaction Field Continuum Dielectric Theory for Short Range HydrogenBonding Effects. J. Phys. Chem. 1996, 100, 11775−11788. (24) Warshel, A.; Florian, J. Computer Simulations of Enzyme Catalysis: Finding out What Has Been Optimized by Evolution. Proc. Natl. Acad. Sci. U. S. A. 1998, 95, 5950−5955. (25) Wiberg, K. B.; Bailey, W. F. Chiral Diamines 4: A Computational Study of the Enantioselective Deprotonation of BocPyrrolidine with an Alkyllithium in the Presence of a Chiral Diamine. J. Am. Chem. Soc. 2001, 123, 8231−8238. (26) Nielsen, R. J.; Keith, J. M.; Stoltz, B. M.; Goddard, W. A., 3rd A Computational Model Relating Structure and Reactivity in Enantioselective Oxidations of Secondary Alcohols by (−)-Sparteine-Pd(II) Complexes. J. Am. Chem. Soc. 2004, 126, 7967−7974. (27) Deubel, D. V.; Lau, J. K. In Silico Evolution of Substrate Selectivity: Comparison of Organometallic Ruthenium Complexes with the Anticancer Drug Cisplatin. Chem. Commun. 2006, 2451− 2453.

Our present energy map of C1 species in aqueous solution, at pH 7 and 25 °C, underscores the stability of amides and carboxylic acids relative to other species when carbon has a + 2 oxidation number. As we consider oxidation and reduction reactions in future work, this will serve as a general guide in generating reaction networks involving a web of molecular species. However, the present work also highlights a weakness in our original protocol in tautomerization reactions involving Y−C=X species (where X and Y are N, O) that lead to anomalously low and spurious negative barriers. Current and future work involves further investigating the limits of our protocol. Our exploration of the oligomerization reactions of HCN and formamide in the present study has been purposefully limited to leave out ionic and zwitterionic species, which are likely to be important in stepwise pathways. As we extend the diversity of molecular and ionic species, we will likely need to make further refinements to our protocol. Thus, current and future work includes investigating the effects of changing the pH and including molecular ions in our calculations.



ASSOCIATED CONTENT

S Supporting Information *

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.jpcb.6b01690. (1) Analysis of the extraction of experimental hydration energy barriers for HCN and formamide from Miller et al.,33 (2) additional data on reactant and product complexes used in rescaling the transition-state entropy factors, and (3) XYZ coordinates of all optimized structures (PDF)



AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. Tel.: (619) 260-7970. Notes

The authors declare no competing financial interest.

■ ■

ACKNOWLEDGMENTS This research was supported by a Camille and Henry Dreyfus Teacher-Scholar award. REFERENCES

(1) Adande, G. R.; Woolf, N. J.; Ziurys, L. M. Observations of Interstellar Formamide: Availability of a Prebiotic Precursor in the Galactic Habitable Zone. Astrobiology 2013, 13, 439−453. (2) Hirota, T.; Yamamoto, S.; Mikami, H.; Ohishi, M. Abundances of HCN and HNC in Dark Cloud Cores. Astrophys. J. 1998, 503, 717− 728. (3) Ziurys, L. M.; Adande, G. R.; Edwards, J. L.; Schmidt, D. R.; Halfen, D. T.; Woolf, N. J. Prebiotic Chemical Evolution in the Astrophysical Context. Origins Life Evol. Biospheres 2015, 45, 275−288. (4) Oberg, K. I.; Guzman, V. V.; Furuya, K.; Qi, C. H.; Aikawa, Y.; Andrews, S. M.; Loomis, R.; Wilner, D. J. The Comet-Like Composition of a Protoplanetary Disk as Revealed by Complex Cyanides. Nature 2015, 520, 198−201. (5) Oro, J.; Kimball, A. P. Synthesis of Purines under Possible Primitive Earth Conditions. I. Adenine from Hydrogen Cyanide. Arch. Biochem. Biophys. 1961, 94, 217−227. (6) Miller, S. L.; Urey, H. C. Organic Compound Synthesis on the Primitive Earth. Science 1959, 130, 245−51. (7) Patel, B. H.; Percivalle, C.; Ritson, D. J.; Duffy, C. D.; Sutherland, J. D. Common Origins of Rna, Protein and Lipid Precursors in a Cyanosulfidic Protometabolism. Nat. Chem. 2015, 7, 301−307. J

DOI: 10.1021/acs.jpcb.6b01690 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B (28) Wertz, D. H. Relationship between the Gas-Phase Entropies of Molecules and Their Entropies of Solvation in Water and 1-Octanol. J. Am. Chem. Soc. 1980, 102, 5316−5322. (29) Abraham, M. H. Relationship between Solution Entropies and Gas Phase Entropies of Nonelectrolytes. J. Am. Chem. Soc. 1981, 103, 6742−6744. (30) Liang, Y.; Liu, S.; Xia, Y.; Li, Y.; Yu, Z. X. Mechanism, Regioselectivity, and the Kinetics of Phosphine-Catalyzed [3 + 2] Cycloaddition Reactions of Allenoates and Electron-Deficient Alkenes. Chem. - Eur. J. 2008, 14, 4361−4373. (31) Huang, F.; Lu, G.; Zhao, L.; Li, H.; Wang, Z. X. The Catalytic Role of N-Heterocyclic Carbene in a Metal-Free Conversion of Carbon Dioxide into Methanol: A Computational Mechanism Study. J. Am. Chem. Soc. 2010, 132, 12388−12396. (32) Liu, W. G.; Sberegaeva, A. V.; Nielsen, R. J.; Goddard, W. A.; Vedernikov, A. N. Mechanism of O-2 Activation and Methanol Production by (Di(2-Pyridyl)Methanesulfonate) (PtMe)-Me-II(OHn)((2-n)-) Complex from Theory with Validation from Experiment. J. Am. Chem. Soc. 2014, 136, 2335−2341. (33) Miyakawa, S.; Cleaves, H. J.; Miller, S. L. The Cold Origin of Life: A. Implications Based on the Hydrolytic Stabilities of Hydrogen Cyanide and Formamide. Origins Life Evol. Biospheres 2002, 32, 195− 208. (34) Gorb, L.; Asensio, A.; Tunon, I.; Ruiz-Lopez, M. F. The Mechanism of Formamide Hydrolysis in Water from Ab Initio Calculations and Simulations. Chem. - Eur. J. 2005, 11, 6743−6753. (35) Guthrie, J. P.; Pitchko, V. Hydration of Carbonyl Compounds, an Analysis in Terms of No Barrier Theory: Prediction of Rates from Equilibrium Constants and Distortion Energies. J. Am. Chem. Soc. 2000, 122, 5520−5528. (36) Slebocka-Tilk, H.; Sauriol, F.; Monette, M.; Brown, R. S. Aspects of the Hydrolysis of Formamide: Revisitation of the Water Reaction and Determination of the Solvent Deuterium Kinetic Isotope Effect in Base. Can. J. Chem. 2002, 80, 1343−1350. (37) Du, D. M.; Fu, A. P.; Zhou, Z. Y. Density Functional Theory Study of Formamide-Formamidic Acid Tautomerization. Int. J. Quantum Chem. 2004, 99, 1−10. (38) Kalia, S.; Sharma, A.; Kaith, B. S. Ab Initio Study of Gas Phase and Water-Assisted Tautomerization of Maleimide and Formamide. Proc. - Indian Acad. Sci., Chem. Sci. 2007, 119, 617−624. (39) Bell, R. L.; Taveras, D. L.; Truong, T. N.; Simons, J. A Direct Ab Initio Dynamics Study of the Water-Assisted Tautomerization of Formamide. Int. J. Quantum Chem. 1997, 63, 861−874.

K

DOI: 10.1021/acs.jpcb.6b01690 J. Phys. Chem. B XXXX, XXX, XXX−XXX