Molecular dynamics

Molecular dynamics Julia M. Goodfellow and Mark A. Williams Birkbeck College, University of London, London, UK As molecular simulations become more sophisticated, their ability to provide a detailed understanding of experimental results, in both structural and energetic terms, is increasing. This review describes some of the growing range of problems that are being tackled by molecular simulations, and the technical developments that will aid future progress. Current Opinion in Structural Biology 1992, 2:211-216

Introduction Computer simulation techniques, including restrained and unrestrained molecular dynamics (MD) and related free energy protocols, are becoming ever more popular, no doubt because of the increasing availability of software packages and the decreasing cost of computers. One of the more reassuring aspects of recent studies is the increasing number and breadth of comparisons being made with experimental data. Free energy calculations will play an important role in providing information for such comparisons. Equally reassuring, therefore, are the number of papers addressing free energy estimations from MD simulation and attempting to pro vide realistic error estimates. Standard simulation protocols treat molecules classically but certain problems, such as electron transfer reactions and catalytic mechanisms, clearly require consideration of quantum effects. Some progress is being made in overcoming this diflqcult problem, which is reviewed separately in this section by Warshel (pp 230--236). Because of the heavy computational cost of a detailed simulation of a macromolecule surrounded by solvent, there continues to be an emphasis on methods to increase the speed of calculations using new hardware, in cluding massively parallel machines and novel approximations. In this review, we have tried not only to sum marize the important technical developments but also to acknowledge the contribution of computation to the wide range of topics that are now being studied using MD and related simulation algorithms.

Protein unfolding Although attempts to simulate protein folding are be coming more sophisticated [see a review in this issue by Jemigan (pp248-256)], all-atom simulations of folding are still some way off. It is proving possible to study protein unfolding, however, with the initial focus being on long simulations of polypeptides with 0~-helical confor-

mation at elevated temperatures. Dynamics simulations of 30-residue alanine and glycine helices both in v a c u o and in solution [1"] show that bending of the glycine helix is intrinsic to the sequence whereas bending of the alanine helix results from the insertion of water molecules into an intramolecular hydrogen bond. This water insertion effect is not seen during a long, 250ps, simulation at 300 K but only when the temperature of the simulation is increased to 350 K for 30 ps. These observations support the suggestion that water insertion is an important part of the mechanism of protein unfolding [2,3]. The stability of another helix, an analogue of the ribonuclease S peptide, has been studied for 300 ps at 278 K and 500 ps at 358 K [4°]. At low temperature, the helix is stable, as measured experimentally. At high temperature, main-chain a-helical hydrogen bonds changed into 310-helical hydrogen bonds, which cannot be distinguished from a type-Ill reverse turn in these calculations. These 310 hydrogen bonds may also break, allowing the main-chain O and NH groups to hydrogen bond to solvent molecules. This study is again related to the unfolding mechanism proposed by Sundralingam and Sekharadu, which involves a transition in hydrogen bonding from a-helical to 310 to turn [2].

Membranes and channels Solvent transport along the membrane-spanning polypeptide gramicidin has aroused much interest because of its small size and the dit~culty of obtaining experimental information on membrane proteins. The free energy changes associated with the movement of Na + ions along a model channel have been calculated [5]. The free energy profile of the channel is found to be controlled by the water-water, water-peptide and peptide-peptide hydrogen bounds and not the large electrostatic energy associated with the ion. Helix flexibility is seen to allow local perturbations of the structure as the Na + ion moves through the channel. In a second paper, Roux and Karplus [6] have shown that, unlike the motion of Na + ,

Abbreviations MD--molecular dynamics; NOE--nuclear Overhausereffect.

(~ Current Biology Ltd ISSN 0959-440X

211

212 Theoryand simulation the motion of H20 and K + is diffusive, with the motion of K + being controlled by the diffusion of water. The importance of polypeptide flexibility for the functioning of the channel has also been demonstrated in a time correlation analysis of water motion [7]. Although individual water molecules in the channel are found to have higher frequentT motions than those in bulk, this is not true for the motion of the centre of mass of the water molecules in the channel (the motion of the centre of mass defines translocation through the channel). In a rigid channel, the mobility of water molecules is reduced "along the long axis of the channel whereas flexibility leads to a reduction in the height of the free energy bartiers to motion along the channel and consequently to increased water mobility. The above simulations do not contain explicit lipid molecules, which would normally surround gramicidin within a membrane. Simulation on ensembles of lipid molecules are progressing, however, and two recent studies have made good use of the fact that membranes are large structures composed of many relatively small building blocks [8°°,9°}. Both studies use some form of average external field and random motion to simu late the environment of a small number of central lipid molecules. One of the methods combines MD of the central molecules with these boundary conditions, and seems to offer a convenient scheme for more precise MD simulations of lipid bilayers and polypeptides in mem branes [8"°].

DNA structure, hydration and drug complexes

[14o]. Mthough clearly of the B form and showing many features of crystal structure, it does not show narrowing of the minor groove, which retains well ordered hydration sites. Comparison with simulation i n v a c u o shows that explicit inclusion of solvent is necessary to properly support the major- and minor-groove structure of the DNA helix. Thus, although the hydration pattern of the minor groove is now seen clearly to depend on the width of the minor groove, the origin of minor-groove narrowing remains unclear. In a recent experimental paper, Dickerson and colleagues [13] noted that opposition of 13ii phosphates across the minor groove may provide a clue to the widening of the groove. This hypothesis could be tested in future simulations. The energetics and conformational properties of the important triple-stranded helices have also been studied during the past year [15]. Given the lack of a crystal structure, there are likely to be many more such stud ies, which clearly must take account of the accumulating experimental data [ 16]. Drug-DNA complexes have received much attention, probably as a result of the growing interest in using MD for rational drug design. Simulations have produced models for a spermine~poly(dG-dC) complex, in which binding occurs along the major groove, distamycin bind ing in the minor groove, a dynemicin-DNA intercalation complex and free energy estimates for the binding of netropsin.

Free energy calculations

The long-standing interest in the effect of hydration on the helical conformation of DNA, particularly B-DNA, continues. In B-DNA, a 'spine' of hydration was seen in the minor groove of the central AT-rich region in the crystal structure of a dodecamer [10]. This spine involves water molecules bridging between thymine and adenine bases on opposite strands in adjacent base pairs in a region where the minor groove is narrower than that in the canonical B-DNA conformation. A Monte Carlo simulation on a dodecamer had provided theoretical evidence fbr a well ordered water structure in the wide minor groove of the canonical B form, which was seen to extend beyond the central AT-rich region [11].

The development of techniques for estimating free en ergy differences between similar systems of molecules was one of the most exciting developments in the use of MD during the 1980s [17]. A variety of these techniques have been incorporated into the popular MD software packages and the number and complexiW of their applications to biological problems is increasing rapidly. The numerous applications during the past year include studies on hydrophobic hydration, stability of dimers of deoxyhemoglobin mutants, several studies of inhibitor binding, hydration of protein cavities, partition coefficients and calculation of pKas.

In a recent study, based on a 40 ps simulation of a de camer of fixed conformation, a second hydration pattern was described in the wider regions of the minor groove [12o]. The water molecules were located in the plane of the bases and individually hydrate each base. In the narrower regions, the water molecules were located between adjacent base pairs and bridged across the strands forming the well known spine, whereas in regions of intermediate width, the water molecules were shifted away from the base plane. These results are in agreement with the latest experimental structural data for B-DNA helices [13].

Not only do free energy estimates give better insight into equilibrium properties of molecular interactions than do potential energy estimates, but they should, in principle, lead to better comparison of simulation results with experimental data. Obtaining precise and accurate estimates for relative free energies, however, is difficult: errors can arise from both the statistical nature of the calculation and systematic deficiencies in both the methodology and the force fields. The fact that errors are difficult to estimate means that the success or failure of a single simulation may be difficult to judge simply by comparing calculated and experimental free energy differences.

A 140 ps MD simulation on the dodecamer in which the nucleotide is allowed to move has also been reported

Mthough it is clearly impractical in many cases to re peat extensive calculations on proteins using different

Molecular dynamicsGoodfellow methodologies or different potentials, it is possible to do this for simpler model systems, and thereby quantify some of these errors and find protocols for minimizing them. Studies using the 'slow-growth' method, on a threonine--+alanine mutation in a dipeptide in solution, only give precise free energy differences for simulations >100ps and cannot be obtained by averaging several shorter simulations [18°]. One probable cause of convergence problems in such calculations is mutations being made more rapidly than the local structure can relax around them [19]. Although free energies are path-independent quantities, it is necessary to find mutation pathways with the shortest relaxation times to allow efficient simulation [20]. Indirect pathways can be used to reduce computer time and obtain better accuracy [21o], A new approach to the thermodynamic integration/slowgrowth method uses several configurations instead of one at each integration step [22]. This procedure leads to a more efficient t~se of computer time and a more accurate way of determining the rate of change of the Hamiltonian. Systematic errors in similar dynamically modified window calculations have also been studied [23].

Restrained dynamics The power of restrained dynamics has been demonstrated by the now routine use of simulated annealing protocols in the refinement of X-ray crystallographic structures against the experimentally determined struc ture factors, and in the determination of solution structures of proteins from two-dimensional nuclear Over hauser effect (NOE) NMR data [24]. We shall look briefly at some developments in techniques but not at the large number of applications. Relaxation matrix refinement has been added to the protocols in order to take account of spin-diffusion effects and thus minimize the differences between the calculated and experimentally determined NOE intensities [25]. The build-up of NOE intensities as a function of time has been calculated from MD trajectories, and it is hoped that comparison of the variations in such results that arise from the use of different models with experimentally measured data will provide a good test of current calculations [26]. The use of time-averaged restraints of internal coordi nates (e.g. applying NOE distance restraints to average rather than instantaneous atom separations) has also been introduced recently [27]. This procedure produces a better picture of inherent flexibility and the fine conformational detail in NMR structures compared with those found in an unrestrained simulation [28]. When applied to X-ray structures, this method has revealed an ensemble of conformations, unlike traditional refinement methods, and leads to better agreement with the experimental data [29]. The extent of disorder in protein structures can be determined rapidly by X-ray restrained dynamics (because of the relatively large radius of convergence) but the actual modelling of these disordered regions can be just as successfully accomplished using conventional least squares refinement (presumably at lower cost) [30].

and Williams

Algorithms and architecture Current technology provides us with computers with diverse architectures. Each type has its own advantages and may be able to offer better performance than others for a specific task. These new architectures promise faster or cheaper simulations but at the cost of developing new algorithms in order to realise this potential. Both the range of different architectures and the extensive rewriting of computer programs present problems which have inhibited the general adoption of novel architecture machines. Accepting these problems, it still seems worthwhile to test the possible benefits of parallel architecture machines which, in principle, offer large increases in speed. With this aim in mind, the systolic loop MD algorithm has been extended to include neighbour listings, where the list is distributed over the nodes [31"]. Fincham and Mitchell [31"] also present an excellent review of MD on multiple instruction multiple data machines. A parallel code with a link cell algorithm has also been developed [32]. Most of these codes do not take account of terms in the potential that require knowledge of three (bond angle) or four (torsion angle) atom coordinates. Raine [33"] has written a parallel algorithm, also on a distributed memory machine, specifically for protein dynamics, in which b o n d and torsion angle terms are included. Algorithms have also been developed for a single instruction multiple data machine [34]. These improvements in algorithms and the apparent movement toward distributed memory machines with flexible architecture should lead to a more stable environment for programmers, and hopefully to parallel versions of the common packages. The use of shared memory parallel machines (e.g. minisupercomputers) is now quite common and programs for these machines using both neighbour lists and linked cell algorithms have been assessed (e.g. [35] ),

Extending the range of applications Using current computers, the most extensive simulations that are readily performed are on systems of several thousand atoms for a time of a few hundred picoseconds. Improvements in algorithms, as well as systematic and reliable simplifications of potential functions, are necessary in order to produce better ensemble averages, to extend the application of MD to processes occurring on longer time scales, and to cope with larger systems. The most time consuming operations in an MD simulation are calculation of the distances and forces between non-bonded atoms. The performance of programs is improved by methods that either reduce the number of distances that need to be calculated or reduce the frequency of their calculation. The first of these routes is represented by a newly proposed set of potentials for peptides, in which the interactions of the atomic centres are represented by fewer 'virtual' centres [36"]. The new potentials have been parameterized against those produced by the program CHARMM and the fit appears to

213

214

Theoryand simulation be quite good. Forthcoming comparisons of MD simulations using both the standard and the new potentials will be very interesting. The second route is represented by a generalization of the verlet algorithm which can be used with distance class schemes to reduce computational costs [37°]. In such schemes, the frequency of the calculation of the interaction between pairs is reduced as their separation increases (on the basis that the rate of change of the interaction of distant pairs is lowered). The application of this scheme to the photosynthetic reaction centre demonstrates the power of the method [37"]. Mimicking the effect of explicit solvent atoms is another route to increasing the speed of simulations, especially as the cost of computing water-water interactions in a highly solvated system can be very high. The incorporation of solvation energy and ion-screening effects, previously calculated using the finite difference Poisson-Boltzmann method, into the MD program XPLOR has been described and tested on an alanine dipeptide [38"]. The results are similar to those generated by an explicit solvent simulation but require much less computing time. Belief in such a method is strengthened by a study which shows that conventional free energy MD protocols and finite difference Poisson-Boltzmann calculations produce similar values for the electrostatic free energy differences [39]. A similar technique has been incorporated into the program CHARMM [40-•]. A popular way of mimicking some of the effects of solvent has been the use of a dielectric constant greater than unity [41]. Microscopic simulations of dielectric constants in solvated trypsin, however, have shown that the value of the dielectric constant varies throughout the protein and that it must be calculated including a solvent reaction field or a sufficiently large number of surrounding water molecules [42"°]. This leads us to believe that the use of electrostatic techniques may replace the use of a distance-dependent dielectric constant. One of the advantages of Monte Carlo simulations has been the ability to perform simulations in the grand canonical ensemble, i.e. one in which the number of particles varies during the simulation. This ensemble is important because it allows the calculation of absolute free energies in a convergent manner, for some problems at least. Now, this monopoly has been broken by the development of a grand canonical MD method [43"]. It is unclear as yet how important such simulations will prove for biomolecular research, although a successful Monte Carlo-based simulation of the ionic environment of DNA [44 •] gives us hope that they will prove extremely useful.

Conclusions Molecular dynamics and related methods are being used in a wide variety of applications. Their ability to provide insight into the causes of experimentally observed phe nomena will encourage their routine use by non-specialism, as has already been seen in the field of structure refinement. This trend will give added impetus to efforts

to describe and quantify sources of error in simulations. Free energy calculations are particularly important and their errors the most difficult to assess, making further work on free energy protocols an important target of future research. The routine use of MD for simple systems and the challenge of studying more complex systems will both require faster and cheaper simulations. The falling cost of conventional computers, the stabilization of parallel algorithms and architectures, and improved approximation schemes will all play important roles in future research and development.

References and recommended reading Papers of particular interest, published within the annual period of re view, have been highlighted as: • of special interest ,,• of outstanding interest 1. •

D1CAPUAFM, SWAMINATHAN S, BEVERIDGE DL." Theoretical Evidence for Water Insertion in ~x-Helix Bending: Molecular Dynamics of Gly30 and Ala30 in vacuo and in S01ution. J Am Chem Soc 199l, 113:6145-6155. Simulations at 27°C and 77°C are used to probe the role of" the stability of intramolecular hydrogen bonding in ~x-helices and the role of solvent in folding intermediates. 2.

SUNDARALINGAM M, SEKHARUDU YC: Water-inserted co-Helical Segments Implicate Reverse Turns as Folding Intermediates. Science 1989, 244:1333-1337.

3.

D~CAPUAFM, SWAMINATHANS, BEVERIDGEDL: Theoretical Evidence for Destabilization of an a-Helix by Water Insertion: Molecular Dynamics of Hydrated Decaalanine. J Am Chem Soc 1991, 112:6768-6771.

4. •

TIRADO-RIVESJ, JORGENSEN WL: Molecular Dynamics Sireulations of the Unfolding of an ~.Hehcal Analogue of Ribonuclease A S-peptide in Water. Biochemistry 1991, 30:3864-3871. Simulations at 5°C and 85°C are used to follow the unfolding of a pepride whose stability is well characterized experimentally. The authors show that possible folding intermediates are seen within 200 ps and are thus amenable to study by computer simulation. 5.

Roux B, KARPLUS M: Ion Transport in a Model Gramicidin Channel: Structure and Thermodynamics. Biophys J 1991, 59:961-981.

6.

Roox B, KARPLUS M: Ion Transport in a Gramicidinlike Channel: Dynamics and Mobility. J Phys Chem 1991, 95:48564868.

7.

CHIU S-W, JAKOBSSONE, SUBRAMANIANS, MCCAMMONJA: Timecorrelation Analysis of Simulated Water Motion in Flexible and Rigid Gramicidin Channels. BiophysJ 1991, 60:273-285.

8. ••

DELOOFH, HARVEY SC, SEGRESTJP, PASTOR RW: Mean Field Stochastic Boundary Molecular Simulations of Phospholipid in a Membrane. Biochemistry 1991, 30:209~2113. This study of a lipid membrane combines atomic level MD of a central molecule with stochastic dynamics to represent the surrounding lipid en~ ronment. 9. •

PASTORPW, VENABLE RM, KARPLUS M: Model for the Structure of the Lipid Bilayer. Proc Natl Acad Sci USA 1991, 88:892-4~96. A lipid bilayer is constructed from the trajectories for a single lipid molecule obtained from brownian dynamics simulations. 10.

DREWHR, DICKERSON RE: Structures of a DNA Dodecamer: lIl. Geometry of Hydration. J Mol Biol 1981, 151:535-556.

11.

SUBRAMANIANPS, RAVISHKANKERG, BEVERIDGE DL: Theoretical Considerations on the 'Spine of Hydration' in the Mi-

Molecular dynamics G o o d f e l l o w a n d W i l l i a m s nor Groove of d(CGCGAATTCGCG).d(GCGCITAAGCGC): Monte Carlo C o m p u t e r Simulation. Proc Natl Acad Sci USA 1988, 85:1836-1840.

28.

P ~ DA, KOLLMAN PA: Are Time-averaged Restraints Necessary for Nuclear Magnetic Resonance Refinement? A Model Study for DNA. J Mol Biol 1991, 220:457-479.

CHUPRIINAVP, HE1NEMA~,~ U, NumstaMOV AA, ZIELENKIEWICZ p, D1CKERSON RE, SAENGERW: Molecular Dynamics Simulation of the Hydration Shell of a B-DNA D e c a m e r Reveals Two Main Types of Minor Groove Hydration D e p e n d i n g on Groove Width. Proc Natl Acad Sci USA 1991, 8:593-597. This standard MD simulation of DNA focuses on the pattern of hydration in the minor groove, which has been the subject of much experimental and theoretical debate.

29.

GROSP, VAN GUNSTERENWF, HOL WGJ: Inclusion of Thermal Motion in Crystallographic Structures by Restrained Molecular Dynamics. Science 1990, 249:1149-1152.

30.

KURWANJ, C)SAPAY K, BURLEY SK, BRUNGER A, HENDRICKSON WA, KARPLUSM: Exploration o f Disorder in Protein Structures by X-ray Restrained Molecular Dynamics. Proteins 1991, 10:34(~358.

13.

31. •o

12. •

GRE.ZESKOWIAKK, YANAGI K, PRIVE GG, DICKERSON RE: The Structure of B-helical C-G-A-T-C-G-A-T-C-G and Comparison with C-C-A-A-C-G-T-T-G-G: the Effect of Base Pair Reversals. j Biol Chem 1991, 266:8861-8883.

14. •

SWAM1NATHANS, RAV1SHANKERG, BEVERIDGE DL: Molecular Dynamics of B-DNA Including Water and Counterions: a 140ps Trajectory for d(CGCGAATrCGCG) Based on the GROMOS Force Field. J Am Chem Soc 1991, 113:5027-5040. A fairly long and detailed simulation of B-DNA, tocussing on the hydration structure in the minor groove.

15.

VAN VUJMEN HWF, RAME GL, PE'ITiTr BM: A Study of Model Energetics and Conformationai Properties of Polynucleotide Triplexes. Biopolymers 1990, 30:517-532.

16.

MACAYARF, GILBERT DE, MALEK S, SINSHEIMER JS, FEIGON J: Structure and Stability of X-G-C Mismatches in t h e Third Strand of Intramolecular Triplexes. Science 1991, 254:270-274.

17.

BEVERIITW, E DL, DICAPUAFM: Free Energy via Molecular Simulation: Applications to Chemical and Biomolecular Systems. Annu Rev Biophys Biophys Chem 1989, 18:431-492.

18. •

MrrCHELLMJ, MCCAMMONJA: Free Energy Difference Calculations by Thermodynamic Integration: Difficulties in Obtaining a Precise Value. J Comput ~ m 1991, 12:271 275. A detailed study of the problems encountered when short-duration ( < 100 ps) simulations are used to estimate free energy differences us ing the now standard slow-growth method. 19.

MAZORM, PETITIT BM: Convergence of t h e Chemical Potential in Aqueous Simulations. Mol Simul 1991, 6:1-4.

20.

BERENDSENHJC: Incomplete Equilibration: a Source of Error in Free Energy Calculation. In Proteins: Structure Dynamics and Design. Edited by Renugopalakrishnan V, Carey PR, Smith ICP, Huang SG, Storer AC. Leiden: ESCOM; 1991:384-392.

21. •

MARKAE, VAN GUNSTERENWF, BERNENDSEN HJC: Calculation of Relative Free Energy via Indirect Pathways. J Chem Pl~ys 1991, 94:380~3816. This study shows that the correct choice of transformation pathway can reduce the simulation time and reduce errors by allowing the system to remain close to equilibrium at all stages. 22.

STRAATSMA TP, MCCAMMONJA: Multiconfiguration T h e r m o d y namic Integration. J Chem P ~ s 1991, 95:1175 1188.

23.

WOODRH, MUHLBAUERWCF, THOMPSON PT: Systematic Errors In Free Energy Perturbation Calculations due to a Finite Sample of Corifigurational Space: Sample Size Hysteresis, J Phys Chem 1991, 95:667045675.

24.

I3RI~INGERAT, KARPLUS M: Molecular Dynamics Simulations with Experimental Restraints. Acc Chem Res 1991, 24:54-61.

25.

NILGESM, HABAZETrLJ, BRI2NGER AT, HOLAK TA: Relaxation Matrix Refinement of t h e Solution Structure of Squash Trypsin Inhibitor. J Mol Biol 1991, 219:499-510.

26.

WrrHKAJM, SWAMINATHANS, BEVERIDGE DL, BOLTON PH: Time D e p e n d e n c e of Nuclear Overhauser Effects of Duplex DNA from MolecuLar Dynamics Trajectories. J A m Chem Soc 1991, 113:5041 5049.

27.

TORDAAE, SCHEEK RM, VAN GUNSTEREN 'q~F: Time-averaged Nuclear Overhauser Effect Distance Restraints Applied to Tendamistat. J Mol Biol 1990, 214:223-235.

FINCHAMD, MITCEIEL1.PJ: Multicomputer Molecular Dynamics Simulations Using Distributed Neighbour Lists. Mol Simul 1991, 7:135-153. This paper describes a parallel MD algorithm, written in Fortran, which includes a conventional neighbour list distributed over the processors. This is well worth reading for the detailed review of current dynamics algorithms for multiple instruction multiple data machines. 32.

PINCHESMRS, TILDESLEY DJ, SMITH W: Large Scale Molecular Dynamics on Parallel C o m p u t e r s Using the Link-cell Algorithm. Mol Simul 1991, 6:51-87.

33. •o

PAINEARC: Systolic Loop Methods for Molecular Dynamics Simulation Generalised for Macromolecules. Mol Simu11991, 7:59-69. This paper describes how ,systolic loop methods for parallel machines may be implemented for the study of protein dynamics. 34.

WINDEMUTHA, SCHULTENK: Molecular Dynamics on t h e Connection Machine. Mol Simul 1991, 6:353 361.

35.

SKEEtRD: Macromolecular Dynamics on a Shared Memory Multiprocessor..1 Comput Chem 1991, 12:175-179.

36. t lEAD-GORDON T, BROOKS CL 11I: Virtual Rigid Body Dynam• ics. Biopolymers 1991, 31:77-100. This details a new set of potentials for proteins which use fewer 'atoms' to represent each molecule. In conjunction with new methods for rep resenting solvent effects, this method could enable very much longer simulations to be performed. 37, •

GRUBMI)LLERH, HELLER H, WINDEMUTH A, SCHULTEN K: Generalized Verlet Algorithm for Efficient Molecular Dynamics Simulations with Long-range Interactions. Mol Simul 1991, 6:121 142. This paper describes a modified verier algorithm for use in conjunction with a distance class scheme, which allows considerable speed up of simulations of very large systems,

38. .,,

SHARP K: Incorporating Solvent and Ion Screening into Molecular Dynamics Using the Finite-difference Poisson-Boltzmann Method. J Comput Chem 1991, 12:454468. The author describes the incorporation into XPLOR of an electrostatic method for estimating the effects of solvation; if this method proves to give results similar to explicit solvation for a range of systems, it will find wide use. 39.

JEAN-CHARLESA, NICHOI.LS A, SHARP K, HONIG H, TEMPCZYK A, HENDR1CKSON TF, STILL WC: Electrostatic Contributions to Solvation Energies: Comparison o f Free Energy Perturbation and C o n t i n u u m Calculations. J Am Chem Soc 1991, 113:1454--1455.

40. ..

ZAUHARRJ: The Incorporation of Hydration Forces Determ i n e d by C o n t i n u u m Electrostatics into Molecular Mechanics Simulations. J Comput Chem 1991, 12:575-583. This paper describes a method similar to that of Sharp [38.e], though, as yet, it has only been incorporated into a minimization routine. 41.

42. •.

FRITSCHV, WESTHOF E: Minimization and Molecular Dynamics Studies of Guanosine and Z-DNA Modified by N-2Acetylaminofluorene. J Comput C.hem 1991, 12:147-166.

KING G, LEE FS, WARSHEL A: Microscopic Simulations of Macroscopic Dielectric Constants of Solvated Proteins. J Chem Phys 1991, 95:43664377. This study demonstrates that the effective dielectric constant of a protein is locally variable. This implies that any scheme that aims to mimic

215

216

Theory and simulation solvent by the use of a distance-dependent dielectric constant is inherently inaccurate. 43.

(~AGINT, PETTITI"BM: Grand Molecular Dynamics: a Method for Open Systems. Mol Simul 1991, 6:5-26. This paper extends the method of MD to systems with variable numbers of particles, and opens the way to the calculation of absolute free energies for some systems. 44. •

JAYARAMB, BeVEmDGE DL: Grand Canonical Monte Carlo Simulations on Aqueous Solutions of NaCI and NaDNA: Excess Chemical Potentials and Sources of Nonideality in Elec-

trolyte and Polyelectrolyte Solutions. J Phys Cbem 1991, 95:2506-2516. This paper considers ionic strength effects in a novel manner and il lustrates the use of grand canonical methods to produce data that can help evaluate the simulation against experimental data.

JM Goodfellow and MA williams, Department of Crystallography, Birkbeck College, Malet Street, University of London, London WClE 7HX, UK.