📄 xmd-5.html
字号:
<!DOCTYPE HTML PUBLIC "-//W3C//DTD HTML 3.2 Final//EN"><HTML><HEAD> <META NAME="GENERATOR" CONTENT="SGML-Tools 1.0.9"> <TITLE>XMD - Molecular Dynamics Program: Theory</TITLE> <LINK HREF="xmd-6.html" REL=next> <LINK HREF="xmd-4.html" REL=previous> <LINK HREF="xmd.html#toc5" REL=contents></HEAD><BODY><A HREF="xmd-6.html">Next</A><A HREF="xmd-4.html">Previous</A><A HREF="xmd.html#toc5">Contents</A><HR><H2><A NAME="s5">5. Theory</A></H2><H2><A NAME="ss5.1">5.1 Molecular Dynamics</A></H2><P><P>In the technique of molecular dynamics the motion for every atom iscalculated. This algorithm proceeds as follows.<P><OL><LI>(1) The initialpositions and velocities of every atom are specified.</LI><LI>(2) Using theinteratomic potentials, the forces on these atoms are calculated.</LI><LI>(3)Using these forces, the atomic positions and velocities are advancedthrough a small time interval (called a time step). These new positionsand velocities become new input to step (2), and when steps 2 and 3 arerepeated, each repetition constitutes an additional time step.</LI></OL>Typically a molecular dynamics simulation will comprise thousands ofsuch time steps, each time step corresponding to fraction of apicosecond (1/10^1^2 seconds). This algorithm is essentially anintegration of the Newtonian equations of motion over time to yield theparticle positions and velocities.<P><H2><A NAME="ss5.2">5.2 Boundary Conditions</A></H2><P>A molecular dynamics simulation is limited in thenumber of atoms one can simualate,typically XMD simulations are between 5,000 and 100,000 atoms.With such a small number of atoms, it is not possible tosimulate bulk material unless one uses<B>repeating boundary conditions</B>.With repeating boundary conditions, one typically conductsa simulation within a box. Atoms that pass out one wall of thebox pass back in to the box throught the opposite wall. Thusat no time is an atom outside the box. This way there is nofree surface, and the system simulates the bulk. While in theory,any parallelpiped (that is, a cell with triclinic symmetry)can be used as the repeating cell, at thistime XMD is restricted to using a rectangular box.The repeating boundary conditions can be turned off in any directionindependently of the others, in this way one can simulation anan infinitely repeating solid (the bulk),an infinitely repeating sheet,an infinitely repeating "wire" ora finite fragment of material.There is more information on this under the <A HREF="xmd-8.html#command-box">BOX command</A>.<P><H2><A NAME="ss5.3">5.3 Adiabatic and Isothermal Simulations</A></H2><P><P>In thermodynamics, systems are divided into two kinds, adiabatic andisothermal. Adiabatic systems are thermally isolated from the externalworld, and no heat will flow in or out. The result is that the totalenergy (the sum of potential and kinetic energies) will always remainconstant. Isothermal systems are systems which maintain a constanttemperature through contact with a heat bath. Both kinds of systems canbe simulated with molecular dynamics. See the section on <A HREF="xmd-6.html#implementation-temp">Temperature Control</A>for more information.<P><H2><A NAME="ss5.4">5.4 Time Step Size</A></H2><P><P>As explained in the above section on molecular dynamics, the forces andvelocities at one step in time are used to calculate the resultingpositions at the next step. The time difference between these twoadjacent steps is called the time step size. The time step size must bespecified by the user. Ideally one would like to use the largest timestep size possible, for this way one would simulate the greatest timepossible. However the molecular dynamics integration algorithm becomesunstable at large time step sizes. This can be understood byconsidering a single particle in a one dimensional harmonic well.Assume that this particle is to one side of the center and that itsinitial velocity is zero. This particle would experience a force towardthe center of the well. The product of the force with the time stepsize squared would yield the displacement of the particle. A small timestep would advance the particle closer to the well center. A too largetime step would overshoot the center of the well, and could even placethe particle on the opposite side higher than its starting position.This would immediately cause a problem. The total energy of theparticle will have increased, and subsequent steps at this large timestep size would only increase the energy more. This clearly violatesthe conservation of energy.<P><P>In practice this violation of energy conservation is exactly the clueused to determine that a time step is too large. One wants to find thelargest time step that will maintain the conservation of energy. Forhints on how to implement a search for the optimum time step see thesection <A HREF="xmd-12.html#techniques">Techniques</A>.<P><H2><A NAME="ss5.5">5.5 Metropolis Algorithm</A></H2><P>The Metropolis algorithm is a type of Monte Carlo technique used for twopurposes (1) to calculate averages of a system at finite temperatures(thermodynamic averages) and (2) to simulate annealing of a system ofparticles. This algorithm varies the atomic positions at random and theproperties of interest for each configuration are calculated. In asimple Monte Carlo (not Metropolis) scheme for calculating thermodynamicaverages random atomic configurations are sampled and each contributionto the average is weighted by the factor exp(-E/kT) where E is thesystem potential energy, k is Boltzmann's constant, and T is thetemperature. However purely random atomic configurations will consistmostly of high energy (and therefore improbable) states, and theresulting weighting factor will be so small as to render insignificantthe contribution to the average. Thus the bulk of the computationaltime will be spent calculating useless numbers. The alternative is touse the Metropolis method for sampling configuration space. This methodsamples phase space via a random walk. Steps in this random walk areonly taken if the energy remains low enough to contribute a significantamount to the average. Steps in this random walk which raise the energywill be rejected unless a randomly generated number uniform on theinterval [0..1] is less than the factor exp(-dE/kT) (where dE is theenergy of the new state minus the energy of the original). Steps whichlower the energy are always accepted. Averages are taken for everyconfiguration (if a new configuration is rejected then the oldconfiguration re- contributes to the average). This selective randomwalk results in the averages taken being weighted by the factorexp(-E/kT). The controlling parameters are the temperature and the jumpsize. The jump size determines how far the atoms move with each randomwalk.<P><H2><A NAME="theory-pressure"></A> <A NAME="ss5.6">5.6 Constant Pressure Algorithms </A></H2><P>Sometimes one wishes to simulate system which can change itsvolume in response to the combination of applied external pressureand the system's own internal stress. XMD has two different methodsfor accomplishing this. The first, Andersen's method [3], isdesigned reproduce the thermodynamics properties of an isobaricadiabatic (NPE) ensembled.Andersen's method is described inmore detail below.<P>The second method is " pressure clamp " method(see the PRESSURE CLAMP command).This method was developed for XMD.It automatically rescales the system volume by the amount neededto balance an applied pressure. It is designed for determiningthe equilbrium volume as a function of applied temperature andpressure. Most typically, it is used to determine the thermalexpansion of a simulated material.<P><P><H3>Andersen's Constant Pressure Method</H3><P>This is a method which maintains an external pressure by imposing anartificial force on every atom. This has the advantage of not creatinga localized interface between the simulation and the external pressuresource. It is implemented with the PRESSURE command. When usingAndersen's constant pressure, you must assign a mass to the externalsystem. This mass determines how quickly the volume of the systemresponds to the pressure. Without such a mass the system volume changeswould be instantaneous - clearly a non-physical situation. The exactvalue of the mass is probably unimportant. A good initial guess wouldbe to set the external system mass to the total internal mass (the massof all the atoms).<P><P>Andersen's constant pressure method is not suitable for determiningthermal expansion. One reason is that is takes much longer toreach equlibrium volume than the " pressure clamp "method. The second reasonis that this methodintroduces a fictitious pressure that is proportional to the systemkinetic energy.The value of the external pressure is given by the simple gas law,<P><BLOCKQUOTE><CODE><PRE> P = nKT/V </PRE></CODE></BLOCKQUOTE><P>where<P>N is the total number of atoms,kB is Boltzmann's constant,T the temperature andV the average system volume.<P><H3>The "Pressure Clamp" method</H3><P>The "Pressure Clamp" works analogously to the temperatureclamp. The volume of the system box is automatically re-sized to maintaina constant pressure. This works as follows. The user specifiesa pressure that she would like to maintain, the default is zero.In addition, the user specifies a Bulk Modulus for the system.This can be the exact Bulk Modulus or an approximation. Then,for each molecular dynamics step, the system stress is calculated.By combining the system stress and Bulk Modulus, a change inthe box size can be estimated which produces the desired pressure,<BLOCKQUOTE><CODE><PRE>(dLi/Li)/Cstep = B / (3*Si)</PRE></CODE></BLOCKQUOTE>where Li is the box size in the i'th direction,dLi is the change in that box size,Si is the stress in the i'th direction (the Voight stress),B is the Bulk Modulus,and Cstep is a parameter which controls how quickly the algorithmconverges on a suitable box size. By default Cstep is 33, butthis is under user control.<P><H3>Stress Calculation and the Virial Term</H3><P>XMD by default does not include the virial term<BLOCKQUOTE><CODE><PRE> Sij = mass * Vi * Vj</PRE></CODE></BLOCKQUOTE>in the stress calculation. This is a controversialissue, but I justify this from the fact that a perfectharmonic solid will not thermally expand, but thevirial term predicts it will.<P>The virial term can be included in the stress calculation.See the STRESS THERMAL command for details.<P><P><H3>Box Geometry</H3><P>When using either of these two constant pressure algorithms, onehas the option of letting each box direction (x,y,z) changeindependently, or forcing the system to change uniformly.When changing uniformly, the ratio of x, y and z box sizes remains thesame, and consequently a system which starts out as a cubiclattice will remain a cubic lattice, and a tetragonal latticewill remain tetragonal and maintain its initial c/a ratio.<P><P><P><P><P><HR><A HREF="xmd-6.html">Next</A><A HREF="xmd-4.html">Previous</A><A HREF="xmd.html#toc5">Contents</A></BODY></HTML>⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -