The Electronic Circular Dichroism of Nicotine in Aqueous Solution: A

Jan 8, 2015 - Test Case for Continuum and Mixed Explicit-Continuum Solvation. Approaches. Franco Egidi. Scuola Normale Superiore, Piazza dei Cavalieri...
0 downloads 0 Views 2MB Size
Article pubs.acs.org/JPCA

The Electronic Circular Dichroism of Nicotine in Aqueous Solution: A Test Case for Continuum and Mixed Explicit-Continuum Solvation Approaches Franco Egidi Scuola Normale Superiore, Piazza dei Cavalieri 7, I-56126 Pisa, Italy

Rosario Russo Dipartimento di Chimica e Chimica Industriale, Università di Pisa, via Giuseppe Moruzzi 3, I-56124 Pisa, Italy

Ivan Carnimeo Compunet, Istituto Italiano di Tecnologia (IIT), via Morego 30, I-16163 Genova, Italy Scuola Normale Superiore, Piazza dei Cavalieri 7, I-56126 Pisa, Italy

Alessandro D’Urso Department of Chemical Sciences, University of Catania, Viale A. Doria 6, I-95125 Catania, Italy

Giordano Mancini Scuola Normale Superiore, Piazza dei Cavalieri 7, I-56126 Pisa, Italy

Chiara Cappelli* Dipartimento di Chimica e Chimica Industriale, Università di Pisa, via Giuseppe Moruzzi 3, I-56124 Pisa, Italy ABSTRACT: In this paper, we extend an integrated QM/MM/polarizable continuum model (PCM) method, which combines a fluctuating charge (FQ) approach to the MM polarization with the PCM, to describe electronic circular dichroism (ECD) spectra of systems in aqueous solution. The main features of the approach are presented, and then applications to the UV and ECD spectra of neutral (S)-nicotine in aqueous solution are reported. The performance of the QM/FQ/PCM is compared with that of the PCM against newly measured UV ECD spectra, which are in agreement with previous findings. The inclusion of specific solvation effects via the QM/FQ/PCM method leads to an improvement in the calculated spectra compared to the experimental findings, though the pure PCM results are still qualitatively correct and are a useful tool for the characterization of the states.

I. INTRODUCTION The simulation of chiroptical properties by means of quantummechanical (QM) methods1−4 has been demonstrated to yield accurate results and to be a powerful methodology to assist the determination of the molecular absolute configuration.5 However, a reliable comparison between experimental data and computations requires a careful consideration of the environment. In fact, chiroptical properties and spectroscopies are usually measured for systems in the condensed phase (i.e., solutions or pure liquids), and the presence of the solvating © XXXX American Chemical Society

environment may affect not only the absolute value but especially the sign of the chiroptical response.6,7 Solvent effects on chiroptical properties and spectroscopies6,8 can originate from both polarization effects and supramolecular Special Issue: Jacopo Tomasi Festschrift Received: October 20, 2014 Revised: January 7, 2015

A

DOI: 10.1021/jp510542x J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A

applied it to (S)-nicotine, which we have already exploited as reference system for the development of protocols rooted in the PCM for calculating OR, infrared, and vibrational circular dichroism (VCD) spectra38 and NMR shielding constants.39 In particular, the QM/FQ/PCM combines a fluctuating charge (FQ) approach to the MM polarization with the PCM,34−37 and it has been shown to give exceptional agreement with experimental findings in the case of the OR of methyloxirane in aqueous solution.16 In this model, polarization effects are considered by endowing each atom with a charge that depends on the atomic electronegativity and hardness.34 The equilibrium between the charge interactions and the differences in the atomic electronegativities are the physical basis of the electronegativity equalization principle,40−42 which defines the FQ model. The QM/FQ/PCM method has recently been extended to the evaluation of analytical first and second derivatives 36 and response equations35 and also to magnetic perturbations with GIAOs37 and excitation energies at the EOM-CCSD level.43 The paper is organized as follows. In the next section, the basic theory of the calculation of ECD spectra within the fully polarizable QM/FQ/PCM approach is outlined. After some details on the computational and experimental segment of this study are provided, the results of the application of PCM44,45 and QM/FQ/PCM to the ECD of (S)-nicotine are discussed and compared with newly measured ECD spectra. The article ends with a brief summary and the relevant conclusions.

clustering patterns, which cannot be neglected whenever strong interactions, such as hydrogen bonding, dominate the solute− solvent interaction. The extent of such effects, however, strongly depends on both the property and the solute−solvent couple, which is of particular relevance in the case of aqueous solutions. A well-known case is the optical rotation (OR) of methyloxirane in aqueous solution,9−15 for which an adequate description of both the absolute value and sign requires careful modeling of the solvation shell around the solute.16 On the contrary, the same property of the same molecule but in different solvents, spanning a large range of polarities, is accurately described by continuum solvation methodologies.9 The OR of nicotine has also been shown to be very much solvent-dependent,17 and given the close relationship between the OR and ECD molecular response, it is worth addressing the problem of modeling solvation effects in the ECD spectrum of this system. In fact, ECD spectra can experience relevant changes upon solvation, especially due to solvent-induced structural changes18 or helicity inversion.19 The OR and other chiroptical spectroscopies of various systems in water (e.g., (S)-N-acetylproline amide)9,20 can be accurately reproduced by implicit solvation models. Therefore, the question of what is the best way to model solvent effects on chiroptical properties of systems in aqueous solution is far from being finally assessed, especially because of the fact that although implicit solvation models have been tested thoroughly for various kinds of chiroptical properties and spectroscopies,7,8,21,22 the same cannot be said about explicit or combined explicit/implicit models.16,21,23−26 Solvation involves different interactions and effects that are not easily modeled by a unique theory. This complexity is reflected in the variety of solvation models that have been proposed. On one hand, continuum solvation models have been shown to represent in many cases an efficient and reliable approach to chiroptical properties and spectroscopies,20,27−29 especially if additional aspects (such as nonequilibrium and/or local-field effects27,30,31) not commonly included in the basic versions of the models are properly taken into account. However, as anticipated above, whenever specific solute− solvent interactions (or solute−solute aggregation effects)27 are relevant, a purely continuum approach may fail. Refinements of continuum models by means of an extension of the definition of the solute have been proposed,32 but they cause an increase in both the computational cost and the complexity of the calculation that is often not accompanied by a definite improvement in the reliability of the calculations, especially whenever chiral imprinting in the solvation shell can be invoked.23,33 In this case, standard QM/molecular mechanics (MM) or supermolecular approaches can also be inadequate because the imprinting is in fact a dynamic and local effect that requires the coupling of accurate electronic methods to evaluate the property with equally accurate dynamic simulations to obtain a reliable description of the fluctuations in the solute− solvent interactions, which in this case do not average out to result in an isotropic effect but instead sustain the chiral response in the solvent cage around the chiral solute. On the basis of these premises and with the aim of continuing our investigation of the performance of various kinds of solvation approaches to describe the chiroptical response of solvated systems, in this contribution we have extended an integrated QM/MM/PCM method we recently developed, named QM/FQ/PCM,34−37 to ECD spectra and

II. THEORY In this section we will recall the main aspects of the TDDFT/ FQ/PCM approach, in view of the application to ECD spectroscopy. The interested reader can find further details in refs 34−37. In the QM/FQ/PCM approach, the total ground state energy of the system E0 can be written as int int E0 = E QM + EPCM + EQM/PCM + EFQ + EQM/FQ int + E FQ/PCM

(1)

The QM/FQ interaction energy, Eint QM/FQ, is the electrostatic interaction between the set of Nq FQ charges qk and the electrostatic potential Vk generated by the nuclear and electronic QM charge densities, i.e. Nq int EQM/FQ = +∑ Vkqk

(2)

k

Such a term is analogous to Eint QM/PCM, where the interaction is between Vk and the set of Nσ PCM charges σk: Nσ int EQM/PCM = +∑ Vkσk

(3)

k

The FQ/PCM interaction energy can be written as Nq int E FQ/PCM



= +∑ ∑ j

k

qjσk |r j − s k |

(4)

where rj and sk are the positions of the FQ and PCM charges, respectively. Then E0 is a functional of the QM density matrix and the FQ and PCM charges. The total energy of the FQ subsystem is B

DOI: 10.1021/jp510542x J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A ,=

1 qJq + qχ + qV 2

(5)

where χ is the array containing the atomic electronegativities, V is the electrostatic potential including the interactions with the PCM and the QM subsystem, and J is the interaction matrix between the charges, which in our implementation is based on the Ohno functional form (see ref 34 and references therein). By exploitation of the electronegativity equalization principle, the optimal set of charges that minimizes such an energy functional is obtained as

q = J − 1 ( − χ − V)

(6)

Thus, at each step of the SCF procedure, the FQ and PCM charges are obtained once the electrostatic potential on the coordinates of the FQ atoms and on the PCM tesserae has been evaluated. The calculation of the excited-state transition energies and amplitudes within the QM/FQ/PCM model is carried out by applying time-dependent linear response theory to the QM/ FQ/PCM Hamiltonian.35 Once the explicit FQ and PCM contributions are added to the response matrices, the TDDFT equations are solved in the Casida formalism46 using standard computational procedures. Then from the excitation energies and the transition amplitudes, the rotatory strength tensor is computed in the velocity form as described in ref 47.

Figure 1. A and B conformers of nicotine. The points on the PES used for the RESP fitting were generated by changing the value of the dihedral angle δ1. The nitrogen atoms N1 and N2 are also displayed.

between conformation i and the global minimum, kB is the Boltzmann constant, and T is the absolute temperature (298.15 K). Dispersion and repulsion terms were modeled according to the general Amber force field (GAFF).56 Classical MD simulations were performed with the GROMACS 4.6.5 software.57 A single nicotine molecule was solvated in a cubic box with a side length of 36.3 Å containing 1545 water molecules modeled using the TIP4P/200558 parameter set. Electrostatic interactions were taken into account by means of the particle mesh Ewald method59 using a cutoff radius of 12.0 Å in real space; a cutoff radius of 12.0 Å was also chosen for the van der Waals interactions. Two short (200 ps) NVT simulations were performed at 150.0 and 298.15 K for thermalization purposes. Then a 10.0 ns trajectory was obtained using the velocity-rescale method60 with a coupling constant of 0.1 ps, saving coordinates every 0.5 ps. A total of 300 snapshots with a constant stride of 30 ps were extracted for FQ calculations. Volumetric maps were produced with VMD.61 We then performed a TDDFT calculation on each snapshot, treating (S)-nicotine at the CAM-B3LYP/aug-cc-pVDZ level of theory and the solvent molecules using the FQ/PCM polarizable embedding scheme.34−36 The FQ calculations were executed using the SPC parameters of Rick et al.41 We extracted the excitation energies and spectroscopic intensities pertaining to the first 12 excited states from each snapshot and then plotted the resulting UV and ECD spectra after convoluting each peak using Gaussian functions with halfwidths at half-maximum (HWHMs) of 10 and 8 nm, respectively. We also measured the experimental UV absorption and ECD spectra of neutral (S)-nicotine in aqueous solution. We diluted 1 μL of neat product in 25 mL of ultrapure water obtained from an Elga Purelab Flex system by Veolia in order to obtain a 250 μM stock solution. Then we diluted 120 μL of (S)-nicotine stock solution with 3 mL of ultrapure water in a Hellma quartz cuvette of 1 cm path length. The concentration of the resulting solution was 10 μM. In order to investigate the neutral form of

III. COMPUTATIONAL AND EXPERIMENTAL DETAILS All of the QM/PCM and QM/FQ/PCM calculations were performed using a locally modified version of the Gaussian 09 quantum chemistry suite of programs,48 employing the CAMB3LYP functional,49 which has been previously shown to perform better than other functionals in the description of the excited states of nicotine in vacuo,50 and the aug-cc-pVDZ basis set.51 UV absorption intensities were calculated from the oscillator strengths computed in the length gauge, while ECD intensities were derived from the velocity-gauge rotatory strengths to ensure origin invariance in the results.47 PCM calculations44,52 were performed by building the molecular cavity in terms of a set of interlocking spheres centered on the atoms, with radii according to the Gaussian 09 defaults (UFF radii). Non-electrostatic effects in the free-energy calculations used to estimate conformer populations were taken into account by using the SMD parametrization.53 To parametrize a reliable charge set for the nicotine molecule, we performed a scan on the potential energy surface along the δ1 coordinate (see Figure 1) at the B3LYP/N07D54 level, which was previously shown to be appropriate to model the ground-state properties of this system,38 relaxing the remaining degrees of freedom at each step. A total of 22 points were calculated, varying the value of δ1 from −180° to 180°, where δ1 = 135° corresponds to the so-called A conformer and δ1 = −45° corresponds to the B conformer (see Figure 1). Atomic charges were then calculated by applying the multiconformation RESP (MultiRESP55) procedure, which takes into account more than one conformation of minimum energy via the following equation: C(n) =

∑i Ci(n)e−ΔEi / kBT ∑i e−ΔEi / kBT

where C(n) is the MultiRESP charge of atom n, Ci(n) is the charge of atom n in conformation i, ΔEi is the energy difference C

DOI: 10.1021/jp510542x J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A

Table I. CAM-B3LYP/aug-cc-pVDZ/PCM Vertical Excitation Wavelengths (λ, in nm), Oscillator Strengths (f), and Rotatory Strengths (R, in 10−40 erg esu cm/gauss) for Nicotine in Water nicotine A

nicotine B

type

λ

f

R

type

λ

f

R

n−πQ* CTQ/La CTQ/La n−π*AQ CTQ/Ryd CTAQ Lb CTAQ CTAQ/Ryd CTAQ/Ryd πQ−Ryd n−Ryd

240 231 228 218 212 209 200 197 193 188 184 179

0.0055 0.0264 0.0577 0.0004 0.0343 0.0113 0.0704 0.0844 0.0170 0.0084 0.0095 0.2363

−15.09 1.64 10.74 0.41 −22.03 6.96 −7.06 0.13 −5.36 5.99 −24.09 39.71

n−πQ* La CTQ/La CTQ/La n−π*AQ Ryd Lb CTAQ CTAQ/Ryd CTAQ/Ryd πQ−Ryd n−Ryd

242 229 228 217 212 211 201 198 194 190 185 180

0.0081 0.0708 0.0069 0.0042 0.0031 0.0222 0.0650 0.1113 0.0211 0.0161 0.0123 0.1738

−23.83 34.64 3.61 −12.33 −1.70 −6.46 15.52 −12.87 −4.20 −24.53 21.89 −10.80

quinoid or antiquinoid) and charge-transfer (CT) transitions between the lone-pair n′ orbital on the pyrrolidine ring and the π* orbitals on the pyridine ring, denoted CTQ and CTAQ. Upon analyzing the results of the TDDFT calculations, we observed * states that result from mixing of the πQ−πQ* and πAQ−πAQ transitions with comparable weights, and we simply denote these excitations with the symbol La, whereas we use the * and symbol Lb to denote excitations of the type πQ−πAQ πAQ−πQ*, which also show some degree of mixing. Table I summarizes the result of the TDDFT/PCM calculations for the first eight transitions of (S)-nicotine in water. The lowest-energy transition is n−πQ* for both conformations, followed by La, which shows some degree of mixing with the quinoid CT transition for both conformations. It should be noted that the La excitation, as expected, has a much higher oscillator strength than n−π*Q, but its rotatory strength is not as dominant and for the B conformation is of opposite sign. At higher excitation energies we find the n−π*AQ transition, which shows negligible oscillator strength (the corresponding transition in pyridine is symmetry-forbidden) and rotatory strength and is therefore is unlikely to contribute significantly to either the absorption or ECD spectrum. We then find the Lb transition, more charge-transfer excitations, and Rydberg excitations. It should be noted that although the two conformations show very similar behavior in terms of the oscillator strengths of similar transitions, this is not true for the rotatory strengths, as the Lb transition, for instance, is negative for the A conformation and strongly positive for the B conformation. In addition, while the B conformation is predicted to be less populated by our calculations, its rotatory strengths are generally higher compared with its sibling, and it is therefore expected to dominate the ECD spectrum nonetheless. Figure 2 shows the convoluted UV−vis spectrum of nicotine. It can be seen that the two conformations show very similar spectra, mainly composed of two bands: the first at 230 nm is assigned to an La transition localized on the pyridine ring, and the band of higher intensity at 200 nm is composed of CT between the two rings and the Lb transition. In agreement with previous investigations,62,63 the experimental spectrum shows only one band at 260 nm, which is much red-shifted with respect to the calculated value. It has in fact been shown50 that CAM-B3LYP greatly overestimates the excitation energies of this transition with respect to the other ones, though it still performs better than other functionals compared with

nicotine, micromolar amounts of 1 M NaOH were added to the cuvette to obtain a solution pH of 10, which was basic enough to ensure that the system was present in its unprotonated form.62,63 UV and ECD spectra were recorded over the wavelength range from 190 to 320 nm using a JASCO UVV530 spectrophotometer and a JASCO J-810 spectropolarimeter, respectively. The ECD spectra were recorded with a wavelength data pitch of 0.5 nm, a response time of 2 s, and a scanning speed of 20 nm/min. The spectra were waterbaseline-corrected and measured at 20 °C.

IV. NUMERICAL RESULTS In this section we will present the results obtained by applying the PCM and QM/FQ/PCM to the calculation of the ECD spectra of (S)-nicotine. We will first show the results obtained by exploiting the PCM, which will be used both as a benchmark and as a reference for the assignment of the excited states and the assessment of their nature. We will then discuss the QM/ FQ/PCM results. In the following we will focus only on the neutral form of nicotine, which is the main species in aqueous solutions at pH greater than 9.62,63 A. PCM Results. The conformational distribution of nicotine has previously been studied both experimentally64 and theoretically.38,39,65−68 As a result, it has been shown that two possible conformations exist for neutral nicotine, denoted as A and B, which are shown in Figure 1. We calculated the populations of the two conformers at the CAM-B3LYP/aug-ccpVDZ/PCM(UFF) level of theory (for both the geometry optimization and free energy calculation) with an estimation of non-electrostatic contributions calculated using the SMD molecular cavity,53 and we obtained a 55.3% population of the A conformation. The possible electronic excitations that can be observed in nicotine originate from the n and π orbitals of the pyridine ring and the nonbonding orbital located on the nitrogen atom of the pyrrolidine ring (hereafter denoted as n′) and end in π* orbitals localized on the pyridine ring or in Rydberg states (Ryd), which are expected to be higher in energy. It is possible to classify the occupied and virtual π orbitals that lie on the pyridine moiety, in analogy to what is done for the Jahn− Teller-distorted benzene radical anion, as either quinoid (Q) or antiquinoid (AQ).69 The possible excitations in nicotine, therefore, include n−π* and π−π* transitions entirely localized on the pyridine ring (where both the occupied and virtual π orbitals can be either D

DOI: 10.1021/jp510542x J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A

indicated by our calculations. In fact, as remarked above, the La transition has a strongly positive rotatory strength for both conformations, and the next transitions with a negative sign appear at much higher energies. Clayton et al.63 measured both the UV and ECD spectra of aqueous nicotine at various pH, and their data indicated that this negative band disappears upon protonation of the pyridine nitrogen, which further supports our assignment. B. QM/FQ/PCM Results. The PCM was shown to yield results that qualitatively agree with the experimental values. Because of its two basic lone pairs, however, nicotine is able to form hydrogen bonds with water molecules. The effects due to specific interactions have been evaluated by exploiting an MD simulation followed by QM/FQ/PCM calculations. 1. Conformational Analysis. Figure 4 shows the value of δ1 sampled along the trajectory. The initial part of the simulation Figure 2. Calculated CAM-B3LYP/aug-cc-pVDZ/PCM and experimental UV spectra of (S)-nicotine in aqueous solution.

correlated wave function methods. At higher energies another strongly allowed transition appears, though the band maximum lies outside the experimental window. The ECD spectra are shown in Figure 3. The spectra of the A and B conformations show analogous sign patterns (also see

Figure 3. CAM-B3LYP/aug-cc-pVDZ/PCM-calculated and experimental ECD spectra of (S)-nicotine in aqueous solution.

Figure 4. (top) Value of the dihedral angle δ1 (see Figure 1) as a function of time. (bottom) Distribution (histogram) of δ1 values yielded by the MD trajectory.

Table I), but they differ much more than the absorption spectra, with the B conformation presenting a strong positive band at 225 nm assigned to the CTQ/La transition. The rotatory strength of the corresponding transition for the A conformation is much lower. The comparison with the experimental findings is fairly good: the sign pattern is reproduced, though there are notable differences in the relative intensities. Both the lowest-energy negative band at 240 nm and the strongly positive band at 225 nm have much higher relative intensities than their experimental counterparts. One possible reason for this behavior is that the latter band is assigned to the La transition, whose energy relative to the other bands is overestimated by DFT; if it were to appear at lower energies, it would partially cancel the negative band, which would also be quenched as a result, giving much better agreement with experiment. It should also be noted that previous experimental works62,63,70,71 have assigned the negative ECD band at 260 nm to the La transition and not the n−π*Q transition as

roughly represents an A conformer (the average angles were 106.4° and −58.1° up to 1500 ps and from 1500 to 10 000 ps, respectively). After a short time (about 1400 ps) a rotation about the δ1 angle places both nitrogen atoms on the same side of the molecule, which now most roughly corresponds to the B conformer. A third peak at about 0° is observed, corresponding to a partial rotation about δ1 that brings the two rings into the same plane. The average lifetime of such a conformation is very short, about 10.0 ps, while the system spends on average 100.0 ps in a “B-like” conformation with an average value of δ1 of −65°. B is therefore the most stable conformation as obtained by the simulation, at variance with QM/PCM findings. Because in the first 10 ns there was only one transition from one conformation to the other, we then extended the simulation to verify that the starting structure (nicotine A) was populated. At longer times we observed additional conformational transitions, though the overall dihedral distribution was not significantly affected (see the bottom panel of Figure 4). E

DOI: 10.1021/jp510542x J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A 2. Hydrogen-Bonding Patterns. Hydration patterns sampled by the MD simulation were analyzed by means of the nitrogen−hydrogen radial distribution function (g(r)) between atoms N1 and N2 in nicotine and H atoms in solvent molecules (Figure 5). Both nitrogen atoms feature a well-

Figure 6. Volumetric maps of the distribution of water hydrogen (green) and oxygen (red) atoms around the nitrogen atoms of nicotine for (top) the A conformer and (bottom) the B conformer. Figure 5. N−H radial distribution functions and corresponding running integrals for atoms N1 (black solid and dot-dashed lines, respectively) and N2 (red dashed and dot-dashed lines, respectively).

the presence of the high number of short-lived δ1 ≈ 0° species, since this is the conformation that minimizes the N1−N2 distance. To investigate this possibility, we analyzed water− water hydrogen bonding in the second hydration shell for these structures (a total of 3080 frames). The analysis was performed as follows: since both nitrogen atoms always have one neighbor solvent molecule, we searched for additional molecules that formed a hydrogen bond (using the same criteria as above) with the water molecule coordinated by N1 and were within a cutoff of 3.0 Å from N2. It turned out that 68% of the frames had two bridging water molecules between N1 and N2, while in 32% of the frames just one molecule was present (Figure 7); completing the rotation about δ1 would completely break down this network, which fits well with the high number of shortlived structures observed in the peak at δ1 ≈ 0°. 3. UV and ECD Spectra. We used the snapshots from the MD simulation to compute the UV absorption and ECD spectra of nicotine in the QM/FQ/PCM scheme. Each snapshot was centered on the center of mass of the solute, around which we cut a solvation sphere of radius 14 Å containing an average of 366 water molecules. A 2.5 Å larger radius was used for the PCM spherical cavity (as measured from the water oxygen atoms). A total of 12 excitations were computed for each snapshot. The computed and experimental spectra are plotted in Figure 8. The computed spectrum presents a strong band around 230 nm that originates mainly from the La transition, and there is a second strong absorption band at higher energies with contribution from the CT excitations. As was observed in the case of the PCM calculations, the computed spectrum is blueshifted with respect to the experimental one. Another difference is the separation between the two peaks, which is much smaller in the case of the QM/PCM calculations. The position of the La band is the same, and therefore, this behavior is to be ascribed to the blue shift of the second band, whose maximum moves from 200 to 180 nm, even though there is a shoulder at 200 nm. Because the structure of the molecule is flexible in the MD simulation, the nature of the states is affected, and therefore, this difference may be due to a conformational effect in addition to the difference in the solvation model. The ECD spectra are reported in Figure 9. The sign pattern of the

defined first peak at 1.60 Å, followed by depletion zones that extend up to 2.5 Å and then by a second peak at 2.91 Å. This indicates the presence of first and second hydration shells separated from the bulk solvent. The g(r) for the aliphatic nitrogen (N1) vanishes in the depletion zone (between 1.90 and 2.5 Å) and the coordination number is constant up to 2.75 Å, indicating that this atom probably forms a very stable hydrogen bond with one water molecule. On the other hand, the hydration structure of the aromatic N2 atom is much less defined compared with that of the aliphatic one: the g(r) does not drop to zero (i.e., there is always a non-negligible probability of finding a water molecule between the first and second hydration shells), and correspondingly, the running integrals N(r) grow linearly (if slowly) from 1.0 to 2.0 before the second peak (Figure 5). On the basis of this analysis of the radial distribution functions, a more detailed analysis of hydrogen-bonding patterns can be carried out: a water molecule was considered to form a hydrogen bond in a given snapshot if the hydrogen− acceptor distance was below 2.75 Å and the acceptor−donor− hydrogen angle was below 30.0°. Averages of 0.99 ± 0.08 and 1.14 ± 0.38 hydrogen bonds were obtained for atoms N1 and N2, respectively, confirming the higher rigidity of the hydration structure of the former. Figure 6 shows the spatial distribution functions of water H and O atoms near N1; frames were divided into A, B, and “intermediate” (δ1 between −30.0° and 30.0°) conformers on the basis of the analysis of Figure 4. N1 shows well-defined first and second hydration shells (note the outermost oxygen probability density maxima), while the N2 hydration structure is less defined. However, when conformer B is considered, it can be observed that the first oxygen isosurfaces for N2 are larger on the side of the pyridine ring nearer to N1; this may indicate the existence of a cooperative effect in which bridging water−water hydrogen bonds stabilize the B conformer, creating less disorder in the water liquid structure. Such an effect may not exist for A conformers but could partially explain F

DOI: 10.1021/jp510542x J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A

Figure 7. Hydrogen-bonding patterns for δ1 ≈ 0° conformers. N−H and O−H distances are indicated. (left) Snapshot with two bridging water molecules. (right) Snapshot with one bridging water molecule.

ECD spectra, respectively, correspond to the same transition (La) and that the blue shift in the band maximum in the ECD spectrum is caused by the partial overlap with the negative band that appears at 235 nm. The observed differences between the QM/PCM and QM/FQ/PCM results can be attributed both to the introduction of specific solvation effects in the TDDFT calculations and to a different conformational sampling from the MD simulations.

V. SUMMARY AND CONCLUSIONS In this paper, we have presented an extension of the QM/FQ/ PCM approach to ECD spectra and its application to the description of the UV and ECD spectra of neutral (S)-nicotine in aqueous solution. The results of the modeling have been compared with those of the PCM approach and to newly measured experimental data. The results obtained show that while both methods correctly reproduce the shape of the UV spectra and the sign pattern of the ECD spectra, the QM/FQ/ PCM approach is a great improvement over the PCM picture in terms of relative intensities, especially in the case of ECD, and an inprovement in the relative positions of the UV and ECD bands is also observed. The differences in the results obtained with the two methods originate from both the more detailed description of the solute−solvent interactions, which include directional effects such as hydrogen bonding, offered by the QM/FQ/PCM approach and the different conformational sampling obtained from the molecular dynamics simulations. The latter could be improved even further by changing the parametrization of the force field and increasing the sampling time, though doing the latter would increase the cost of the simulation and the number of snapshots to be considered in the subsequent QM calculations. Both the parametrization and simulation time strongly affect the accuracy and reliability of mixed explicit/implicit approaches, whereas purely implicit approaches are based on only a very limited number of parameters (the solvent dielectric constants and cavity shape). On the basis of such considerations, mixed explicit/implicit models are physically more adequate to treat flexible systems, for which a QM treatment will require the evaluation of a large number of minima, as well as solute−solvent couples dominated by specific interactions, which are only partially considered in implicit strategies. However, continuum modes are still able to provide qualitatively correct results at a much lower cost, in addition to providing a simpler way to characterize the electronic excitations of the system and thus ease the spectral assignment.

Figure 8. CAM-B3LYP/aug-cc-pVDZ/FQ/PCM-calculated and experimental UV spectra of (S)-nicotine in water. The spectra were obtained by averaging spectra from single snapshots.

Figure 9. CAM-B3LYP/aug-cc-pVDZ/FQ/PCM-calculated and experimental ECD spectra of (S)-nicotine in aqueous solution. The spectra were obtained by averaging spectra from single snapshots.

computed spectrum agrees with the experimental one, even though the relative intensity of the positive band with respect to the negative bands is higher in comparison with the experimental data. This difference, however, is not as large as in the QM/PCM results of Figure 3, demonstrating a significant improvement in the results. It should be noted that the peaks at about 230 and 220 nm in the absorption and G

DOI: 10.1021/jp510542x J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A

(15) Tam, M. C.; Russ, N. J.; Crawford, T. D. Coupled Cluster Calculations of Optical Rotatory Dispersion of (S)-Methyloxirane. J. Chem. Phys. 2004, 121, 3550−3557. (16) Lipparini, F.; Egidi, F.; Cappelli, C.; Barone, V. The Optical Rotation of Methyloxirane in Aqueous Solution: A Never Ending Story? J. Chem. Theory Comput. 2013, 9, 1880−1884. (17) Eliel, E. L.; Wilen, S. H.; Doyle, M. P. Basic Organic Stereochemistry; Wiley Interscience: New York, 2001. (18) Macleod, N. A.; Butz, P.; Simons, J. P.; Grant, G. H.; Baker, C. M.; Tranter, G. E. Structure, Electronic Circular Dichroism and Raman Optical Activity in the Gas Phase and in Solution: A Computational and Experimental Investigation. Phys. Chem. Chem. Phys. 2005, 7, 1432−1440. (19) Gregoliński, J.; Ślepokura, K.; Lisowski, J. Lanthanide Complexes of the Chiral Hexaaza Macrocycle and Its meso-Type Isomer: Solvent-Controlled Helicity Inversion. Inorg. Chem. 2007, 46, 7923−7934. (20) Cappelli, C.; Mennucci, B. Modeling the Solvation of Peptides. The Case of (S)-N-Acetylproline Amide in Liquid Water. J. Phys. Chem. B 2008, 112, 3441−3450. (21) Autschbach, J.; Nitsch-Velasquez, L.; Rudolph, M. TimeDependent Density Functional Response Theory for Electronic Chiroptical Properties of Chiral Molecules. Top. Curr. Chem. 2011, 298, 1−98. (22) Lahiri, P.; Wiberg, K. B.; Vaccaro, P. H.; Caricato, M.; Crawford, T. D. Large Solvation Effect in the Optical Rotatory Dispersion of Norbornenone. Angew. Chem., Int. Ed. 2014, 53, 1386−1389. (23) Mukhopadhyay, P.; Zuber, G.; Wipf, P.; Beratan, D. N. Contribution of a Solute’s Chiral Solvent Imprint to Optical Rotation. Angew. Chem., Int. Ed. 2007, 46, 6450−6452. (24) Cheeseman, J. R.; Shaik, M. S.; Popelier, P. L. A.; Blanch, E. W. Calculation of Raman Optical Activity Spectra of Methyl-D-Glucose Incorporating a Full Molecular Dynamics Simulation of Hydration Effects. J. Am. Chem. Soc. 2011, 133, 4991−4997. (25) Hopmann, K. H.; Ruud, K.; Pecul, M.; Kudelski, A.; Dracinsky, M.; Bour, P. Explicit versus Implicit Solvent Modeling of Raman Optical Activity Spectra. J. Phys. Chem. B 2011, 115, 4128−4137. (26) Jurinovich, S.; Pescitelli, G.; Di Bari, L.; Mennucci, B. A TDDFT/MMPol/PCM Model for the Simulation of Exciton-Coupled Circular Dichroism Spectra. Phys. Chem. Chem. Phys. 2014, 16, 16407−16418. (27) Cappelli, C.; Corni, S.; Mennucci, B.; Cammi, R.; Tomasi, J. Vibrational Circular Dichroism within the Polarizable Continuum Model: A Theoretical Evidence of Conformation Effects and Hydrogen Bonding for (S)-(−)-3-Butyn-2-ol in CCl4 Solution. J. Phys. Chem. A 2002, 106, 12331−12339. (28) Kongsted, J.; Ruud, K. Solvent Effects on Zero-Point Vibrational Corrections to Optical Rotations and Nuclear Magnetic Resonance Shielding Constants. Chem. Phys. Lett. 2008, 451, 226−232. (29) Wilson, S. M.; Wiberg, K. B.; Murphy, M. J.; Vaccaro, P. H. The Effects of Conformation and Solvation on Optical Rotation: Substituted Epoxides. Chirality 2008, 20, 357−369. (30) Cappelli, C. In Continuum Solvation Models in Chemical Physics; Cammi, R., Mennucci, B., Eds.; John Wiley & Sons: Chichester, U.K., 2012; pp 167−179. (31) Cappelli, C.; Lipparini, F.; Bloino, J.; Barone, V. Towards an Accurate Description of Anharmonic Infrared Spectra in Solution within the Polarizable Continuum Model: Reaction Field, Cavity Field and Nonequilibrium Effects. J. Chem. Phys. 2011, 135, No. 104505. (32) Cappelli, C.; Mennucci, B.; Monti, S. Environmental Effects on the Spectroscopic Properties of Gallic Acid: A Combined Classical and Quantum Mechanical Study. J. Phys. Chem. A 2005, 109, 1933−1943. (33) Neugebauer, J. Induced Chirality in Achiral Media: How Theory Unravels Mysterious Solvent Effects. Angew. Chem., Int. Ed. 2007, 46, 7738−7740. (34) Lipparini, F.; Barone, V. Polarizable Force Fields and Polarizable Continuum Model: A Fluctuating Charges/PCM Approach. 1. Theory and Implementation. J. Chem. Theory Comput. 2011, 7, 3711−3724.

More extended testing of the performance of the QM/FQ/ PCM method is surely needed, but from the experience gained to date, QM/PCM and QM/FQ/PCM are to be seen as complementary approaches for the investigation of the spectroscopic properties of systems in aqueous solution.



AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. Notes

The authors declare no competing financial interest.



ACKNOWLEDGMENTS C.C. acknowledges support from the Italian MIUR PRIN 2012 NB3KLK_002 and FIRB 2010-Futuro in Ricerca: protocollo: RBFR10Y5VW. The high-performance computing facilities of the DREAMS center (http://dreamslab.sns.it) are acknowledged for providing computer resources. The support of COST CMST-Action CM1002 COnvergent Distributed Environment for Computational Spectroscopy (CODECS) is also acknowledged. I.C. acknowledges Compunet of IIT for the financial support.



REFERENCES

(1) Crawford, T. Ab Initio Calculation of Molecular Chiroptical Properties. Theor. Chem. Acc. 2006, 115, 227−245. (2) Autschbach, J. In Comprehensive Chiroptical Spectroscopy, Vol. 1; Berova, N., Polavarapu, P., Nakanishi, K., Woody, R. W., Eds.; Wiley: New York, 2012; pp 593−642. (3) Crawford, T. D. In Comprehensive Chiroptical Spectroscopy, Vol. 1; Berova, N., Polavarapu, P., Nakanishi, K., Woody, R. W., Eds.; Wiley: New York, 2012; pp 675−697. (4) Goerigk, L.; Kruse, H.; Grimme, S. In Comprehensive Chiroptical Spectroscopy, Vol. 1; Berova, N., Polavarapu, P., Nakanishi, K., Woody, R. W., Eds.; Wiley: New York, 2012; pp 643−674. (5) Superchi, S.; Rosini, C.; Mazzeo, G.; Giorgio, E. In Comprehensive Chiroptical Spectroscopy, Vol. 2; Berova, N., Polavarapu, P., Nakanishi, K., Woody, R. W., Eds.; Wiley: New York, 2012; pp 421−448. (6) Comprehensive Chiroptical Spectroscopy; Berova, N., Polavarapu, P. L., Nakanishi, K., Woody, R. W., Eds.; Wiley: New York, 2012. (7) Pecul, M.; Ruud, K. In Continuum Solvation Models in Chemical Physics; Mennucci, B., Cammi, R., Eds.; Wiley: New York, 2007; Chapter 1.4, pp 64−81. (8) Mennucci, B.; Cappelli, C.; Cammi, R.; Tomasi, J. Modeling Solvent Effects on Chiroptical Properties. Chirality 2011, 23, 717− 729. (9) Egidi, F.; Barone, V.; Bloino, J.; Cappelli, C. Toward an Accurate Modeling of Optical Rotation for Solvated Systems: Anharmonic Vibrational Contributions Coupled to the Polarizable Continuum Model. J. Chem. Theory Comput. 2012, 8, 585−597. (10) Mukhopadhyay, P.; Zuber, G.; Goldsmith, M.-R.; Wipf, P.; Beratan, D. N. Solvent Effect on Optical Rotation: A Case Study of Methyloxirane in Water. ChemPhysChem 2006, 7, 2483−2486. (11) Wiberg, K. B.; Vaccaro, P. H.; Cheeseman, J. R. Conformational Effects on Optical Rotation. 3-Substituted 1-Butenes. J. Am. Chem. Soc. 2003, 125, 1888−1896. (12) Ruud, K.; Zanasi, R. The Importance of Molecular Vibrations: The Sign Change of the Optical Rotation of Methyloxirane. Angew. Chem., Int. Ed. 2005, 44, 3594−3596. (13) Mort, B. C.; Autschbach, J. A Pragmatic Recipe for the Treatment of Hindered Rotations in the Vibrational Averaging of Molecular Properties. ChemPhysChem 2008, 9, 159−170. (14) Kongsted, J.; Pedersen, T. B.; Jensen, L.; Hansen, A. E.; Mikkelsen, K. V. Coupled Cluster and Density Functional Theory Studies of the Vibrational Contribution to the Optical Rotation of (S)Propylene Oxide. J. Am. Chem. Soc. 2006, 128, 976−982. H

DOI: 10.1021/jp510542x J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A (35) Lipparini, F.; Cappelli, C.; Barone, V. Linear Response Theory and Electronic Transition Energies for a Fully Polarizable QM/ Classical Hamiltonian. J. Chem. Theory Comput. 2012, 8, 4153−4165. (36) Lipparini, F.; Cappelli, C.; Scalmani, G.; De Mitri, N.; Barone, V. Analytical First and Second Derivatives for a Fully Polarizable QM/ Classical Hamiltonian. J. Chem. Theory Comput. 2012, 8, 4270−4278. (37) Lipparini, F.; Cappelli, C.; Barone, V. A Gauge Invariant Multiscale Approach to Magnetic Spectroscopies in Condensed Phase: General Three-Layer Model, Computational Implementation and Pilot Applications. J. Chem. Phys. 2013, 138, No. 234108. (38) Egidi, F.; Bloino, J.; Cappelli, C.; Barone, V. Development of a Virtual Spectrometer for Chiroptical Spectroscopies: The Case of Nicotine. Chirality 2013, 25, 701−708. (39) Egidi, F.; Bloino, J.; Cappelli, C.; Barone, V.; Tomasi, J. Tuning of NMR and EPR Parameters by Vibrational Averaging and Environmental Effects: An Integrated Computational Approach. Mol. Phys. 2013, 111, 1345−1354. (40) Mortier, W. J.; Van Genechten, K.; Gasteiger, J. Electronegativity Equalization: Application and Parametrization. J. Am. Chem. Soc. 1985, 107, 829−835. (41) Rick, S. W.; Stuart, S. J.; Berne, B. J. Dynamical Fluctuating Charge Force Fields: Application to Liquid Water. J. Chem. Phys. 1994, 101, 6141−6156. (42) Rappé, A. K.; Goddard, W. A., III Charge Equilibration for Molecular Dynamics Simulations. J. Phys. Chem. 1991, 95, 3358−3363. (43) Caricato, M.; Lipparini, F.; Scalmani, G.; Cappelli, C.; Barone, V. Vertical Electronic Excitations in Solution with the EOM-CCSD Method Combined with a Polarizable Explicit/Implicit Solvent Model. J. Chem. Theory Comput. 2013, 9, 3035−3042. (44) Tomasi, J.; Mennucci, B.; Cammi, R. Quantum Mechanical Continuum Solvation Models. Chem. Rev. 2005, 105, 2999−3093. (45) Mennucci, B. Continuum Solvation Models: What Else Can We Learn from Them? J. Phys. Chem. Lett. 2010, 1, 1666−1674. (46) Casida, M. In Recent Advances in Density Functional Methods (Part I); Chong, D., Ed.; World Scientific: Singapore, 1997; pp 155− 192. (47) Pedersen, T. B.; Hansen, A. E. Ab Initio Calculation and Display of the Rotatory Strength Tensor in the Random Phase Approximation. Method and Model Studies. Chem. Phys. Lett. 1995, 246, 1−8. (48) Frisch, M. J.; Trucks, G. W.; Schlegel, H. B.; Scuseria, G. E.; Robb, M. A.; Cheeseman, J. R.; Scalmani, G.; Barone, V.; Mennucci, B.; Petersson, G. A.; et al. Gaussian 09, revision D.01; Gaussian, Inc.: Wallingford, CT, 2009. (49) Yanai, T.; Tew, D. P.; Handy, N. C. A New Hybrid Exchange− Correlation Functional Using the Coulomb-Attenuating Method (CAM-B3LYP). Chem. Phys. Lett. 2004, 393, 51−57. (50) Egidi, F.; Segado, M.; Koch, H.; Cappelli, C.; Barone, V. A Benchmark Study of Electronic Excitation Energies, Transition Moments, and Excited-State Energy Gradients on the Nicotine Molecule. J. Chem. Phys. 2014, 141, No. 224114. (51) Kendall, R. A.; Dunning, T. H., Jr.; Harrison, R. J. Electron Affinities of the First-Row Atoms Revisited. Systematic Basis Sets and Wave Functions. J. Chem. Phys. 1992, 96, 6796−6806. (52) Mennucci, B. Polarizable Continuum Model. Wiley Interdiscip. Rev.: Comput. Mol. Sci. 2012, 2, 386−404. (53) Marenich, A. V.; Cramer, C. J.; Truhlar, D. G. Universal Solvation Model Based on Solute Electron Density and a Continuum Model of the Solvent Defined by the Bulk Dielectric Constant and Atomic Surface Tensions. J. Phys. Chem. B 2009, 113, 6378−6396. (54) Barone, V.; Cimino, P.; Stendardo, E. Development and Validation of the B3LYP/N07D Computational Model for Structural Parameter and Magnetic Tensors of Large Free Radicals. J. Chem. Theory Comput. 2008, 4, 751−764. (55) Cieplak, P.; Cornell, W. D.; Bayly, C.; Kollman, P. A. Application of the Multimolecule and Multiconformational RESP Methodology to Biopolymers: Charge Derivation for DNA, RNA, and Proteins. J. Comput. Chem. 1995, 16, 1357−1377.

(56) Wang, J.; Wolf, R. M.; Caldwell, J. W.; Kollman, P. A.; Case, D. A. Development and testing of a general amber force field. J. Comput. Chem. 2004, 25, 1157−1174. (57) Pronk, S.; Pall, S.; Schulz, R.; Larsson, P.; Bjelkmar, P.; Apostolov, R.; Shirts, M. R.; Smith, J. C.; Kasson, P. M.; van der Spoel, D.; et al. GROMACS 4.5: A High-Throughput and Highly Parallel Open Source Molecular Simulation Toolkit. Bioinformatics 2013, 29, 845−854. (58) Abascal, J. L. F.; Vega, C. A General Purpose Model for the Condensed Phases of Water: TIP4P/2005. J. Chem. Phys. 2005, 123, No. 234505. (59) Darden, T.; York, D.; Pedersen, L. Particle mesh Ewald: An N log(N) Method for Ewald Sums in Large Systems. J. Chem. Phys. 1993, 98, 10089−10092. (60) Bussi, G.; Donadio, D.; Parrinello, M. Canonical Sampling through Velocity Rescaling. J. Chem. Phys. 2007, 126, No. 014101. (61) Humphrey, W.; Dalke, A.; Schulten, K. VMD: Visual Molecular Dynamics. J. Mol. Graphics 1996, 14, 33−38. (62) Clayton, P. M.; Vas, C. A.; Bui, T. T. T.; Drake, A. F.; McAdam, K. Spectroscopic Studies on Nicotine and Nornicotine in the UV Region. Chirality 2013, 5, 288−293. (63) Clayton, P. M.; Vas, C. A.; Bui, T. T. T.; Drake, A. F.; McAdam, K. Spectroscopic Investigations into the Acid−Base Properties of Nicotine at Different Temperatures. Anal. Methods 2013, 5, 81−88. (64) Takeshima, T.; Fukumoto, R.; Egawa, T.; Konaka, S. Molecular Structure of Nicotine As Studied by Gas Electron Diffraction Combined with Theoretical Calculations. J. Phys. Chem. A 2002, 106, 8734−8740. (65) Elmore, D. E.; Dougherty, D. A. A Computational Study of Nicotine Conformations in the Gas Phase and in Water. J. Org. Chem. 2000, 65, 742−747. (66) Mora, M.; Castro, M. E.; Nino, A.; Melendez, F. J.; MuñozCaro, C. Analysis of B3LYP and MP2 Conformational Population Distributions in trans-Nicotine, Acetylcholine, and ABT-594. Int. J. Quantum Chem. 2005, 103, 25−33. (67) Muñoz-Caro, C.; Nino, A.; Mora, M.; Reyes, S.; Melendez, F. J.; Castro, M. E. Conformational Population Distribution of Acetylcholine, Nicotine and Muscarine in Vacuum and Solution. J. Mol. Struct.: THEOCHEM 2005, 726, 115−124. (68) Ortega, P. G. R.; Montejo, M.; González, J. J. L. Vibrational Circular Dichroism and Theoretical Study of the Conformational Equilibrium in (−)-S-Nicotine. ChemPhysChem 2014, DOI: 10.1002/ cphc.201402652. (69) Cogan, S.; Zilberg, S.; Haas, Y. The Electronic Origin of the Dual Fluorescence in Donor−Acceptor Substituted Benzene Derivatives. J. Am. Chem. Soc. 2006, 128, 3335−3345. (70) Tomita, A.; Ochiai, E.-I.; Hirai, H.; Makishi, S. The Effect of Hydrogen Chloride on the Optical Rotation of Nicotine. J. Org. Chem. 1966, 31, 4307−4308. (71) Testa, B.; Jenner, P. Circular Dichroic Determination of the Preferred Conformation of Nicotine and Related Chiral Alkaloids in Aqueous Solution. Mol. Pharmacol. 1973, 9, 10−16.

I

DOI: 10.1021/jp510542x J. Phys. Chem. A XXXX, XXX, XXX−XXX