pubs.acs.org/JPCL
Thermodynamic Insights into the Dynamic Switching of a Cyclodextrin in a Bistable Molecular Shuttle Peng Liu,† Wensheng Cai,*,† Christophe Chipot,‡,§ and Xueguang Shao† †
College of Chemistry, Nankai University, Tianjin, 300071, People's Republic of China, and ‡Theoretical and Computational Biophysics Group, Beckman Institute, University of Illinois at Urbana-Champaign, Urbana, Illinois 61801
ABSTRACT A prototypical molecular shuttle formed by dodecamethylene and 4,40 -bipyridinium units and an R-cyclodextrin (CD) was investigated employing molecular dynamics simulations. The free-energy profile characterizing the shuttling process of the R-CD along the molecular thread was determined using the adaptive biasing force method, revealing two thermodynamically stable states separated by a pronounced energy barrier. The free-energy barrier with respect to the stable states is calculated to be þ20.6 kcal/mol, in excellent agreement with the experimentally measured quantity. Partitioning of the free energy into contributions of different nature indicates that the shuttling process is primarily controlled by the favorable inclusion of the dodecamethylene chain and the unfavorable inclusion of the cationic bipyridinium group by R-CD. The predominant contribution to the energy barrier stems from the disruption of the solvation shell of the charged group. Deciphering the molecular mechanism of the shuttling process is expected to help design new controllable molecular shuttles. SECTION Statistical Mechanics, Thermodynamics, Medium Effects
D
An efficient approach, coined the adaptive biasing force (ABF),19-23 was developed to calculate the free energy along a chosen order parameter from unconstrained MD simulations. This scheme has been used previously to explore the recognition and association processes of CDs with steroid drugs.24,25 In ABF, a biasing force is rapidly estimated and refined to erase the ruggedness of the free-energy surface and, hence, allows the order parameter to be sampled uniformly. It can, therefore, effectively accelerate barrier crossing and improve the accuracy of the calculated free energies. In the present contribution, this method combined with unconstrained MD simulations was employed to investigate the shuttling process of a CD in a rotaxane designed by the Harada research group.11 This molecular shuttle is formed by an R-CD and a dumbbell-shaped thread composed of two dodecamethylene units (stations), three 4,40 -bipyridinium units (linkers), and two 2,4-dinitrophenyl groups (stoppers), as shown in Scheme 1. The shuttling ability of R-CD was found to be solvent- and temperature-sensitive and was confirmed by 1H NMR spectra to occur in dimethylsulfoxide (DMSO) at 130 C. At coalescence temperature, the rate constant characterizing the shuttling process is determined from the difference in the chemical shifts for the relevant target split signals. The free energy of activation, ΔG‡, is subsequently inferred based on Eyring's law to be on the order of þ20 kcal/mol.11,26 The primary thrust of the present study is the demonstration that
esign, synthesis, and characterization of molecular devices have recently emerged as promising research areas on account of their potential applications in miniaturization of the components of electronic devices. Of particular interest, rotaxanes, mechanically interlocked molecular architectures, have been widely employed as molecular frameworks to construct such supramolecular assemblies. Molecular shuttles provide a cogent illustration of the characteristic features shared by these interlocked assemblies.1,2 A cyclic molecule, like cyclodextrin (CD), threaded on a dumbbell-shaped strand, shuttles between two or more thermodynamically stable states separated by linkers. This translocation can be triggered by external stimuli, such as a redox process,3,4 a pH change,5,6 light,7-9 metal-ion coordination,10 as well as the temperature and the solvent.9,11 The shuttling process forms the basic function of molecular switches. Consequently, considerable effort on the experimental front has been invested for their study.11-13 On the computational front, however, only a handful of investigations have tackled this problem hitherto, limiting their scope to the energy of discrete states.14-17 Exploration of the freeenergy landscape that characterizes the shuttling process has been seldom reported. The only investigation undertaken for electrochemically switchable rotaxanes is that recently proposed by Goddard et al.18 using constrained molecular dynamics (MD) simulations. This scarcity of theoretical studies can probably be ascribed to the existence of pronounced energy barriers, which are anticipated to be overcome only with great difficulty on the nanosecond time scale commonly amenable to MD.
r 2010 American Chemical Society
Received Date: March 30, 2010 Accepted Date: May 18, 2010 Published on Web Date: May 24, 2010
1776
DOI: 10.1021/jz100414j |J. Phys. Chem. Lett. 2010, 1, 1776–1780
pubs.acs.org/JPCL
Scheme 1. [2]Rotaxane Molecule Formed by an R-CD, Two Dodecamethylene Units (Stations), Three 4,40 -Bipyridinium Units (Linkers), and Two 2,4-Dinitrophenyl Groups (Stoppers)
shutting processes can be tackled at the theoretical level, employing the appropriate methodology. In addition, new light is being shed on the atomic-level detail of the phenomenon and its underlying molecular mechanism. Improved knowledge of molecular shuttles is envisioned to help design novel rotaxane structures with a target rate of site exchange. All of the MD simulations reported here rely upon the CHARMM27 force filed.27 A number of parameters describing the linkers and the stoppers of the chain are absent from this potential energy function. The missing parameters were optimized using the strategies suggested by MacKerell et al.28,29 Details of the parametrization are supplied in the Supporting Information. The carbohydrate solution force field (CSFF)30 parameters for CD were employed. A box of DMSO molecules of suitable size was constructed and equilibrated employing the all-atom DMSO model proposed by Strader and Feller.31 This additive model offers a reasonable reproduction of the experimental dielectric constant of the liquid. The initial supramolecular assemblies featuring the chain and an R-CD were immersed in this solvent bath. All MD simulations were carried out with the NAMD 2.6 program.32 The reaction coordinate, ξ, was chosen as the distance separating the center of mass of the CD from that of the central linker. In order to distinguish the position of the former with respect to that of the latter, ξ was chosen arbitrarily to be negative or positive depending on whether the CD was located on the left- or on the right-hand side of the central linker, respectively (see Scheme 1). The pathway, which extends approximately from -15 to þ15 Å, was divided into 30 nonoverlapping windows. Instantaneous values of the force were accrued in bins 0.1 Å wide. The variation of the free energy, ΔG(ξ), was determined by integrating the average force acting on ξ.19,21 To ensure that the free energy associated with the shuttling process is converged over the time scale explored, the backbone of the dumbbell thread was restrained by means of a harmonic potential, excluding the middle linker, which was left free to rotate. The entire shuttling process was explored in two independent simulations, wherein the CD diffuses freely from the left station to the central linker and from the right station to this point. The orientation of the CD is illustrated in Scheme 1. The simulation time was extended incrementally to probe the convergence of the free energy in regions featuring barriers and minima, corresponding to an overall simulation time of 240 ns. Block-average regression was utilized to estimate the standard error associated with the free-energy differences.20,21,33 The potential of mean force (PMF) characterizing the motion of R-CD along ξ is depicted in Figure 1B. As can be
r 2010 American Chemical Society
Figure 1. (A) Representative configurations of R-CD shuttled along the polymeric chain of the rotaxane. (B) Free-energy profile delineating the shuttling process along ξ. The error bars represent the standard error of the free-energy difference at each ξ. Insets: Number of accrued force samples and average force acting. Image rendered with VMD.34
seen, the curve is essentially symmetrical with respect to ξ = 0 Å. Two relatively even regions, namely, -15 e ξ e -7 and þ7 e ξ e þ15 Å, are separated by a broad and high energy barrier spanning ∼-7 < ξ < þ7 Å. These two even regions correspond to the stable states wherein the R-CD is found on the dodecamethylene section of the thread. The barrier reflects the metastable state of the supramolecular assembly, wherein the R-CD coincides with the central bipyridinium moiety. Three representative configurations are gathered in Figure 1A. The slight asymmetry of the PMF can be related to the structural asymmetry of the R-CD. A closer look at Figure 1B reveals that the free energy between -12 and -10 Å remains roughly constant at a value of ∼-1.9 kcal/mol. It abruptly increases as R-CD approaches the central linker, leading to the first maximum at around -1.5 Å, with a free energy of ∼17.9 kcal/mol. In the structure corresponding to this maximum, the pyridine ring A is found to coincide with the middle plane of the R-CD. Sequentially, a local minimum emerges at 0 Å, where the center of mass of the R-CD is located at the midpoint of the bond connecting the A and A0 rings. Symmetrically, a second maximum arises at around þ1.5 Å with a free energy of ∼16.9 kcal/mol, corresponding
1777
DOI: 10.1021/jz100414j |J. Phys. Chem. Lett. 2010, 1, 1776–1780
pubs.acs.org/JPCL
Figure 2. (A) Decomposition of the total free energy into van der Waals CD-thread, electrostatic CD-thread, and CD-solvent contributions. (B) Fluctuation of the area of the primary rim of the R-CD as a function of ξ. Partitioning of the CD-thread interactions into (C) CD OHs-thread electrostatic and van der Waals contributions, and (D) CD cavity-thread electrostatic and van der Waals contributions.
to a similar structure, wherein the pyridine ring A0 is docked near the center of the cavity. As the R-CD moves forth, the free energy sharply decreases and reaches a plateau at ∼-3.5 kcal/mol found between þ11 and þ13 Å. The freeenergy difference between the central linker and the station for the forward and the backward transfer processes is, therefore, estimated to be þ19.8 ( 0.1 and þ21.4 ( 0.2 kcal/mol for the left- and the right-hand side pathway, respectively, the difference of which falls within chemical accuracy. Its mean value, equal to þ20.6 ( 0.2 kcal/mol, agrees nicely with the free energy of activation for the siteexchange process determined experimentally.11 On the basis of this measure of ΔG‡, the reaction time required to overcome the free-energy barrier is estimated, using the Eyring equation,26 to be approximately 1.6 10-2 s at 130 C and 98.6 s at 30 C. This result rationalizes the experimental observation that the R-CD shuttles quickly between the two stations at higher temperature but that site exchange does not occur on the NMR time scale at room temperature. To understand how the forces at play vary in the course of the shuttle process, a number of intermolecular interactions were monitored as a function of time. This was achieved by (i) partitioning the total average force acting along ξ into two components, namely, CD-thread and CD-solvent, and (ii) integrating these forces separately. The resulting free-energy contributions are depicted in Figure 2A. As can be observed, the shape of the CD-chain van der Waals term is akin to that of the complete PMF. The two stable states stem from favorable van der Waals interactions of the CD with the dodecamethylene chain. The free-energy barrier for this component arises at the position of the central linker and can be ascribed to steric hindrances. To delve further into this observation, the geometric change of the R-CD cavity was examined by measuring the average area of the hexagon formed by the six C6 atoms in the primary rim of
r 2010 American Chemical Society
the CD, as shown in Figure 2B. Expansion of the cavity can be visualized in the pronounced peak found at around 0 Å, which corresponds to the R-CD passing through the central linker. Yet, the effect of such steric hindrances does not apparently constitute the main contribution to the barrier seen in the PMF. At variance with the repulsive interactions, the free energy due to CD-thread electrostatic interactions abruptly decreases as ξ approaches zero, indicating a positive contribution to the shuttling process. Additional analysis was carried out by partitioning the CD-thread interactions into contributions from the hydroxyl groups and the atoms lining the cavity. The corresponding free-energy components are depicted in Figure 2C and D, suggesting that attractive electrostatic interactions arise primarily from the hydroxyl groups of the CD, whereas the electrostatic interaction of atoms lining the cavity with the positively charged linker is, as expected, repulsive. The latter observation is in line with several previous reports that the interior of the CD exerts repulsive interaction with positively charged moieties.35,36 In sharp contrast with steric hindrances, CD-solvent interactions are found to constitute the main contribution to the barrier of the PMF (see Figure 2A). To understand the effect of the solvent on the shuttling process, the distribution of the DMSO molecules around the central linker was monitored. This was achieved by computing the radial distribution functions, g(r;ξ), of the sulfur atom of DMSO with respect to the nitrogen atoms of the linker for each value of ξ. The distributions of DMSO around the A and A0 rings are shown in Figure 3A and B. The maximum of the density found in the region spanning þ4 e r e þ6 Å observed in both figures corresponds to the first solvation shell of the central linker. This solvation shell is, however, discontinuous along ξ, in the range of -8 e ξ e 0 Å in Figure 3A and 0 e ξ e þ8 Å in Figure 3B, which implies that solvation of the nitrogen atoms
1778
DOI: 10.1021/jz100414j |J. Phys. Chem. Lett. 2010, 1, 1776–1780
pubs.acs.org/JPCL
CD is extremely disfavored. This, in turn, modulates considerably the free-energy barrier due to CD-chain van der Waals interactions. To conclude, the free-energy profile delineating the motion of an R-CD along the polymeric chain of the rotaxane has been determined as a paradigm for the investigation of the shuttling process of a CD between two stable states. The free-energy barrier for switching is estimated to be þ20.6 ( 0.2 kcal/mol, which agrees well with the experimental value of about 20 kcal/mol.11 The interactions that underly the shuttling process were disentangled by partitioning the PMF into different components. The present results reveal that the sliding motion of the R-CD over the two dodecamethylene chains corresponds to two thermodynamically stable states separated by a pronounced barrier, which arises from steric hindrances and unfavorable CD-solvent interactions. The latter represent the predominant contribution to the free-energy barrier and stem from the disruption of the solvation shell of the charged 4,40 -bipyridinium group. The present work paves the way for further investigation aimed at understanding the influence of the solvent on the shuttling process.
SUPPORTING INFORMATION AVAILABLE Parameterization and simulation details. This material is available free of charge via the Internet at http://pubs.acs.org.
AUTHOR INFORMATION Corresponding Author: *To whom correspondence should be addressed. Tel: þ86-2223503430. Fax: þ86-22-23502458. E-mail:
[email protected].
Notes Figure 3. (A) Evolution along ξ of the radial distribution function of the N(central linker)-S(DMSO) pair; r is the distance of the atomic pair. The N atom pertains to the A ring of the 4,40 -bipyridinium moiety. (B) Same as (A) with the A0 ring of the central linker. (C) Volumetric map of the DMSO density around the rotaxane in its free-energy minimum, ξ ≈ þ11 Å. (D) Same as (C) at the transition state, ξ ≈ þ2 Å. (E) Orientational anisotropy of DMSO near the chain measured in terms of the average of cos θ as a function of |R|, where θ is the angle formed between the dipole moment of DMSO and vector R that connects the COM of DMSO and the central bipyridinium moiety. |R| = 3.5-5.5 Å corresponds to the first solvation shell.
On leave from Equipe de Dynamique des Assemblages Membranaires, UMR 7565, Nancy Universit e, BP 239, 54506 Vandoeuvre-l esNancy cedex, France. §
ACKNOWLEDGMENT The study was supported by the National
Natural Science Foundation of China (Nos. 20873066 and 20835002). C.C. acknowledges the French Embassy in Beijing for travel support.
pertaining to the A and A0 rings is disrupted as the CD moves toward the latter. Figure 3C and D depicts this evolution by means of three-dimensional density maps of DMSO around the rotaxane calculated for two representative states. In addition, the orientational anisotropy of the solvent near the central bipyridinium unit was measured, as shown in Figure 3E. The peak in the distribution of Æcos θ æ in a radial range characterizing the first solvation shell indicates that the electrostatic interaction of the positively charged moiety and its surrounding DMSO molecules is appreciable. Breaking the solvation shell is energetically unfavorable, resulting in a significant increase of CD-solvent interactions. Furthermore, this energetic increment cannot be entirely compensated for by favorable CD-chain electrostatic interactions, as can be seen in the peak of Figure 2A, suggestive that inclusion of cationic groups in the
r 2010 American Chemical Society
REFERENCES (1) (2) (3)
(4)
(5)
1779
Raymo, F. M.; Stoddart, J. F. Interlocked Macromolecules. Chem. Rev. 1999, 99, 1643–1663. Harada, A. Cyclodextrin-Based Molecular Machines. Acc. Chem. Res. 2001, 34, 456–464. Bissell, R. A.; C ordova, E.; Kaifer, A. E.; Stoddart, J. F. A Chemically and Electrochemically Switchable Molecular Shuttle. Nature 1994, 369, 133–137. Zhao, Y.-L.; Dichtel, W. R.; Trabolsi, A.; Saha, S.; Aprahamian, I.; Stoddart, J. F. A Redox-Switchable R-Cyclodextrin-Based [2]Rotaxane. J. Am. Chem. Soc. 2008, 130, 11294–11296. Bern a, J.; Goldup, S. M.; Lee, A.-L.; Leigh, D. A.; Symes, M. D.; Teobaldi, G.; Zerbetto, F. Cadiot-Chodkiewicz Active Template Synthesis of Rotaxanes and Switchable Molecular Shuttles with Weak Intercomponent Interactions. Angew. Chem., Int. Ed. 2008, 47, 4392–4396.
DOI: 10.1021/jz100414j |J. Phys. Chem. Lett. 2010, 1, 1776–1780
pubs.acs.org/JPCL
(6)
(7)
(8)
(9)
(10)
(11)
(12)
(13)
(14)
(15)
(16)
(17)
(18)
(19) (20)
(21)
(22)
(23)
(24)
Coutrot, F.; Romuald, C.; Busseron, E. A New pH-Switchable Dimannosyl[c2]Daisy Chain Molecular Machine. Org. Lett. 2008, 10, 3741–3744. Coskun, A.; Friedman, D. C.; Li, H.; Patel, K.; Khatib, H. A.; Stoddart, J. F. A Light-Gated STOP-GO Molecular Shuttle. J. Am. Chem. Soc. 2009, 131, 2493–2495. Balzani, V.; Clemente-Le on, M.; Credi, A.; Ferrer, B.; Venturi, M.; Flood, A. H.; Stoddart, J. F. Autonomous Artificial Nanomotor Powered by Sunlight. Proc. Natl. Acad. Sci. U.S.A. 2006, 103, 1178–1183. Murakami, H.; Kawabuchi, A.; Matsumoto, R.; Ido, T.; Nakashima, N. A Multi-Mode-Driven Molecular Shuttle: Photochemically and Thermally Reactive Azobenzene Rotaxanes. J. Am. Chem. Soc. 2005, 127, 15891–15899. Marlin, D. S.; Cabrera, D. G.; Leigh, D. A.; Slawin, A. M. Z. Complexation-Induced Translational Isomerism: Shuttling through Stepwise Competitive Binding. Angew. Chem., Int. Ed. 2006, 45, 77–83. Kawaguchi, Y.; Harada, A. A Cyclodextrin-Based Molecular Shuttle Containing Energetically Favored and Disfavored Portions in Its Dumbbell Component. Org. Lett. 2000, 2, 1353–1356. Belohradsky, M.; Elizarov, A. M.; Stoddart, J. F. Speed-Controlled Molecular Shuttles. Collect. Czech. Chem. Commun. 2002, 67, 1719–1728. Brough, B.; Northrop, B. H.; Schmidt, J. J.; Tseng, H. R.; Houk, K. N.; Stoddart, J. F.; Ho, C. M. Evaluation of Synthetic Linear Motor-Molecule Actuation Energetics. Proc. Natl. Acad. Sci. U. S.A. 2006, 103, 8583–8588. Phoa, K.; Neaton, J. B.; Subramanian, V. First-Principles Studies of the Dynamics of [2]Rotaxane Molecular Switches. Nano Lett. 2009, 9, 3225–3229. Jang, Y. H.; Hwang, S. G.; Kim, Y. H.; Jang, S. S.; Goddard, W. A. Density Functional Theory Studies of the [2]Rotaxane Component of the Stoddart-Heath Molecular Switch. J. Am. Chem. Soc. 2004, 126, 12636–12645. Deng, W. Q.; Muller, R. P.; Goddard, W. A. Mechanism of the Stoddart-Heath Bistable Rotaxane Molecular Switch. J. Am. Chem. Soc. 2004, 126, 13562–13563. Jang, Y. H.; Goddard, W. A. Stoddart-Heath [2]Rotaxane Molecular Switch Made Simple: A Density Functional Theory Study on Model Junction Devices. J. Phys. Chem. C 2010, 114, 4611–4616. Kim, H.; Goddard, W. A.; Jang, S. S.; Dichtel, W. R.; Heath, J. R.; Stoddart, J. F. Free Energy Barrier for Molecular Motions in Bistable [2]Rotaxane Molecular Electronic Devices. J. Phys. Chem. A 2009, 113, 2136–2143. Darve, E.; Pohorille, A. Calculating Free Energies Using Average Force. J. Chem. Phys. 2001, 115, 9169–9183. Henin, J.; Chipot, C. Overcoming Free Energy Barriers Using Unconstrained Molecular Dynamics Simulations. J. Chem. Phys. 2004, 121, 2904–2914. Rodriguez-Gomez, D.; Darve, E.; Pohorille, A. Assessing the Efficiency of Free Energy Calculation Methods. J. Chem. Phys. 2004, 120, 3563–3578. Chipot, C.; Henin, J. Exploring the Free-Energy Landscape of a Short Peptide Using an Average Force. J. Chem. Phys. 2005, 123, 244906. Darve, E.; Rodríguez-G omez, D.; Pohorille, A. Adaptive Biasing Force Method for Scalar and Vector Free Energy Calculations. J. Chem. Phys. 2008, 128, 144120. Yu, Y. M.; Chipot, C.; Cai, W. S.; Shao, X. G. Molecular Dynamics Study of the Inclusion of Cholesterol into Cyclodextrins. J. Phys. Chem. B 2006, 110, 6372–6378.
r 2010 American Chemical Society
(25)
(26) (27)
(28)
(29)
(30)
(31)
(32)
(33) (34) (35)
(36)
1780
Yu, Y. M.; Cai, W. S.; Chipot, C.; Sun, T. T.; Shao, X. G. Spatial Arrangement of R-Cyclodextrins in a Rotaxane. Insights from Free-Energy Calculations. J. Phys. Chem. B 2008, 112, 5268–5271. Eyring, H. The Activated Complex in Chemical Reactions. J. Chem. Phys. 1935, 3, 107–115. MacKerell, A. D.; Bashford, D.; Bellott, M.; Dunbrack, R. L.; Evanseck, J. D.; Field, M. J.; Fischer, S.; Gao, J.; Guo, H.; Ha, S.; et al. All-Atom Empirical Potential for Molecular Modeling and Dynamics Studies of Proteins. J. Phys. Chem. B 1998, 102, 3586–3616. Pavelites, J. J.; Gao, J. L.; Bash, P. A.; Mackerell, A. D. A Molecular Mechanics Force Field for NADþ, NADH, and the Pyrophosphate Groups of Nucleotides. J. Comput. Chem. 1997, 18, 221–239. Foloppe, N.; MacKerell, A. D. All-Atom Empirical Force Field for Nucleic Acids: I. Parameter Optimization Based on Small Molecule and Condensed Phase Macromolecular Target Data. J. Comput. Chem. 2000, 21, 86–104. Kuttel, M.; Brady, J. W.; Naidoo, K. J. Carbohydrate Solution Simulations: Producing a Force Field with Experimentally Consistent Primary Alcohol Rotational Frequencies and Populations. J. Comput. Chem. 2002, 23, 1236–1243. Strader, M. L.; Feller, S. E. A Flexible All-Atom Model of Dimethyl Sulfoxide for Molecular Dynamics Simulations. J. Phys. Chem. A 2002, 106, 1074–1080. Phillips, J. C.; Braun, R.; Wang, W.; Gumbart, J.; Tajkhorshid, E.; Villa, E.; Chipot, C.; Skeel, R. D.; Kal e, L.; Schulten, K. Scalable Molecular Dynamics with NAMD. J. Comput. Chem. 2005, 26, 1781–1802. Flyvbjerg, H.; Petersen, H. G. Error Estimates on Averages of Correlated Data. J. Chem. Phys. 1989, 91, 461–466. Humphrey, W.; Dalke, A.; Schulten, K. VMD: Visual Molecular Dynamics. J. Mol. Graphics 1996, 14, 33–38. Bergeron, R. J.; Channing, M. A.; Gibeily, G. J.; Pillor, D. M. Disposition Requirements for Binding in Aqueous Solution of Polar Substrates in the Cyclohexaamylose Cavity. J. Am. Chem. Soc. 1977, 99, 5146–5151. Li, J.; Harada, A.; Kamachi, M. Formation of Inclusion Complexes of Oligoethylene and Its Derivatives with R-Cyclodextrin. Bull. Chem. Soc. Jpn. 1994, 67, 2808–2818.
DOI: 10.1021/jz100414j |J. Phys. Chem. Lett. 2010, 1, 1776–1780