How solvent determines the molecular reactive conformation and the selectivity: Solvation spheres and energy

In solution, the solvent determines the molecular conformation and the chemical reaction viability and selectivity. When solvent-solute and solvent-solvent interactions present similar strengths, explicit salvation is the best way to describe a system. The problem to solve is how big the explicit shell should be. In this paper, we want to answer one of the fundamental questions in the implementation of explicit solvation, exactly how many solvent molecules should be added and where they should be placed. Here we determine the first solvent sphere around a molecule and describe how it controls the conformation and selectivity of a selected reaction. NMR experiments were carried out to identify the number of solvent molecules around the solute that constitutes the first solvent sphere, and the interaction between this solvent sphere and the solute was detected using DFT and QTAIM calculations. A new approach to the solvation energy is presented. Finally, we established the role of solvent molecules in the conformation of the solute and in the transition states that produce the two possible products of the reaction.

1 Introduction

Several factors determine the path a chemical reaction follows; some depend on the structural and electronic features of the reactants, but the reaction environment controls others. In solution, the solvent determines the molecular conformation and, therefore, the reaction viability and selectivity. There are several examples of the influence of the solvent on the conformation and the reactivity of organic molecules. (Reichardt and Welton, 2011) In this work, we present an NMR approach to determine the first and the external solvation spheres of a molecule and a theoretical approach to calculate the explicit solvation energy. We also present the analysis of the specific interactions between the reactant and its solvent sphere, which controls the conformation and selectivity of a reaction.

Usually, the solvent effect is modeled by continuous models, where the system is composed of a single solute molecule immersed in an infinite solvent reservoir. The solvent is considered isotropic, and the solute-solvent interactions are limited to those of electrostatic origin, and usually, no dynamic effects are considered. (Tomasi et al., 2005) Implicit solvation uses a molecular cavity following the geometry of the system and a surface charge distribution to represent the polarization of the environment. (Floris et al., 1991) This approach ignores the short-range effects of hydrogen bonding and other weak interactions. In the case of nonpolar solvents, with low dielectric constant and small dipole moment, like hexane, cyclohexane, benzene, toluene, or chloroform, the implicit models do not provide higher accuracy than gas-phase calculations. The main shortcoming of continuum models is their failure to account for chemical interactions between solute and solvent molecules in the first solvation shell, especially for solutes with hydrophobic groups, (da Silva et al., 2009) where the solvent-solvent and solvent-solute interactions have similar strength. In some cases, despite using an accurate implicit solvent model, specific solute-solvent interactions play an important role in the solution, where the addition of explicit solvent molecules (micro-solvation) significantly improves the accuracy of the properties prediction. This is the case of pKa, where the presence of the hydrogen bonding of a single water molecule is necessary to obtain a correlation with the experiments. (Kelly et al., 2006) Another example is the calculation of the electron paramagnetic resonance spectrum, where the addition of several solvent molecules helps to reduce the difference with the experimental data. (Grande-Aztatzi et al., 2013) To approximate the bulk water’s role in the solution-phase proton-transfer processes, Bayse and others have included clusters of explicit water molecules in the gas phase to provide an indirect pathway for the proton exchange assisted by the solvent. (Bayse, 2007; Bayse and Antony, 2009; Bayse, 2010; Bayse, 2011).

Chemical reactions may involve species more complex than just molecules of reagents. The chemical reactions occur in non-ideal solutions and not between the isolated molecules of reagents. The presence of various supramolecular species, including solution domains, solute clusters, solute aggregates, and supramers, along with non-aggregated solute molecules should be considered (Kononov, 2015). The solute aggregates may differ dramatically in their chemical properties from those of isolated molecules, the latter being commonly considered by traditional organic chemistry textbooks as the genuine chemical properties of the corresponding compound. The trade-off between the inaccuracies of implicit solvent models and the computational cost of explicit models is an important issue that produces several methods which provide a reasonable accuracy at an appropriate cost. Combined implicit–explicit hybrid models present a central part of the system that contains explicit solute and a few explicit solvent molecules, and the rest of the system is treated as a dielectric continuum, also called the cluster-continuum model (Sunoj and Anand, 2012; Skyner et al., 2015). On the other hand, a solute can be embedded in a solvent box and optimized with hybrid methods such as ONIOM (Chung et al., 2015) and QM/MM (Bulo et al., 2009), where some solvent molecules are included in the higher level. In these cases, the question remains valid, how many solvent molecules include within the solute level. For example, Boereboom et al. employed two spherical QM regions encompassing the two reactive moieties in the reactant, while the micro-solvation model confined five QM water molecules near a reactive site at all times, where the dual-sphere explicit solvation allows a superior description of solvent rearrangement along the entire reaction path. (Boereboom et al., 2018) Recently Hasen et al. investigated how the solvation influences the competition between the S_N2 and E2 pathways of the model F⁻ + C₂H₅Cl reaction. In this case, the system was solvated in a stepwise manner by going from the gas phase, then via micro-solvation of one to three explicit solvent molecules, then last to bulk solvation. (Hansen et al., 2022) Karmeling et al. showed that when the solvent conformations are chosen at random mixed explicit/implicit solvation models, the phosphate hydrolysis reaction barriers are much larger than the experimental one. To obtain the acceptable values, it would be necessary to sample over many solvent configurations (Kamerlin et al., 2009).

One of the fundamental questions in implementing explicit solvation is exactly how many solvent molecules should be added and where they should be placed (Tantillo, 2018). Another important question is to determine the role of the explicit solvent molecules with active participation during the reaction or just as a spectator. In this work, we present the experimental determination of the first (FSS) and external (ESS) solvent spheres around a molecule, a theoretical approach to the explicit solvation energy and the analysis of the specific interaction within it, the role of the subsequent solvent spheres, and the description of how FSS controls the conformation and selectivity of a reaction. FSS can be defined as the set of solvent molecules that are the first neighbors of the solute. ESS can be defined by the group containing all solvent molecules which are bonded to the FSS fixing their positions around the solute (Schienbein et al., 2017).

To explore these ideas, we selected a known chemical reaction whose selectivity depends on solvents with low dielectric constant and small dipole moment, where the implicit models do not account for the difference in the product ratios. This is the case of the epoxy-nitrile cyclization, a reaction developed by Gilbert Stork, which gives one of two possible rings after an intramolecular nucleophilic attack by the α-nitrile-carbanion to an oxirane ring (Stork et al., 1974). Baldwin’s rules for opening this kind of three-membered ring to form cyclic structures seem to lie between those for tetrahedral and trigonal systems, generally preferring Exo modes. (Baldwin, 1976) Stork found that the δ-epoxy-nitrile cyclization prefers to give cyclobutanes to cyclopentanes, as expected by Baldwin’s rules, leading to a general and non-photochemical synthetic method of substituted cyclobutanes as the (±)-Grandisol (Stork and Cohen, 1974). Stork considered that the oxirane ring imposes a rigidity on the system, which makes it difficult to achieve the proper collinear arrangement for an SN₂ attack to give a five-membered ring, whereas the transition state for the formation of a four-membered ring seems to allow easy attainment of collinearity. However, there are examples where one can observe the formation of a cyclopentane rather than a cyclobutane ring, as in the cases of the total syntheses of (±)-α-Cuparenone and (±)-Epilaurene (Avila-Zárraga and Maldonado, 2000). Luján-Montelongo et al. found that the cyclizations of δ-epoxy-nitrile are also controlled through the choice of cation within the alkyl metal base and the solvent (Fleming et al., 2008; Lujan-Montelongo et al., 2009). For example, using lithium hexamethyldisilylamide (LiHMDS) as the base, the reaction prefers a four-membered ring with benzene as the solvent, whereas toluene changes the preference to a five-membered ring. This work shows how the specific solute-solvent interactions determine the molecular conformation of the solute and then the selectivity of a reaction, in this case, the cyclization of δ-epoxy nitriles.

2 Methodology

2.1 Experimental methods

2.1.1 Materials

Most chemicals were obtained from commercial sources and used as received. Solvents such as ethyl acetate, hexane, acetone, dichloromethane, were distilled before being used. Benzene, toluene, and tetrahydrofuran (THF) were used freshly distilled from metallic sodium. Deuterated solvents were obtained from Cambridge Isotope Laboratories Inc. Lithium hexamethyldisilylamide (LiHMDS) was used as a standard commercial solution in toluene from Sigma-Aldrich. Laboratories Flash chromatography was carried out on hand-packed silica gel 60 (230–400 mesh) columns. Merck TLC silica gel 60 F254 aluminum plates were used for thin-layer chromatography. NMR spectra were acquired using a Bruker Avance III 750 MHz (17.6 Tesla) equipped with a TCI triple resonance cryoprobe and in a Varian-Inova 500 MHz. Chemical shifts are reported in ppm (δ) and referenced to TMS (δH = 0 ppm). The spectral data required for the study comprise 1D (¹H–NMR, ¹³C–NMR and titrations) and 2D (COSY, HSQC, HMBC, NOESY). Calorimetric data were obtained using a MicroCal VP-ITC.

2.1.2 Synthesis

The synthesis of the 2-phenyl-6,6-dimethyl-5-epoxynitrile 1) was performed by the alkylation-oxidation procedure reported by Lujan-Montelongo et al. report (Lujan-Montelongo et al., 2009). The alkylation reaction of the phenylacetonitrile with 5-iodo-2-methyl-pent-2-ene was performed in THF at -78°C, using lithium diisopropylamide (LDA) as base to give 2-phenyl-6-methylhept-5-enonitrile. Epoxidation of the alkylated nitrile was then done with meta-chloroperoxybenzoic acid (MCPBA) as oxidant to obtain 4-(3,3-dimethyloxiran-2-yl)-2-phenylbutanenitrile (1). The intermolecular cyclization reactions were performed in benzene or toluene using (LiHMDS) as base.

2.1.3 Characterization

The characterization of racemic 2-phenyl-6-methylhept-5-enonitrile, 4-(3,3-dimethyloxiran -2-yl)-2-phenylbutanenitrile (1), 2-(2-hydroxypropan-2-yl) 1-phenylcyclobutane-1- carbonitrile, (4MR) and 3-hydroxy-2,2- dimethyl-1-phenylcyclopentane-1-carbonitrile, (5MR) were performed by NMR techniques and compared with a previous report of these compounds (Lujan-Montelongo et al., 2009).

2-phenyl-6-methylhept-5-enonitrile. ¹H–NMR (500 MHz, CDCl₃): δ ppm 7.67–7.16 (m, 5H), 5.07 (tdt, J = 5.7, 2.8, 1.4 Hz, 1H), 3.78 (dd, J = 8.8, 6.2 Hz, 1H), 2.28–2.08 (m, 2H), 2.08–1.78 (m, 2H), 1.71 (q, J = 1.3 Hz, 3H), 1.61 (d, J = 1.3 Hz, 3H). ¹³C–NMR (125 MHz, CDCl₃): δ ppm 134.96, 133.80, 129.01(2C), 128.22, 127.31 (2C), 122.01, 119.98, 36.98, 36.01, 25.6, 25.78, 16.99.

4-(3,3-dimethyloxiran-2-yl)-2-phenylbutanenitrile (1). ¹H–NMR (750 MHz, CDCl₃): δ ppm 7.33–7.24 (m, 10H), 3.86–3.80 (m, 2H), 2.66 (ddd, J = 14.7, 7.8, 4.7 Hz, 2H), 2.15–1.90 (m, 4H), 1.83–1.73 (m, 1H), 1.66 (dddd, J = 14.2, 9.5, 6.1, 4.8 Hz, 1H), 1.62–1.54 (m, 1H), 1.54–1.45 (m, 1H), 1.23 (s, 6H), 1.20 (s, 3H), 1.15 (s, 3H). ¹³C–NMR (187.5 MHz, CDCl₃): δ ppm 137.00, 136.68, 130.52 (2C), 130.49 (2C), 129.58, 129.55, 128.76 (2C), 128.53(2C), 121.96, 121.84, 64.71, 64.37, 59.78, 59.67, 38.66, 38.05, 34.62, 34.14, 27.94, 27.23, 26.08, 26.06, 20.09, 20.07.

2-(2-hydroxypropan-2-yl) 1-phenylcyclobutane-1- carbonitrile, (4MR). ¹H–NMR (750 MHz, CDCl₃): δ ppm 7.61 (dd, J = 7.6, 1.7 Hz, 2H), 7.34 (dt, J = 53.1, 7.6 Hz, 3H), 3.06 (dd, J = 11.0, 8.9 Hz, 1H), 2.73 (ddd, J = 12.4, 10.7, 8.9 Hz, 1H), 2.64 (ddd, J = 12.1, 9.1, 2.6 Hz, 1H), 2.47 (qd, J = 11.0, 9.1 Hz, 1H), 2.12 (dtd, J = 11.5, 8.9, 2.6 Hz, 1H), 1.71–1.41 (m, 1H), 1.04 (s, 3H), 0.86 (s, 3H). ¹³C–NMR (187.5 MHz, CDCl₃): δ ppm 133.58, 128.07 (2C), 127.52, 126.98 (2C), 123.30, 70.55, 54.66, 41.57, 29.02, 26.24, 26.04, 19.21.

3-hydroxy-2,2-dimethyl-1-phenylcyclopentane-1-carbonitrile, (5MR). ¹H–NMR (750 MHz, CDCl₃): δ ppm 7.49 (td, J = 7.2, 1.3 Hz, 1H), 7.43–7.36 (m, 4H), 4.30 (t, J = 8.4 Hz, 1H), 2.67 (ddd, J = 13.9, 11.9, 6.5 Hz, 1H), 2.38 (dtd, J = 13.6, 9.3, 6.4 Hz, 1H), 2.21 (ddd, J = 13.9, 9.6, 4.2 Hz, 1H), 1.73–1.68 (m, 1H), 1.13 (s, 4H), 0.51 (s, 3H). ¹³C–NMR (187.5 MHz, CDCl₃): δ ppm 134.99, 128.30(2C), 128.28, 127.68(2C), 123.93, 79.56, 54.22, 49.21, 29.71, 28.98, 21.22, 14.61.

The assignment of the 1, 4MR, and 5MR NMR signals using the 2D techniques are detailed in the supplementary information.

2.1.4 Determination of the solvation spheres

2.1.4.1 Anisotropic titration

The ¹H–NMR titration curve was carried out in CDCl₃ with benzene-d₆ and toluene-d₈ as titration solvents. 1 (15 mg, 0.06 mmol) was dissolved in CDCl₃ (0.75 ml), then additions of aromatic deuterated solvent (benzene-d₆ or toluene-d₈, 0 to 7 equivalent) were done before each ¹H–NMR measure. Chemical shifts were used to establish the first solvent sphere. The data is available in supporting information. The number of molecules within the first solvation sphere was achieved by the second derivative of the chemical shifts with respect to the equivalents of an added aromatic solvent. The slope of a titration curve reaches its maximum value at the inflection point. The first derivative of a titration curve shows a separate peak for each end point. The first derivative is approximated as Δδ/ΔV, where Δδ is the change in the chemical shift between successive additions of titrant solvent. The second derivative of a titration curve, Δ²δ/ΔV², is more useful than the first derivative since the endpoint is indicated by its intersection with the volume axis. (Harvey, 1999; Denney et al., 2000). The same procedure was carried out for 4MR and 5MR.

2.1.4.2 Calorimetric measurements

Using the technique of Isothermal Titration Calorimetry (ITC) it is possible to determine the energy evolved (Δ_transH) when transferring the solute dissolved in chloroform to benzene or toluene. Since experimentally, the transfer of the solute also involves the mixing process of chloroform with benzene or toluene (Δ_obsH), we performed the titration of pure chloroform into benzene and toluene in separate experiments from where the energy of mixing of the solvents (Δ_solvH) is determined and can be subtracted from the titration of toluene or benzene with the solute dissolved in chloroform to calculate the Δ_transH:

${\mathrm{\Delta }}_{\mathit{trans}}H={\mathrm{\Delta }}_{\mathit{obs}}H-{\mathrm{\Delta }}_{\mathit{solv}}H(1)$

Experimentally the cell of the calorimeter is filled either with benzene or toluene while the syringe of the calorimeter is filled with either pure chloroform (for the experiments from which Δ_solvH will be determined) or with a 32 mM solution of the solute dissolved in chloroform (for the experiments from which Δ_transH will the calculated). Additions of 1 μL of the solution placed in the syringe are introduced into the cell. The additions produce an evolution of heat and therefore a change in temperature in the cell which the calorimeter compensates with a change in power to keep the temperature constant. The power signals (peaks) are integrated to obtain the corresponding heats associated with each of the solution additions placed in the syringe into the calorimeter cell with the Origin 7.0 software macros supplied with the calorimeter. All the experiments were done at a temperature of 298.15 K.

2.1.5 Computational methods

All structures were optimized at M062X (Zhao and Truhlar, 2008)/6-31(d,p) level with Gaussian 09 (Frisch et al., 2009). To obtain NMR spectroscopy data, the GIAO method (Ditchfield, 1974) was used at the same level. Three molecular environments utilized gas phase, implicit and hybrid solvation, to understand the role of specific solute-solvent interactions. Conductor-like polarisable continuum model, CPCM, was the implicit method chosen as a continuous medium for benzene and toluene solvent (Cossi et al., 2003). Vibrational frequency calculations were performed to characterize the nature of the stationary points. As reported, the initial geometries were established using experimental observation, where the lithium cations link oxygen and nitrogen atoms of 1. For mixed solvation (explicit-implicit), we considered conditions to obtain a model that adequately described the specific solute-solvent interactions and whose calculation level was not computationally expensive: 1) The reactant, products, and transition states were first obtained in the gas phase, and then they were surrounded by solvent molecules. 2) The data from the anisotropic titrations were used to locate the position of the first solvation sphere. 3) Due to the diffusion of the solvent molecules in the medium, it was decided to involve a greater number of solvent molecules until the NMR spectra were reproduced. The capability of the systems to reproduce the experimental chemical shifts was monitored by comparing the experimental with the theoretical data of 1, 4MR and 5MR. To increase the solvent molecules, a cluster of 30 benzene molecules was used as reported by Llanio-Trujillo et al. (2011) 4) We included two lithium cations linked to the nitrogen and oxygen atoms due to the experimental conditions including two equivalents of LDA base. 5) Finally, the set of Kohn–Sham orbitals obtained for each system in the different solvation situations was used to calculate the atomic properties using the AIMAll program package to analyze the solute-solvent interactions. (Keith, 2019).

2.1.6 Explicit solvation energy

The solvation energy reflects the environmental change around a solute from vacuum to a solvent, accounting for the cavity formation, the dipole orientation, electrostatic and specific interactions. (Pure and Chemistry, 2009) In this work, we propose a new theoretical approach to the solvation energy (E_solv) based on the knowledge of the solvation spheres. The calculation of E_solv two steps: 1) determination of the energy change associated with the formation of the solvent cavity from an optimized solvent cluster, and 2) the energy change of the solute from the gas phase to the explicit solvation, as shown in Eq. (2). In this approach, all the factors involved in the solvation are obtained from the explicit solute-solvent interactions and the cavity is not just an electrostatic space but a real hole within a solvent cluster.

$\begin{array}{ll}\hfill & E_{\mathit{solv}}=\mathrm{\Delta }{E}_{\mathit{cavity}}+\mathrm{\Delta }{E}_{\mathit{solute}}\hfill \\ \hfill & \mathrm{\Delta }{E}_{\mathit{cavity}}=E_{\mathit{solvent-cavity}}-E_{\mathit{optimized-solvent}}\hfill \\ \hfill & \mathrm{\Delta }{E}_{\mathit{solute}}=E_{\mathit{solvated-solute}}-E_{\mathit{gas-phase-solute}}\hfill \end{array}(2)$

3 Results and discussion

3.1 First solvation sphere determination

To establish the composition of the first solvation sphere of each structure in the aromatic solvents, the ¹H–NMR signals of 1, 4MR and 5MR in CDCl₃, benzene, and toluene were compared to find the useful chemical shifts to understand the solvation changes (see S.I. for details). The chemical shift of protons H₇, H₈, H₉, and H₁₀ gives more positive values than the other protons of the three structures (see Table 1 and Figure 1 for labels). Then we performed NMR anisotropy titrations using benzene-d₆ or toluene-d₈, to determine the number of solvent molecules that provoke chemical shift changes, which constitute the first solvation spheres. The NMR titrations have been used for monitoring intermolecular aromatic solvent–sugar interactions. (Vandenbussche et al., 2008) A sample of an anisotropy titrations is presented in Figure 2. The number of solvent molecules was achieved by the second derivative of the chemical shifts concerning the equivalents of an added aromatic solvent (see S.I. for details). In this way, it was found that 1 is surrounded by five benzene molecules or four toluene molecules in their first solvation spheres (Figure 2D). The NMR information was useful for modeling the explicit solvation around the molecules, locating the faces of the benzene rings in the direction of the protons with the greatest change chemical shift. The size of these spheres agrees with those da Silva et al. reported for organic anions. (da Silva et al., 2009).

TABLE 1

TABLE 1. Chemical shift changes (ppm) provoked by the aromatic solvent on 1

FIGURE 1

FIGURE 1. Cyclization of the δ-epoxy nitriles 1 to give 4MR and 5MR. Atomic labels are shown for the three structures.

FIGURE 2

FIGURE 2. ¹H–NMR of 1 using benzene-d₆ as titration solvent to determine the number of solvent molecules that provoke chemical shift changes. (A) and (B) H₁₀ Chemical shift evolution during the anisotropy titration, (C) the second derivative of the chemical shift with respect to the solvent equivalents, (D) solvent sphere model extracted from the NMR titration.

However, it was not possible to reproduce the experimental chemical shifts in the presence of just the first solvation sphere due to the mobility of solvent molecules on the surface of the solutes. The correlation between experimental and theoretical NMR data in the implicit solvent is R² = 0.98 but when the first solvation sphere is included, the correlation decreases to R² = 0.79. As we will show in the next sections, in spite of the implicit solvent reproducing the NMR shifts, it does not reproduce the energy profiles that agree with the experimental results. Then subsequent solvation spheres need to be included to locate the first sphere in the proper positions on the solute, to reproduce the NMR spectra. We surrounded the solute with its first solvation sphere, with a cluster of solvent molecules. In the case of benzene, we used the cluster reported by Llanio-Trujillo et al., which included 30 benzene molecules. (Llanio-Trujillo et al., 2011) Afterward, we removed benzene molecules until we reached the minimum number of solvent molecules to fit the chemical shift, finding that 22 molecules were enough to achieve a linear fit between the theoretical and experimental data (see S.I. for details). In the case of toluene, we surrounded the solute with the same number of toluene molecules (Frisch et al., 2009), using the Vega ZZ program (Pedretti et al., 2002) to create a solvent cluster. In this case, the correlation between experimental and theoretical NMR data is R² = 0.93; if the implicit solvent is also included, R² = 0.95. da Silva et al. suggested that when enough solvent molecules are added to saturate the solvation shell, the calculated solvation energies should converge. (da Silva et al., 2009) Unlike implicit methods to calculate the cavitation, where a volume is determined, the explicit methods must find the first solvent sphere which directly interacts with the solute and the external shell that maintains the FSS in the right place.

In fact, the number of molecules in the external solvent sphere can not be determined exactly but we have calorimetric evidence to support a similar number of solvent molecules in these solvation spheres. From the calorimetric experiments (data summarized in Table 2), we found that the amount of energy associated with the titration of chloroform over toluene or benzene is small and exothermic. However, when 1 is solved in chloroform and it is titrated over toluene or benzene, the reaction is very exothermic. Using Eq. 1, Δ_transH = -206.712 kJ/mol for the transfer of 1, from chloroform to benzene whereas from the theoretical calculations with 22 benzene molecules, Δ_transH is -194.177 kJ/mol, which is very similar considering the size of the system. For the transfer to toluene, Δ_transH is -539.157 kJ/mol, which involves a larger number of toluene molecules than in the case of benzene. However, we have to note that our interest is not to characterize the external solvent sphere but just to find the solvent environment that allows us to fix the first solvent sphere in the proper positions around the solute to understand the reaction process. From this information, we extrapolated the solvation of the anionic intermediate, the transition states, and cyclic products under the hypothesis that the solvent molecules in the first solvation sphere would remain approximately in the same regions throughout the reaction process.

TABLE 2

TABLE 2. Heats of mixing, Δ_transH, when transferring 1 dissolved in chloroform to benzene or toluene.

3.2 Solvation energy

Experimentally, the thermodynamic determination of the heat of solvation takes into account the heat of the broken solute intramolecular interactions and the heat released when each molecule of solute gets in the bulk of the solution and is surrounded by a shell of solvent molecules. (Ramirez-Gualito et al., 2009) Based on the experimental procedure, our approach to the explicit solvation energy accounts for the energies differences associated with the cavity formation and the presence of the solvent interactions with the solute, as shown in Figure 3 and Eq. 2. In this way, the E_solv of 1 in benzene is -379.65 kJ/mol and -516.01 kJ/mol in toluene.

FIGURE 3

FIGURE 3. Approach to the solvation energy, which includes the solvent cavity formation and the solute energy change from the gas phase to explicit solvent, counting the specific solvent-solute interactions.

3.3 Solvent effect on the conformation

Usually, the molecular conformation is explained in terms of the stereoelectronic effects (Alabugin, 2016) or by interactions between the molecular and solvent dipoles, (Wang et al., 2014) but there are only a few examples where the specific solute-solvent interactions are taken into account. (Schmitt et al., 2005; Cainelli et al., 2009) In our case, the conformation of 1 changes as the reaction system becomes more complex and realistic. The deprotonation of the α-carbon produces a delocalized carbanion with the nitrile group, (Corset et al., 2003; Strzalko et al., 2012; Yang and Fleming, 2017) which attacks the epoxy group in the next step of the reaction. If the carbanion is modeled without a counterion, the conformation allows the negative nitrogen atom and the phenyl ring to be as far as possible. When the counterions are included, one lithium-ion linked to the N-terminus of the nitrile group and the other in the oxygen of the oxirane, the conformation change dramatically by an intramolecular cation-π interaction, producing a conformer 104.35 kJ/mol more stable than the open conformer, from which it is not possible to obtain the products, in both gas phase and implicit solvent. The explicit aromatic solvent molecules rectify the conformation (Figure 5), interacting with the lithium atoms and allowing the molecule to assume the orientation where the delocalized carbanion can attack the epoxy group. Figure 4 shows the dihedral angle values of the 4MR and 5MR transition states in the gas phase, implicit and explicit solvents. The solvent presence provokes specific interactions which modify the dihedral angle changing the position of the carbon atom C₇ and the orientation of C₁₄ and C₁₀ to give five-membered ring (5MR) and four-membered ring (4MR), respectively. In the gas phase, TS(5MR) depicts a dihedral of 10.13°, but the explicit solvation increase it to 13.6° and 18.5° in benzene and toluene, respectively, however, in the implicit solvent, the dihedral change is smaller. In the case of TS(4MR), the gas phase dihedral is -30.6° and the implicit solvent increases the value to -28° and - 29.5° in benzene and toluene respectively. The C(7)–C(14) distance is shorter in explicit toluene, whereas the C(7)–C(10) distance is in explicit benzene. These observations correlate with the preference of the reaction in each implicit solvent.

FIGURE 4

FIGURE 4. Conformation Newman projections of the transition states in different solvent scenarios. Upper and bottom lines show TS(5MR) and TS(4MR) conformers respectively. In these projections, C₁₄ is the proximal atom, while C₇ is the distal one. (A) and (B) without cation in implicit benzene; (C) and (D) without cation in implicit toluene; (E) and (F) two possible conformation of 1 with two lithium cations in gas phase; (G) and (H) in gas phase; (I) and (J) in implicit benzene; (K) and (L) in implicit toluene; (M) and (N) in explicit benzene; (O) and (P) in explicit toluene.

FIGURE 5

FIGURE 5. Intermolecular interactions within the benzene solvation sphere of TS(4MR). Four solvent molecules of the second solvation sphere are included.

3.4 Cyclization reactions within the solvent models

In Figure 6, we present the internal reaction coordinate (IRC) electronic energy profiles of the cyclization of the deprotonated 1 with two lithium-ions linked to the N-extreme of the nitrile group and other in the oxygen of the oxirane. The free energy and enthalpic profiles are presented in the supplementary information. We focused on electronic energy because in the next section this energy is partitioned into solute and solvent contributions, to understand the role of the solvent in the selectivity. This figure shows that in the gas phase, 5MR is thermodynamically more favored than 4MR, whereas both products present similar reaction barriers, with a difference of less than 4.18 kJ/mol. We found the same behavior when the implicit benzene or toluene model is presented. Based on these results, the experimental observation should show just the presence of 5MR. However, the main product is 4MR when benzene is used as a solvent and 5MR with toluene. This contradiction makes clear the necessity to include a better solvent model to improve the description of the cyclization reaction. We found that the inclusion of an explicit solvent is a determinant factor in reproducing the experimental stabilities and barriers. In this way, the observed energetic differences between 5MR and 4MR are linked to specific interactions between them and the solvent.

FIGURE 6

FIGURE 6. Reaction energy profiles in gas phase, implicit and explicit solvent. Benzene right and Toluene left. In the gas phase and implicit solvent, 5MR is thermodynamically favored and both products present similar reaction barriers. Including explicit solvent reduces the 5MR stability and increases its reaction barrier.

Figure 6 shows the effect of the explicit solvation on the energetic profile of both cyclization paths. In benzene, the barrier for 5MR is 25.73 kJ/mol larger than for 4MR, whereas in toluene the height of both barriers is almost equal but 5MR is more stable. In this way, 4MR is kinetically favored in benzene, while 5MR is thermodynamically favored in toluene. The presence of the explicit solvent modifies the energetic profile to reproduce the experimental results.

3.5 Solvent and solute contribution to the energy

We used the Quantum Theory of Atoms in Molecules (QTAIM) (Bader, 1994) to obtain the solvent (E_solute) and solute (E_solvent) contributions to the total electronic energy. Each contribution comprises the sum of the atomic energy of every atom in the solvent and the solute. The atomic energies were calculated using the partition of the molecular space defined by the gradient of electron density. Eq. 3 describes the relationship between each contribution and the total energy.

$E_{\mathit{total}}=E_{\mathit{solute}}+E_{\mathit{solvent}}=\sum _{A\in Solute}E\left(A\right)+\sum _{B\in Solvent}E\left(B\right)(3)$

Figure 7 shows the solute and solvent contributions to the energy of the reactant, transition state, and product. In this way, in benzene, the barrier to 4MR, is originated by the destabilization of 1, by 90.58 kJ/mol, along with a stabilization of the solvent, –6.36 kJ/mol. In the case of toluene for the same product, both contributions are destabilizing, the solvent by 68.49 kJ/mol, and 1 by 25.19 kJ/mol. On the other hand, 5MR destabilizes in both solvents whereas the solvent stabilizes. This behavior is more notorious in the case of toluene. One can observe that the solute in explicit solvent behaves opposite to the molecule in the gas phase or in implicit solvent. The contribution of the solvent modulates the energetic behavior of the solute, which is a determinant to reproduce the experimental observations.

FIGURE 7

FIGURE 7. Solute and solvent energetic contributions to the reaction profile in the explicit solvent environment. Benzene right and Toluene left. The solvent modulates the energetic behavior of the solute. In both solvents T5(5MR) is determined by the solute contribution, whereas TS(4MR) is determined by the solute in benzene and by the solvent in toluene.

3.6 Specific solute-solvent interactions

The analysis of the specific solute-solvent interactions allows us to determine the role of each one in the energy profiles. The interaction energies were estimated with the Espinosa–Molins–Lecomte (EML) approximation. (Espinosa et al., 1998) The EML approximation is based on the value of the potential energy density (V(r)) at the bond critical points (BCP) that links atoms of the solute with atoms of the solvent molecules. The value of V(r) is proportional to the interatomic interaction energy. The sum of interaction energies associated with every solute-solvent BCP gives the total interaction energy in which k is a constant that depends on each interaction; for hydrogen bonds, the value is 0.5, which was determined by empirical correlations between the interaction energies and BCP properties of a set of hydrogen bonds. In a previous job, our group found, studying the ligand-DNA interaction, that k = 0.433 describes non-hydrogen bond weak interactions. (Díaz-Gómez et al., 2017) The energetic contributions of the specific interactions are around 4–29 kJ/mol, table 3 shows the type of weak interactions we found in the systems and the total contribution of each type in each solvent. The experimental facts are supported by the interaction energy contribution since the solute-solvent interactions give a larger stabilization in one of the two solvents, 5MR is more stabilized in toluene, while 4MR is in benzene.

TABLE 3

TABLE 3. Energetic contribution (kJ/mole) of the specific interactions to the intermolecular interaction energy between transition state structures with the explicit solvent. The times the interaction occurs are in parentheses.

The breakdown of solute-solvent contributions indicates an important contribution of the specific interactions on the barrier to reaching a transition state. In all cases, the C–H⋯π interactions are present in greater numbers, which is to be expected due to the aromatic nature of the solvents, so they have a significant stabilizing contribution. Subsequently, the hydrogen-hydrogen bonds are presented in an average of 16 interactions and represent a net stabilizing contribution. (Matta et al., 2003) In benzene, Li⋯C interactions favor TS(5MR), whereas O⋯H does for TS(4MR).

4 Conclusion

In this work, we answered one fundamental question in chemical reactivity, how many solvent molecules should be added and where they should be placed to understand the role of the explicit solvent in the reaction selectivity, in the specific case of the cyclization of δ-epoxy nitriles. For this purpose, we performed the experimental determination of the first solvent sphere around a molecule using NMR methods. We found that an epoxy nitrile is surrounded by five benzene or four toluene molecules. With this experimental information, we determined the role of the subsequent solvent spheres ( ∼20 solvent molecules) in locating the first sphere around the solute to reproduce the NMR spectra. We analyzed the specific interactions between the solute with the solvent spheres and how they control the conformation of the reactant and reaction barriers. We showed that explicit solvation is necessary to reproduce experimental results of the cyclization reaction. Also, a new approach to solvation energy is presented, which includes solvent cavity formation and specific solute-solvent interactions from explicit solvation.

Data availability statement

The original contributions presented in the study are included in the article/Supplementary Material, further inquiries can be directed to the corresponding author.

Author contributions

JH-L: Synthesis, NMR analysis, DFT Calculations, and Visualization, Writing; KR-G: NMR titration and solvent sphere determination; BQ-G: NMR Characterization and analysis; ALS-P: Synthesis and characterization; EC-N: Calorometric characterization of solvent spheres; FC-G: Conceptualization, Methodology, Formal analysis and Writing.

Funding

Authors acknowledge CONACyT-México (Grants CB 220392, 293294, and 1561802), and PAPIIT/DGAPA/UNAM (Grants IN202717, IN212520, and IN207822) for financial support, and DGTIC-UNAM (LANCAD-UNAM-DGTIC-194) for computer time.

Acknowledgments

The authors thank B. Meza Gonzalez, G. E. Cortés Romero, J. D. Vázquez Cuevas, M. M. Aguilar Araiza and Y. A. García Hernández for their help with computational systems.

Conflict of interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher’s note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Supplementary material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fchem.2022.1012769/full#supplementary-material

References

Alabugin, I. (2016). Stereoelectronic effects: A bridge between structure and reactivity. Chichester, UK: John Wiley & Sons. doi 10, 9781118906378.

Google Scholar

Avila-Zárraga, J. G., and Maldonado, L. A. (2000). 2-aryl-5-epoxynitrile anion cyclizations: Total syntheses of (±)-α-cuparenone and (±)-epilaurene. Chem. Lett. 29, 512–513. doi:10.1246/cl.2000.512

CrossRef Full Text | Google Scholar

Bader, R. (1994). Atoms in molecules: A Quantum theory. International ser. Of monogr. On chem. Clarendon Press.

Google Scholar

Baldwin, J. E. (1976). Rules for ring closure. J. Chem. Soc. Chem. Commun. 1976, 734–736. doi:10.1039/c39760000734

CrossRef Full Text | Google Scholar

Bayse, C. A., and Antony, S. (2009). Modeling the oxidation of ebselen and other organoselenium compounds using explicit solvent networks. J. Phys. Chem. A 113, 5780–5785. doi:10.1021/jp901880n

PubMed Abstract | CrossRef Full Text | Google Scholar

Bayse, C. A. (2007). Dft study of the glutathione peroxidase-like activity of phenylselenol incorporating solvent-assisted proton exchange. J. Phys. Chem. A 111, 9070–9075. doi:10.1021/jp072297u

PubMed Abstract | CrossRef Full Text | Google Scholar

Bayse, C. A. (2010). Model mechanisms of sulfhydryl oxidation by methyl-and benzeneseleninic acid, inhibitors of zinc-finger transcription factors. J. Inorg. Biochem. 104, 1–8. doi:10.1016/j.jinorgbio.2009.09.007

PubMed Abstract | CrossRef Full Text | Google Scholar

Bayse, C. A. (2011). Transition states for cysteine redox processes modeled by dft and solvent-assisted proton exchange. Org. Biomol. Chem. 9, 4748–4751. doi:10.1039/c1ob05497j

PubMed Abstract | CrossRef Full Text | Google Scholar

Boereboom, J. M., Fleurat-Lessard, P., and Bulo, R. E. (2018). Explicit solvation matters: Performance of qm/mm solvation models in nucleophilic addition. J. Chem. Theory Comput. 14, 1841–1852. doi:10.1021/acs.jctc.7b01206

PubMed Abstract | CrossRef Full Text | Google Scholar

Bulo, R. E., Ensing, B., Sikkema, J., and Visscher, L. (2009). Toward a practical method for adaptive qm/mm simulations. J. Chem. Theory Comput. 5, 2212–2221. doi:10.1021/ct900148e

PubMed Abstract | CrossRef Full Text | Google Scholar

Cainelli, G., Galletti, P., and Giacomini, D. (2009). Solvent effects on stereoselectivity: More than just an environment. Chem. Soc. Rev. 38, 990–1001. doi:10.1039/b802815j

PubMed Abstract | CrossRef Full Text | Google Scholar

Chung, L. W., Sameera, W., Ramozzi, R., Page, A. J., Hatanaka, M., Petrova, G. P., et al. (2015). The oniom method and its applications. Chem. Rev. 115, 5678–5796. doi:10.1021/cr5004419

PubMed Abstract | CrossRef Full Text | Google Scholar

Corset, J., Castellà-Ventura, M., Froment, F., Strzalko, T., and Wartski, L. (2003). Study of the lithiated phenylacetonitrile monoanions and dianions formed according to the lithiated base used (lhmds, lda, or n-buli). 1. evidence of heterodimer (“quadac”) or dianion formation by vibrational spectroscopy. J. Org. Chem. 68, 3902–3911. doi:10.1021/jo020492t

PubMed Abstract | CrossRef Full Text | Google Scholar

Cossi, M., Rega, N., Scalmani, G., and Barone, V. (2003). Energies, structures, and electronic properties of molecules in solution with the c-pcm solvation model. J. Comput. Chem. 24, 669–681. doi:10.1002/jcc.10189

PubMed Abstract | CrossRef Full Text | Google Scholar

da Silva, E. F., Svendsen, H. F., and Merz, K. M. (2009). Explicitly representing the solvation shell in continuum solvent calculations. J. Phys. Chem. A 113, 6404–6409. doi:10.1021/jp809712y

PubMed Abstract | CrossRef Full Text | Google Scholar

Díaz-Gómez, D. G., Galindo-Murillo, R., and Cortés-Guzmán, F. (2017). The role of the dna backbone in minor-groove ligand binding. ChemPhysChem 18, 1909–1915. doi:10.1002/cphc.201700260

PubMed Abstract | CrossRef Full Text | Google Scholar

Ditchfield, R. (1974). Self-consistent perturbation theory of diamagnetism: I. A gauge-invariant lcao method for nmr chemical shifts. Mol. Phys. 27, 789–807. doi:10.1080/00268977400100711

CrossRef Full Text | Google Scholar

Espinosa, E., Molins, E., and Lecomte, C. (1998). Hydrogen bond strengths revealed by topological analyses of experimentally observed electron densities. Chem. Phys. Lett. 285, 170–173. doi:10.1016/s0009-2614(98)00036-0

CrossRef Full Text | Google Scholar

Fleming, F. F., Wei, Y., Liu, W., and Zhang, Z. (2008). Metalated nitriles: Stereodivergent cation-controlled cyclizations. Tetrahedron 64, 7477–7488. doi:10.1016/j.tet.2008.05.110

PubMed Abstract | CrossRef Full Text | Google Scholar

Floris, F. M., Tomasi, J., and Ahuir, J. P. (1991). Dispersion and repulsion contributions to the solvation energy: Refinements to a simple computational model in the continuum approximation. J. Comput. Chem. 12, 784–791. doi:10.1002/jcc.540120703

CrossRef Full Text | Google Scholar

Frisch, M. J., Trucks, G. W., Schlegel, H. B., Scuseria, G. E., Robb, M. A., Cheeseman, J. R., et al. (2009). Gaussian09 revision E.01. Wallingford CT: Gaussian Inc. [Dataset].

Google Scholar

Grande-Aztatzi, R., Rivillas-Acevedo, L., Quintanar, L., and Vela, A. (2013). Structural models for cu (ii) bound to the fragment 92–96 of the human prion protein. J. Phys. Chem. B 117, 789–799. doi:10.1021/jp310000h

PubMed Abstract | CrossRef Full Text | Google Scholar

Hansen, T., Roozee, J. C., Bickelhaupt, F. M., and Hamlin, T. A. (2022). How solvation influences the sn2 versus e2 competition. J. Org. Chem. 87, 1805–1813. doi:10.1021/acs.joc.1c02354

PubMed Abstract | CrossRef Full Text | Google Scholar

Harvey, D. T. (1999). Modern analytical chemistry. 1 edn. McGraw-Hill.

Google Scholar

Kamerlin, S. C., Haranczyk, M., and Warshel, A. (2009). Are mixed explicit/implicit solvation models reliable for studying phosphate hydrolysis? A comparative study of continuum, explicit and mixed solvation models. ChemPhysChem 10, 1125–1134. doi:10.1002/cphc.200800753

PubMed Abstract | CrossRef Full Text | Google Scholar

Keith, T. A. (2019). Aimall (version 19.10.12. Overland Park KS, USA: TK Gristmill Software. (aim.tkgristmill.com) [Dataset].

Google Scholar

Kelly, C. P., Cramer, C. J., and Truhlar, D. G. (2006). Adding explicit solvent molecules to continuum solvent calculations for the calculation of aqueous acid dissociation constants. J. Phys. Chem. A 110, 2493–2499. doi:10.1021/jp055336f

PubMed Abstract | CrossRef Full Text | Google Scholar

Kononov, L. (2015). Chemical reactivity and solution structure: On the way to a paradigm shift? RSC Adv. 5, 46718–46734. doi:10.1039/c4ra17257d

CrossRef Full Text | Google Scholar

Llanio-Trujillo, J., Marques, J., and Pereira, F. (2011). An evolutionary algorithm for the global optimization of molecular clusters: Application to water, benzene, and benzene cation. J. Phys. Chem. A 115, 2130–2138. doi:10.1021/jp1117695

PubMed Abstract | CrossRef Full Text | Google Scholar

Lujan-Montelongo, J. A., Vazquez-Sanchez, A., and Avila-Zarraga, J. G. (2009). Regio and stereo selectivity on the α-aryl-δ-epoxynitrile anionic cyclization reactions: Modulation towards a 5-endo or 4-exo process. Heterocycles 78, 1955–1976. doi:10.3987/com-09-11669

CrossRef Full Text | Google Scholar

Matta, C. F., Hernández-Trujillo, J., Tang, T.-H., and Bader, R. F. (2003). Hydrogen–hydrogen bonding: A stabilizing interaction in molecules and crystals. Chem. Eur. J. 9, 1940–1951. doi:10.1002/chem.200204626

PubMed Abstract | CrossRef Full Text | Google Scholar

Pedretti, A., Villa, L., and Vistoli, G. (2002). Vega: A versatile program to convert, handle and visualize molecular structure on windows-based pcs. J. Mol. Graph. Model. 21, 47–49. doi:10.1016/s1093-3263(02)00123-7

PubMed Abstract | CrossRef Full Text | Google Scholar

Pure, I. U., and Chemistry, A. (2009). Iupac compendium of chemical terminology – the gold book. [Dataset].

Google Scholar

Ramirez-Gualito, K., Alonso-Rios, R., Quiroz-Garcia, B., Rojas-Aguilar, A., Diaz, D., Jimenez-Barbero, J., et al. (2009). Enthalpic nature of the ch/π interaction involved in the recognition of carbohydrates by aromatic compounds, confirmed by a novel interplay of nmr, calorimetry, and theoretical calculations. J. Am. Chem. Soc. 131, 18129–18138. doi:10.1021/ja903950t

PubMed Abstract | CrossRef Full Text | Google Scholar

Reichardt, C., and Welton, T. (2011). Solvents and solvent effects in organic chemistry. John Wiley & Sons.

Google Scholar

Schienbein, P., Schwaab, G., Forbert, H., Havenith, M., and Marx, D. (2017). Correlations in the solute–solvent dynamics reach beyond the first hydration shell of ions. J. Phys. Chem. Lett. 8, 2373–2380. doi:10.1021/acs.jpclett.7b00713

PubMed Abstract | CrossRef Full Text | Google Scholar

Schmitt, M., Böhm, M., Ratzer, C., Vu, C., Kalkman, I., and Meerts, W. L. (2005). Structural selection by microsolvation: Conformational locking of tryptamine. J. Am. Chem. Soc. 127, 10356–10364. doi:10.1021/ja0522377

PubMed Abstract | CrossRef Full Text | Google Scholar

Skyner, R. E., Groom, C. R., Van Mourik, T., McDonagh, J. L., and Mitchell, J. B. (2015). A review of methods for the calculation of solution free energies and the modelling of systems in solution. Phys. Chem. Chem. Phys. 17, 6174–6191. doi:10.1039/c5cp00288e

PubMed Abstract | CrossRef Full Text | Google Scholar

Stork, G., Cama, L. D., and Coulson, D. (1974). Epoxynitrile cyclization. general method of ring formation. J. Am. Chem. Soc. 96, 5268–5270. doi:10.1021/ja00823a052

CrossRef Full Text | Google Scholar

Stork, G., and Cohen, J. F. (1974). Ring size in epoxynitrile cyclization. general synthesis of functionally substituted cyclobutanes. application to (+)grandisol. J. Am. Chem. Soc. 96, 5270–5272. doi:10.1021/ja00823a053

CrossRef Full Text | Google Scholar

Strzalko, T., Wartski, L., Corset, J., Castellà-Ventura, M., and Froment, F. (2012). Study of the lithiated phenylacetonitrile monoanions and dianions formed according to the lithiated base used (lhmds, lda, or n-buli). 2. alkylation and deuteriation mechanism study by vibrational and nmr spectroscopy and quantum chemistry calculations. J. Org. Chem. 77, 6431–6442. doi:10.1021/jo300758g

PubMed Abstract | CrossRef Full Text | Google Scholar

Sunoj, R. B., and Anand, M. (2012). Microsolvated transition state models for improved insight into chemical properties and reaction mechanisms. Phys. Chem. Chem. Phys. 14, 12715–12736. doi:10.1039/c2cp41719g

PubMed Abstract | CrossRef Full Text | Google Scholar

Tantillo, D. J. (2018). Applied Theoretical Organic Chemistry. doi:10.1142/q0119

CrossRef Full Text

Tomasi, J., Mennucci, B., and Cammi, R. (2005). Quantum mechanical continuum solvation models. Chem. Rev. 105, 2999–3094. doi:10.1021/cr9904009

PubMed Abstract | CrossRef Full Text | Google Scholar

Vandenbussche, S., Díaz, D., Fernández-Alonso, M. C., Pan, W., Vincent, S. P., Cuevas, G., et al. (2008). Aromatic–carbohydrate interactions: An nmr and computational study of model systems. Chem. Eur. J. 14, 7570–7578. doi:10.1002/chem.200800247

PubMed Abstract | CrossRef Full Text | Google Scholar

Vogel, A. I., Mendham, J., Denney, R. C., and Barnes, J. D. (2000). Vogel’s Textbook of Quantitative Chemical Analysis. Prentice Hall.

Google Scholar

Wang, C., Ying, F., Wu, W., and Mo, Y. (2014). How solvent influences the anomeric effect: Roles of hyperconjugative versus steric interactions on the conformational preference. J. Org. Chem. 79, 1571–1581. doi:10.1021/jo402306e

PubMed Abstract | CrossRef Full Text | Google Scholar

Yang, X., and Fleming, F. F. (2017). C- and n-metalated nitriles: The relationship between structure and selectivity. Acc. Chem. Res. 50, 2556–2568. doi:10.1021/acs.accounts.7b00329

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhao, Y., and Truhlar, D. G. (2008). The m06 suite of density functionals for main group thermochemistry, thermochemical kinetics, noncovalent interactions, excited states, and transition elements: Two new functionals and systematic testing of four m06-class functionals and 12 other functionals. Theor. Chem. Acc. 120, 215–241. doi:10.1007/s00214-007-0310-x

CrossRef Full Text | Google Scholar