《几种常见的优化方法ppt课件》由会员分享，可在线阅读，更多相关《几种常见的优化方法ppt课件（29页珍藏版）》请在金锄头文库上搜索。
1、几种常见的优化方法电子结构几何机构函数稳定点最小点Taylor 展开：V(x) = V(xk) + (x-xk)V(xk) +1/2 (x-xk)2 V(xk)+. 当x是3N个变量的时候， V(xk)成为3Nx1的向量，而V(xk)成为3Nx3N的矩阵，矩阵元如：Hessian11.一阶梯度法2. a. Steepest descendentSk = -gk/|gk|directiongradient知道了方向，如何确定步长呢？最常用的是先选择任意步长l，然后在计算中调节用体系的能量作为外界衡量标准，能量升高了则逐步减小步长。robust, but slow最速下降法 2最陡下降法（SD）
2、3b. Conjugate Gradient （CG) 共轭梯度第k步的方向标量Usually more efficient than SD, also robust不需要外界能量等作为衡量量利用了上一步的信息42。二阶梯度方法这类方法很多，最简单的称为Newton-Raphson方法，而最常用的是Quasi-Newton方法。Quasi-Newton方法： use an approximation of the inverse Hessian.Form of approximation differs among methods牛顿-拉夫逊法 BFGS methodBroyden-Fletc
3、her-Golfarb-ShannoDFP methodDavidon-Fletcher-Powell5Molecular dynamics 分子动力学分子动力学HistoryIt was not until 1964 that MD was used to study a realistic molecular system, in which the atoms interacted via a Lennard-Jones potential.After this point, MD techniques developed rapidly toencompass diatomic spe
4、cies, water (which is still thesubject of current research today!), small rigid molecules,flexible hydrocarbons and now even macromoleculessuch as proteins and DNA.These are all examples of continuous dynamicalsimulations, and the way in which the atomic motion iscalculated is quite different from t
5、hat in impulsivesimulations containing hard-core repulsions.6What can we do with MD Calculate equilibrium configurational properties in a similarfashion to MC. Study transport properties (e.g. mean-squared displacement anddiffusion coefficients). MD in the NVT, NpT and NpH ensembles The united atom
6、approximation Constraint dynamics and SHAKE Rigid body dynamics Multiple time step algorithmsExtend the basic MD algorithm7Impulsive molecular dynamics1.Dynamics of perfectly hard particles can be solved exactly, but process becomes involved for many part (N-body problem). 2.Can use a numerical sche
7、me that advances the system forward in time until a collision occurs. 3.Velocities of colliding particles (usually a pair!) then recalculated and system put into motion again. 4. Simulation proceeds by fits and starts, with a mean time between collisions related to the average kinetic energy of the
8、particles. 5.Potentially very efficient algorithm, but collisions between particles of complex shape are not easy to solve, and cannot be generalised to continuous potentials. 8Continuous time molecular dynamics1. By calculating the derivative of a macromolecular forcefield, we can find the forces o
9、n each atomas a function of its position.2. Require a method of evolving the positions of the particles in space and time to produce a true dynamical trajectory.3. Standard technique is to solve Newtons equations ofmotion numerically, using some finite difference scheme,which is known as integration
10、.4. This means that we advance the system by some smalltime step t, recalculate the forces and velocities, and thenrepeat the process iteratively.5. Provided t is small enough, this produces an acceptable approximate solution to the continuous equations of motion.9Example of integrator for MD simula
11、tion One of the most popular and widely used integrators isthe Verlet leapfrog method: positions and velocities ofparticles are successively leap-frogged over each otherusing accelerations calculated from force field. The Verlet scheme has the advantage of high precision(of order t4), which means th
12、at a longer time step canbe used for a given level of fluctuations. The method also enjoys very low drift, provided anappropriate time step and force cut-off are used.r(t+D Dt)=r(t)+v(t+D Dt/2)D Dtv(t+D Dt/2)=v(t-D Dt/2)+a(t+D Dt/2)D Dt10Other integrators for MD simulations Although the Verlet leapf
13、rog method is not particularlyfast, this is relatively unimportant because the timerequired for integration is usually trivial in comparison tothe time required for the force calculations. The most important concern for an integrator is that itexhibits low drift, i.e. that the total energy fluctuate
14、s aboutsome constant value. A necessary (but not sufficient)condition for this is that it is symplectic. Crudely speaking, this means that it should be timereversible (like Newtons equations), i.e. if we reverse themomenta of all particles at a given instant, the systemshould trace back along its pr
15、evious trajectory.11Other integrators for MD simulations The Verlet method is symplectic, but methods such aspredictor-corrector schemes are not. Non-symplectic methods generally have problems withlong term energy conservation. Having achieved low drift, would also like the energyfluctuations for a
16、given time step to be as low as possible. Always desirable to use the largest time step possible. In general, the trajectories produced by integration willdiverge exponentially from their true continuous paths dueto the Lyapunov instability. However, this does not concern us greatly, as the thermals
17、ampling is unaffected expectation values unchanged.12Choosing the correct time step 1. The choice of time step is crucial: too short and phase space is sampled inefficiently, too long and the energy will fluctuate wildly and the simulation may become catastrophically unstable (“blow up”). 2. The ins
18、tabilities are caused by the motion of atoms being extrapolated into regions where the potential energy is prohibitively high (e.g. atoms overlapping). 3. A good rule of thumb is that when simulating an atomic fluid, the time step should be comparable to the mean time between collisions (about 5 fs
19、for Ar at 298K). 4. For flexible molecules, the time step should be an order of magnitude less than the period of the fastest motion (usually bond stretching: CH around 10 fs so use 1 fs).13For classic MD, there could be many tricks to speed up calculations, all centering around reducing the effort
20、involved in the calculation of the interatomic forces, as this is generally much more time-consuming than integration.For example1.Truncate the long-range forces: charge-charge, charge-dipole2.Look-up tablesFor first principles MD, as forces are evaluated from quantum mechanics, we are only concerne
21、d with the time-step.14Because the interactions are completely elastic andpairwise acting, both energy and momentum areconserved. Therefore, MD naturally samples from themicrocanonical or NVE ensemble.As mentioned previously, the NVE ensemble is not very useful for studying real systems. We would li
22、ke to be able to simulate systems at constant temperature or constant pressure.The simplest MD, like verlet method, is a deterministic simulation technique for evolving systems to equilibrium by solving Newtons lawsnumerically.15MD in different thermodynamic ensembles In this lecture, we will discus
23、s ways of using MD tosample from different thermodynamic ensembles, whichare identified by their conserved quantities. Canonical (NVT) Fixed number of particles, total volume and temperature.Requires the particles to interact with a thermostat. Isobaric-isothermal (NpT) Fixed number of particles, pr
24、essure and temperature. Requiresparticles to interact with a thermostat and barostat. Isobaric-isenthalpic (NpH) Fixed number of particles, pressure and enthalpy. Unusual, butrequires particles to interact with a barostat only.16Advanced applications of MD We will then study some more advanced MD me
25、thodsthat are designed specifically to speed up, or makepossible, the simulation of large scale macromolecularsystems. All these methods share a common principle: they freezeout, or decouple, the high frequency degrees of freedom. This enables the use of a larger time step withoutnumerical instabili
26、ty. These methods include: United atom approximation Constraint dynamics and SHAKE Rigid body dynamics Multiple time step algorithms17Revision of NVE MD Lets start by revising how to do NVE MD. Recall that we calculated the forces on all atoms from thederivative of the force field, then integrated t
27、he e.o.m.using a finite difference scheme with some time step t. We then recalculated the forces on the atoms, andrepeated the process to generate a dynamical trajectoryin the NVE ensemble. Because the mean kinetic energy is constant, the averagekinetic temperature TK is also constant. However, in t
28、hermal equilibrium, we know thatinstantaneous TK will fluctuate. If we want to sample fromthe NVT ensemble, we should keep the statisticaltemperature constant.18Extended Lagrangians There are essentially two ways to keep the statisticaltemperature constant, and therefore sample from the trueNVT ense
29、mble. Stochastically, using hybrid MC/MD methods Dynamically, via an extended LagrangianWe will describe the latter method in this lecture An extended Lagrangian is simply a way of including adegree of freedom which represents the reservoir, andthen carrying out a simulation on this extended system.
30、 Energy can flow dynamically back and forth from thereservoir, which has a certain thermal inertia associatedwith it. All we have to do is add some terms to Newtonsequations of motion for the system.19Extended Lagrangians The standard Lagrangian L is written as the difference ofthe kinetic and poten
31、tial energies: Newtons laws then follow by substituting this into theEuler-Lagrange equation: Newtons equations and Lagrangian formalism areequivalent, but the latter uses generalised coordinates.20Canonical MD So, our extended Lagrangian includes an extracoordinate , which is a frictional coefficie
32、nt that evolves in time so as to minimise the difference between the instantaneous kinetic and statistical temperatures. The modified equations of motion are:The conserved quantity is the Helmholtz free energy.(modified form of Newton II)21Canonical MD By adjusting the thermostat relaxation time tT
33、(usuallyin the range 0.5 to 2 ps) the simulation will reach anequilibrium state with constant statistical temperature TS. TS is now a parameter of our system, as opposed to themeasured instantaneous value of TK which fluctuatesaccording to the amount of thermal energy in the systemat any particular
34、time. Too high a value of tT and energy will flow very slowlybetween the system and the reservoir (overdamped). Too low a value of tT and temperature will oscillate aboutits equilibrium value (underdamped). This is the Nos-Hoover thermostat method.22Canonical MD There are many other methods for achi
35、eving constanttemperature, but not all of them sample from the trueNVT ensemble due to a lack of microscopic reversibility.We call these pseudo-NVT methods, and they include: Berendsen methodVelocities are rescaled deterministically after each step so that the system is forced towards the desired te
36、mperature Gaussian constraintsMakes the kinetic energy a constant of the motion by minimising the least squares difference between the Newtonian and constrained trajectories These methods are often faster, but only converge on thetrue canonical average properties as O(1/N).23Isothermal-isobaric MD W
37、e can apply the extended Lagrangian approach to simulations at constant pressure by simply adding yet another coordinate to our system. We use , which is a frictional coefficient that evolves intime to minimise the difference between the instantaneous pressure p(t), measured by a virial expression,
38、and the pressure of an external reservoir pext. The equations of motion for the system can then be obtained by substituting the modified Lagrangian into the Euler-Lagrange equations. These now include two relaxation times: one for thethermostat tT, and one for the barostat tp.24Isothermal-isobaric M
39、DThe is known as the Nos-Hoover method (Melchionnatype) and the equations of motion are:25Isothermal-isobaric MD26Constraint dynamicsfreeze the bond stretching motions of the hydrogens (or any other bond, in principle). We apply a set of holonomic constraints to the system, which are just algebraic
40、equations connecting the coordinates of the individual atoms. These constraints contain information about the bond lengths. The most widely used algorithm for implementing these constraints is called SHAKE. This works on the principle of applying fictitious forces to all the particles such that they
41、 follow appropriately constrained trajectories. There is also a variant of SHAKE called RATTLE.27We can calculate what the constraint forces are (obviously they must be directed parallel to the bonds) by introducing Lagrange multipliers lij:These multipliers can be determined by substituting the mod
42、ified expressions for Newton II into our Verlet integration scheme and imposing the constraints that:28 This leads to a set of simultaneous equations for lij whichcan be recast into a matrix form and solved by invertingthe matrix. For larger molecules, this is very time consuming and soan iterative
43、procedure is used (hence the “shaking” of themolecule as each constraint is enforced in turn). Angle and higher order constraints can easily beincorporated into the algorithm by enforcing fictitiousconstraint bonds between atoms such that the anglebetween them is held constant. However, there are problems with convergence of thealgorithm if the molecule is overconstrained.29
