分子模拟论文PPT课件
《分子模拟教程》PPT课件_OK
17
GEMC 的配分函数
对于原子系统，位型（构型）的配分函数
Q(N,V1,V2 ,T )
N
V N1 1
(V
V1 ) N N1
复杂流体如：胶体悬浮液、高分子溶液、 表面活性剂溶胶等。
❖ 超临界过程研究中的应用
37
❖ 在多相催化研究中的应用：
➢ 对催化剂进行表征 ➢ 表面吸附与脱附过程及表面性质的模拟 ➢ 催化剂表面反应机理的模拟
38
5.3 介观层次材料的计算机模拟
结构是多层次、多尺度的，并且不仅要研究平衡结构， 还要研究结构随时间的演变。所谓结构，就是粒子在空 间有规律的分布。（胡英院士）
N
exp[U (s N )]
14
Metropolis GCMC algorithm
产生巨正则系综的马尔可夫链的过程涉及到典型的、 三种不同的随机移动：
❖ Attempt to move a particle (just like canonical MC)
❖ Attempt to create a particle ❖ Attempt to delete a particle
19
各种随机移动的概率：
Pmove (m n) min{1, exp[Unm / kT]}
Pvolume (m n) min{1,
change
(V1n ) N1 (V - Vቤተ መጻሕፍቲ ባይዱn ) N-N1 (V1m ) N1 (V - V1m ) N-N1
exp[Unm / kT]}
分子模拟PPT—第二章 力场
2. 蒙特卡罗随机采样法
分子内能
分子的能量
Etot Eele Evib Erot Etrans
分子的振动自由度 = 3N-6 = 3N-6 分子的振动能 非线性 线性
Evib Ebonds Eangles Etorsions Ecross
分子间相互作用能
1. 范德华能
AMBER
• / • "Amber" refers to two things: a set of molecular mechanical force fields for the simulation of biomolecules ;and a package of molecular simulation programs which includes source code and demos. The current version of the code is Amber version 10, which is distributed by UCSF subject to a licensing agreement described below. • Amber is now distributed in two parts: AmberTools and Amber10. AmberTools 1.2 & Amber 10 is now available! • Amber was originally developed under the leadership of Peter Kollman, and Version 9 is dedicated to his memory.
第二代力场
• 第二代力场的形式远较上述的经典力场复杂，需 要大量的力常数。其设计的目的为能精确地计算 分子的各种性质、结构、光谱、热力学特性、晶 体特性等资料。其力常数的推导除引用大量的实 验数据外，还参照精确的量子计算的结果。尤其 适用于有机分子或不含过渡金属元素的分子系统。 • 第二代力场因其参数的不同，包括CFF91、 CFF95、PCFF与MMFF93等。
分子模拟PPT—第四章 分子动力学模拟原理
2
rb
or
S2
(r
)
1
r b
2
rb
b
模拟中的控制
• 温度控制 简单方法 – 每一步调节动量都使动能逼近期望值
• 更新位置和力 • 迭代方法使得 v v f
f T Tcalc
– 缺点
• “运行方程” 是不可逆的，系统不遵循细致平衡 • 不属于任何有明确定义的系综
复杂要求 (多时间尺度):
当出现多时间尺度 e.g., 不同质量的混合粒子, 溶剂 聚合体, 柔性和刚性的共存分子体系等等, t 的选取 必须依照体系中动力学变化最快的成分或模型。
积分步长的选取
• 积分步长应小于系统中最快运动周期的1/10。 • 以氩原子的分子动力学计算为例：
U (r) d U (r) d [4 (12 12 / r13 6 6 / r7 )] 4 (156 12 / r14 42 6 / r8)
S(r) 0 r b
b
cutoff 方案
2. Switching
1
ra
S(r) 1 y(r)22 y(r) 3 a r b
0
rb
这里
y(r)
r2 b2
a2 a2
ab
3. Shifting
S1(r)
1
r b
2
dr
dr
U (rmin ) 57.14 / 2
k 4 2 2
mAr mAr 19.981.662 1022 (g)
mAr mAr
57.12 0.24 6.9446 1014 erg
《分子模拟方法》课件
加速研发进程
分子模拟可以大大缩短药 物研发、材料合成等领域 的实验周期，降低研发成 本。
揭示微观机制
通过模拟，可以揭示分子 间的相互作用机制和反应 过程，有助于深入理解物 质的性质和行为。
分子模拟的发展历程
经典力学模拟
基于牛顿力学，适用于 较大分子体系，但精度
较低。
量子力学模拟
适用于小分子体系，精 度高，但计算量大，需
详细描述
利用分子模拟方法，模拟小分子药物与生物大分子（如蛋白质、核酸等）的相 互作用过程，探究药物的作用机制和药效，为新药研发提供理论支持。
高分子材料的模拟研究
总结词
研究高分子材料的结构和性能，优化 材料的设计和制备。
详细描述
通过模拟高分子材料的结构和性能， 探究高分子材料的物理和化学性质， 优化材料的设计和制备过程，为新材 料的研发提供理论指导。
分子动力学方法需要较高的计算资源和 精度，但可以获得较为准确的结果，因 此在计算化学、生物学、材料科学等领
域得到广泛应用。
介观模拟的原理
介观模拟是一种介于微观和宏观之间的模拟方 法，通过模拟一定数量的粒子的相互作用和演 化来研究介观尺度的结构和性质。
介观模拟方法通常采用格子波尔兹曼方法、粒 子流体动力学等方法，适用于模拟流体、表面 、界面等介观尺度的问题。
分子模拟基于量子力学、经典力 学、蒙特卡洛等理论，通过建立 数学模型来描述分子间的相互作
用和运动。
分子模拟可以用于药物研发、材 料科学、环境科学等领域，为实 验研究和工业应用提供重要支持
。
分子模拟的重要性
01
02
03
预测分子性质
通过模拟，可以预测分子 的性质，如稳定性、溶解 度、光谱等，为实验设计 和优化提供指导。
分子模拟课件
ME346–Introduction to Molecular Simulations–Stanford University–Winter2007 Handout1.An Overview of Molecular SimulationWei Caic All rights reservedSeptember26,2005Contents1Examples of molecular simulations1 2What consists of a molecular dynamics?2 3Newton’s equation of motion2 4A simple numerical integrator:Verlet algorithm5 1Examples of molecular simulationsPictures and movies of modern large scale molecular simulations can be downloaded from the course web site /∼me346(click“Course Materials”).Exam-ples include plastic deformation in metals containing a crack,phase transformation induced by shock waves,nano-indentation of metalfilm,laser ablation on metal surface,etc.In each of these examples,molecular simulations not only produce beautiful pictures,but they also offer new understanding on the behavior of materials at the nano-scale.Many times such understanding provides new insights on how the materials would behave at the micro-scale and nano-scale.With unique capabilities to offer,molecular simulation has become a widely used tool for scientific investigation that complements our traditional way of doing science: experiment and theory.The purpose of this course is to provide thefirst introduction to molecular level simulations—how do they work,what questions can they answer,and how much can we trust their results—so that the simulation tools are no longer black boxes to us.Students will develop a set of simulation codes in Matlab by themselves.(A simple Matlab tutorial written by Sigmon[5]may help you refresh your memory.)We will also give tutorial and exercise problems based on MD++,a molecular simulation package used in my research group.MD++will be used for illustration in class as well as tools to solve homework problems.2What consists of a molecular dynamics?In a molecular simulation,we view materials as a collection of discrete atoms.The atoms interact by exerting forces on each other and they follow the Newton’s equation of motion. Based on the interaction model,a simulation compute the atoms’trajectories numerically. Thus a molecular simulation necessarily contains the following ingredients.•The model that describes the interaction between atoms.This is usually called the interatomic potential:V({r i}),where{r i}represent the position of all atoms.•Numerical integrator that follows the atoms equation of motion.This is the heart of the ually,we also need auxiliary algorithms to set up the initial and bound-ary conditions and to monitor and control the systems state(such as temperature) during the simulation.•Extract useful data from the raw atomic trajectory pute materials properties of interest.Visualization.We will start our discussion with a simple but complete example that illustrates a lot of issues in simulation,although it is not on the molecular scale.We will simulate the orbit of the Earth around the Sun.After we have understood how this simulation works(in a couple of lectures),we will then discuss how molecular simulations differ from this type of “macro”-scale simulations.3Newton’s equation of motionThroughout this course,we limit our discussions to classical dynamics(no quantum mechan-ics).The dynamics of classical objects follow the three laws of Newton.Let’s review the Newton’s laws.Figure1:Sir Isaac Newton(1643-1727England).First law:Every object in a state of uniform motion tends to remain in that state of motion unless an external force is applied to it.This is also called the Law of inertia.Second law:An object’s mass m,its acceleration a,and the applied force F are related byF=m a.The bold font here emphasizes that the acceleration and force are vectors.Here they point to the same direction.According to Newton,a force causes only a change in velocity(v)but is not needed to maintain the velocity as Aristotle claimed(erroneously).Third law:For every action there is an equal and opposite reaction.The second law gives us the equation of motion for classical particles.Consider a particle (be it an atom or the Earth!depending on which scale we look at the world),its position is described by a vector r=(x,y,z).The velocity is how fast r changes with time and is also a vector:v=d r/dt=(v x,v y,v z).In component form,v x=dx/dt,v y=dy/dt,v z=dz/dt. The acceleration vector is then time derivative of velocity,i.e.,a=d v/dt.Therefore the particle’s equation of motion can be written asd2r dt2=Fm(1)To complete the equation of motion,we need to know how does the force F in turn depends on position r.As an example,consider the Earth orbiting around the Sun.Assume the Sun isfixed at origin(0,0,0).The gravitational force on the Earth is,F=−GmM|r|2·ˆe r(2)whereˆe r=r/|r|is the unit vector along r direction.|r|=√x2+y2+z2is the magnitudeof r.m is the mass of the Earth(5.9736×1024kg),M is mass of the Sun(1.9891×1030kg), and G=6.67259×10−11N·m2/kg2is gravitational constant.If we introduce a potential field(the gravitationalfield of the Sun),V(r)=−GmM|r|(3)Then the gravitational force can be written as minus the spatial derivative of the potential,F=−dV(r)d r(4)In component form,this equation should be interpreted as,F x=−dV(x,y,z)dxF y=−dV(x,y,z)dyF z=−dV(x,y,z)dzNow we have the complete equation of motion for the Earth:d2r dt2=−1mdV(r)d r(5)in component form,this reads,d2x dt2=−GMx|r|3d2y dt2=−GMy|r|3d2z dt2=−GMz|r|3In fact,most of the system we will study has a interaction potential V(r),from which the equation of motion can be specified through Eq.(5).V(r)is the potential energy of the system.On the other hand,the kinetic energy is given by the velocity,E kin=12m|v|2(6)The total energy is the sum of potential and kinetic energy contributions,E tot=E kin+V(7) When we express the total energy as a function of particle position r and momentum p=m v, it is called the Hamiltonian of the system,H(r,p)=|p|22m+V(r)(8)The Hamiltonian(i.e.the total energy)is a conserved quantity as the particle moves.To see this,let us compute its time derivative,dE tot dt =m v·d vdt+dV(r)d r·d rdt=md rdt·−1mdV(r)d r+dV(r)d r·d rdt=0(9)Therefore the total energy is conserved if the particle follows the Newton’s equation of motion in a conservative forcefield(when force can be written in terms of spatial derivative of a potentialfield),while the kinetic energy and potential energy can interchange with each other.The Newton’s equation of motion can also be written in the Hamiltonian form,d r dt =∂H(r,p)∂p(10)d p dt =−∂H(r,p)∂r(11)Hence,dH dt =∂H(r,p)∂r·d rdt+∂H(r,p)∂p·d pdt=0(12)Figure2:Sir William Rowan Hamilton(1805-1865Ireland).4A simple numerical integrator:Verlet algorithmA dynamical simulation computes the particle position as a function of time(i.e.its trajec-tory)given its initial position and velocity.This is called the initial value problem(IVP). Because the Newton’s equation of motion is second order in r,the initial condition needs to specify both particle position and velocity.To solve the equation of motion on a computer,thefirst thing we need is to discretize the time.In other words,we will solve for r(t)on a series of time instances t ually the time axis is discretized uniformly:t i=i·∆t,where∆t is called the time step.Again,the task of the simulation algorithm is tofind r(t i)for i=1,2,3,···given r(0)and v(0).The Verlet algorithm begins by approximating d2r(t)/dt2,d2r(t) dt2≈r(t+∆t)−2r(t)+r(t−∆t)∆t2(13)Thusr(t+∆t)−2r(t)+r(t−∆t)∆t2=−1mdV(r(t))dr(14)r(t+∆t)=2r(t)−r(t−∆t)−∆t·1m·dV(r(t))dr(15)Therefore,we can solve r(t+∆t)as along as we know r(t)and r(t−∆t).In other words,if we know r(0)and r(−∆t),we can solve for all r(t i)for i=1,3,4,···.Question:Notice that to start the Verlet algorithm we need r(0)and r(−∆t).However, the initial condition of a simulation is usually given in r(0)and v(0).How do we start the simulation when this initial condition is specified?References1.Alan and Tildesley,Computer Simulation of Liquids,(Oxford University Press,1987)pp.71-80.2.Frenkel and Smit,Understanding Molecular Simulation:From Algorithms to Applica-tions,2nd ed.(Academic Press,2002).3.Bulatov and Cai,Computer Simulation of Dislocations,(Oxford University Press,2006)pp.49-51.ndau and Lifshitz,Mechanics,(Butterworth-Heinemann,1976).5.Kermit Sigmon,Matlab Primer,Third Edition,/Matlab/matlab-primer.pdfME346–Introduction to Molecular Simulations–Stanford University–Winter2007 Handout2.Numerical Integrators IWei Caic All rights reservedJanuary12,2007Contents1Which integrator to use?1 2Several versions of Verlet algorithm2 3Order of accuracy3 4Other integrators4 1Which integrator to use?In the previous lecture we looked at a very simple simulation.The Newton’s equation of motion is integrated numerically by the original Verlet algorithm.There are many algorithms that can be used to integrate an ordinary differential equation(ODE).Each has a its own advantages and disadvantages.In this lecture,we will look into more details at the numerical integrators.Our goal is to be able to make sensible choices when we need to select a numerical integrator for a specific simulation.What questions should we ask when choosing an integrator?Obviously,we would like something that is both fast and accurate.But usually there are more issues we should care about.Here are a list of the most important ones(Allen and Tildesley1987):•Order of accuracy•Range of stability•Reversibility•Speed(computation time per step)•Memory requirement(as low as possible)•Satisfying conservation law•Simple form and easy to implementItem1and2(especially2)are important to allow the use of long time steps,which makes the simulation more efficient.2Several versions of Verlet algorithmThe Verlet algorithm used in the previous lecture can be written as,r(t+∆t)=2r(t)−r(t−∆t)+a(t)·∆t2(1)wherea(t)=−1mdV(r(t))d r(t)(2)While this is enough to compute the entire trajectory of the particle,it does not specify the velocity of the particle.Many times,we would like to know the particle velocity as a function of time as well.For example,the velocity is needed to compute kinetic energy and hence the total energy of the system.One way to compute velocity is to use the followingapproximation,v(t)≡d r(t)dt≈r(t+∆t)−r(t−∆t)2∆t(3)Question:if the initial condition is given in terms of r(0)and v(0),how do we use the Verlet algorithm?Answer:The Verlet algorithm can be used to compute r(n∆t)for n=1,2,···provided r(−∆t)and r(0)are known.The task is to solve for r(−∆t)given r(0)and v(0).This is achieved from the following two equations,r(∆t)−r(−∆t)2∆t=v(0)(4)r(∆t)−2r(0)+r(−∆t)∆t2=−1mdV(r(t))d r(t)t=0(5)Exercise:express r(−∆t)and r(∆t)in terms of r(0)and v(0).Therefore,r(−∆t)and r(∆t)are solved together.During the simulation,we only know v(t)after we have computed(t+∆t).These are the(slight)inconveniences of the(original) Verlet algorithm.As we will see later,in this algorithm the accuracy of velocity calculation is not very high(not as high as particle positions).Several modifications of the Verlet algorithms were introduced to give a better treatment of velocities.The leap-frog algorithm works as follows,v(t+∆t/2)=v(t−∆t/2)+a(t)·∆t(6)r(t+∆t)=r(t)+v(t+∆t/2)·∆t(7) Notice that the velocity and position are not stored at the same time slices.Their values are updated in alternation,hence the name“leap-frog”.For a nice illustration see(Allen and Tildesley1987).This algorithm also has the draw back that v and r are not known simultaneously.To get velocity at time t,we can use the following approximation,v(t)≈v(t−∆t/2)+v(t+∆t/2)2(8)The velocity Verlet algorithm solves this problem in a better way.The algorithm reads,r(t+∆t)=r(t)+v(t)∆t+12a(t)∆t2(9)v(t+∆t/2)=v(t)+12a(t)∆t(10)a(t+∆t)=−1mdV(r(t+∆t))d r(t+∆t)(11)v(t+∆t)=v(t+∆t/2)+12a(t+∆t)∆t(12)Exercise:show that both leap-frog and velocity Verlet algorithms gives identical results as the original Verlet algorithm(except forfloating point truncation error).3Order of accuracyHow accurate are the Verlet-family algorithms?Let’s take a look at the original Verletfirst. Consider the trajectory of the particle as a continuous function of time,i.e.r(t).Let us Taylor expand this function around a given time t.r(t+∆t)=r(t)+d r(t)dt∆t+12d2r(t)dt2∆t2+13!d3r(t)dt3∆t3+O(∆t4),(13)where O(∆t4)means that error of this expression is on the order of∆t4.In other words, (when∆t is sufficiently small)if we reduce∆t by10times,the error term would reduce by 104times.Replace∆t by−∆t in the above equation,we have,r(t−∆t)=r(t)−d r(t)dt∆t+12d2r(t)dt2∆t2−13!d3r(t)dt3∆t3+O(∆t4),(14)Summing the two equations together,we have,r(t+∆t)+r(t−∆t)=2r(t)+d2r(t)dt2∆t2+O(∆t4)≡2r(t)+a(t)∆t2+O(∆t4)(15)r(t+∆t)=2r(t)−r(t−∆t)+a(t)∆t2+O(∆t4)(16) Therefore,the local accuracy for the Verlet-family algorithms is to the order of O(∆t4). This is quite accurate for such a simple algorithm!The trick here is that,by adding t+∆t and t−∆t terms,the third order term cancels out.Notice that the accuracy analysis here deals with the local error.That is,assuming r(t)is exact,how much error we make when computing r(t+∆t).In practice,the total error of computed r at time t is the accumulated (and possibly amplified)effect of all time steps before time t.Because for a smaller timestep ∆t more steps is required to reach a pre-specified time t,the global order of accuracy is always local order of accuracy minus one.For example,the global order of accuracy for Verlet algorithms is3.Since all Verlet-family algorithms generate the same trajectory,they are all O(∆t4)in terms of particle positions.But they differ in velocity calculations.How accurate is the velocity calculation in the original Verlet algorithm?Let us subtract Eq.(13)and(14),r(t+∆t)−r(t−∆t)=2v(t)∆t+13d3r(t)dt3∆t3+O(∆t4),(17)v(t)=12∆tr(t+∆t)−r(t−∆t)−16d3r(t)dt2∆t3+O(∆t3)(18)v(t)=r(t+∆t)−r(t−∆t)2∆t+O(∆t2)(19)Thus the velocity is only accurate to the order of O(∆t2)!This is much worse than O(∆t4) of the positions.A reason for the lose of accuracy is that the velocities are computed by taking the difference between two large numbers,i.e.r(t+∆t)and r(t−∆t),and then divide by a small number.1Intuitively,higher order algorithms allows the use of bigger time steps to achieve the same level of accuracy.However,the more important determining factor for the choice of time steps is usually not the order of accuracy,but instead is the algorithm’s numerical stability (see next section).4Other integratorsBefore closing this lecture,let us briefly look at two other well known integrators.Gear’s predictor-corrector method is a high-order integrator.It keeps track of several time deriva-tives of the particle position.Definer0(t)≡r(t)(20) 1Strictly speaking,Eq.(17)is valid when we have exact values for r(t+∆t),r(t)and r(t−∆t).But in doing the order-of-accuracy analysis,we can only assume that we have exact values for r(t)and r(t−∆t) and we know that we are making an error in r(t+∆t).Fortunately,the error we make in r(t+∆t)is also O(∆t4).Thus Eq.(17)still holds.r1(t)≡∆t d r(t)dt(21)r2(t)≡∆t22!d2r(t)dt2(22)r3(t)≡∆t33!d3r(t)dt3(23)···Sincer(t+∆t)=r(t)+∆t d r(t)dt+∆t22!d2r(t)dt2+∆t33!d3r(t)dt3+O(∆t4)=r(t)+r1(t)+r2(t)+r3(t)+O(∆t4)r1(t+∆t)=r1(t)+2r2(t)+3r3(t)+O(∆t4)r2(t+∆t)=r2(t)+3r3(t)+O(∆t4)r3(t+∆t)=r3(t)+O(∆t4)we can predict the values of r and r i at time t+∆t simply based on their values at time t,r P(t+∆t)=r(t)+r1(t)+r2(t)+r3(t)r P1(t+∆t)=r1(t)+2r2(t)+3r3(t)r P2(t+∆t)=r2(t)+3r3(t)r P3(t+∆t)=r3(t)Notice that we have not even calculated the force at all!If we compute the force at time t+∆t based on the predicted position r P(t+∆t),we will obtain the second derivative of r at time t+∆t—let it be r C(t+∆t).Mostly likely,r C(t+∆t)will be different from the predicted value r P(t+∆t).Thefinal results are obtained by adding a correction term to each predictor,r n(t+∆t)=r Pn (t+∆t)+C n[r C2(t+∆t)−r P2(t+∆t)](24)The Gear’s6th order scheme(uses r0through r5)takes the following C n values,C0=3 16C1=251 360C2=1C3=11 18C4=1 6C5=1 60The global accuracy is of the order O(∆t)6,hence the method is called the Gear’s6th order algorithm.Exercise:show that the local truncation error of the above scheme is of the order of O (∆t 7).In a way,the Runge-Kutta method is similar to predictor-corrector methods in that the final answer is produced by improving upon previous predictions.However,the Runge-Kutta method requires multiple evaluation of forces.The Runge-Kutta method is most conveniently represented when the equation of motion is written as a first order differential equation,˙η=g (η)≡ω∂H (η)∂η(25)whereη≡(r ,p )(26)andω= 0I−I 0 (27)and I is the 3×3identity matrix and H is the Hamiltonian.The fourth-order Runge-Kutta method goes as (from ),k 1=∆t ·g (η(t ))k 2=∆t ·g η(t )+12k 1 k 3=∆t ·g η(t )+12k 2 k 4=∆t ·g (η(t )+k 3)η(t +∆t )=η(t )+16k 1+13k 2+13k 3+16k 4(28)Exercise:show that the local truncation error of the above scheme is of the order of O (∆t 5)and the global accuracy is of the order of O (∆t 4).Suggested reading1.(Alan and Tildesley 1987)pp.71-84.2.J.C.Butcher,The numerical analysis of ordinary differential equations:Runge-Kutta and general linear methods ,(Wiley,1987),pp.105-151.ME346–Introduction to Molecular Simulations–Stanford University–Winter2007 Handout3.Symplectic ConditionWei Caic All rights reservedJanuary17,2007Contents1Conservation laws of Hamiltonian systems1 2Symplectic methods4 3Reversibility and chaos5 4Energy conservation7 1Conservation laws of Hamiltonian systemsThe numerical integrators for simulating the Newton’s equation of motion are solvers of ordinary differential equations(ODE).However,the Newton’s equation of motion that can be derived from a Hamiltonian is a special type of ODE(not every ODE can be derived from a Hamiltonian).For example,a Hamiltonian system conserves it total energy,while an arbitrary ODE may not have any conserved quantity.Therefore,special attention is needed to select an ODE solver for Molecular Dynamics simulations.We need to pick those ODE solvers that obeys the conservation laws of a Hamiltonian system.For example,if an ODE solver predicts that the Earth will drop into the Sun after100cycles around the Sun,it is certainly the ODE solver to blame and we should not immediately start worrying about the fate of our planet based on this prediction.A conservation law is a(time invariant)symmetry.A Hamiltonian system possess more symmetries than just energy conservation.Ideally the numerical integration scheme we pick should satisfy all the symmetries of the true dynamics of the system.We will look at these symmetries in this and following sections.When discussing a Hamiltonian system,it is customary to use letter q to represent the particle position(instead of r).Recall that the equation of motion for a system with Hamiltonian H(q,p)is,˙q=∂H(q,p)∂p(1)˙p=−∂H(q,p)∂q,(2)where˙q≡d q/dt,˙p≡d p/dt.If the system has N particles,then q is a3N-dimensional vector specifying its position and p is a3N-dimensional vector specifying its momentum. (In the above Earth orbiting example,N=1.)Notice that in terms of q and p,the equation of motion is a set of(coupled)first-order differential equations(whereas F=m a is second order).The equation can be reduced to an even simpler form if we define a6N-dimensional vectorη≡(q,p)(3) Then the equation of motion becomes,˙η=ω∂H(η)∂η(4)whereω=0I−I0(5)and I is the3×3identity matrix.The6N-dimensional space whichη≡(q,p)belongs to is called the phase space.The evolution of system in time can be regarded as the motion of a pointηin the6N-dimensional phase space following thefirst-order differential equation(4).(a)(b)Figure1:(a)Motion of a point(q,p)in phase space is confined to a shell with constant energy.(b)Area of a phase space element(shaded)is conserved.Because the Newton’s equation of motion conserves total energy,the motion of a point in phase space must be confined to a subspace(or hyperplane)with constant energy,i.e.a constant energy shell,as illustrated in Fig.1(a).Think about all the points in the constant energy shell.As time evolves,the entire volume that these points occupy remains the same —the points simply rearrange within the constant energy shell.Let us now consider the evolution of a small element in phase space over time,as illus-trated in Fig.1(b).Think of this as the ensemble of trajectories of many Molecular Dynamics simulations with very similar initial conditions.We will show that the area enclosed by this element remains constant,even though the element inevitably experiences translation anddistortion.1Let the element at time t be a hypercube around pointη,whose area is,|dη|=|d q|·|d p|=dq1···dq3N dp1···dp3N(6) For a pointηat time t,letξbe its new position at time t+δt.In the limit ofδt→0,wehave,ξ=η+˙ηδt=η+ω∂H(η)∂ηδt(7)Because the above equation is true only in the limit ofδt→0,we will ignore any higher order terms(e.g.δt2)in the following discussions.Let M be the Jacobian matrix of the transformation fromηtoξ,i.e.,M≡∂ξ∂η=1+ω∂2H(η)∂η∂ηδt(8)For clarity,we write M explicitly in terms of q and p,M=1+0I−I0·∂2H∂q∂q∂2H∂q∂p∂2H∂p∂q∂2H∂p∂p·δt(9)=1+∂2Hδt∂2Hδt−∂2H∂q∂qδt1−∂2H∂q∂pδt(10)The area of the new element is related to the the determinant of the Jacobian matrix,|dξ|=|dη|·|det M|=|dη|·(1+O(δt2))(11) Because thefirst order term ofδt in det M vanishes,we can show that the area of the element remains constant for an arbitrary period of time(i.e.forever).2This is an important property of a Hamiltonian system.Because the area of any element in phase space always remains constant,the evolution of phase space points is analogous to that in an incompressiblefluid.Exercise:Write down the explicit expression for the Jacobian matrix M for problem of the Earth orbiting around the Sun.The Hamiltonian dynamics has even more symmetries.Notice that the transpose of M is,M T=1−∂2H∂η∂ηωδt(12)1It then follows that the area enclosed by an arbitrarily shaped element in the phase space also remains constant.2To show that area remains constant after afinite period∆t,we only need to divide the time interval into N subintervals,each withδt=∆t/N.The area change per subinterval is of the order of1/N2.The accumulated area change over time period∆t is then of the order of1/N,which vanishes as N→∞.we haveM ωM T = 1+ω∂2H ∂η∂ηδt ω 1−∂2H ∂η∂ηωδt =ω+ω∂2H ∂η∂ηωδt −ω∂2H ∂η∂ηωδt +O (δt 2)(13)=ω+O (δt 2)(14)(15)Therefore the Jacobian matrix M thus satisfies the so called symplectic condition (up to O (δt 2)),M ωM T =ω(16)Again,we can show that the symplectic condition is satisfied for an arbitrary period of time.To be specific,let ξbe the point in phase space at time t ;at time 0this point was at η.ξcan then be regarded as a function of t and η(initial condition),i.e.ξ=ξ(η,t ).Define M as the Jacobian matrix between ξand η,we can show that M satisfies the symplectic condition Eq.(16).3Obviously,the Hamiltonian dynamics is also time reversible.This means that if ηevolves to ξafter time t .Another phase space point η that has the same q as ηbut with reversed p (momentum)will exactly come to point ξafter time t .To summarize,the dynamics of a system that has a Hamiltonian has the following symmetry properties:•Conserves total energy•Conserves phase space area (incompressible flow in phase space)•Satisfies symplectic condition Eq.(16)•Is time reversibleIdeally,the numerical integrator we choose to simulate Hamiltonian system should satisfy all of these symmetries.In the following,let us take a look at how specific integrators are doing in this regard.In fact,if an algorithm satisfies the symplectic condition,it automatically conserves phase space area 4,is reversible,and should have good long-term energy conserva-tion properties.This is why symplectic methods are good choices for simulating Hamiltonian systems.2Symplectic methodsWhat do we mean if we ask whether a numerical integrator satisfies the symplectic condition or not?For a numerical integrator,the system moves forward in time by discrete jumps.Let3Again,we devide the time period t into N subintervals,each with δt =∆t/N .Let M i be the Jacobian matrix from time (i −1)δt to iδt .We have M = N i =1M i .We can show that M T ωM =ωin the limit of N →∞.4Since M T ωM =ω,det ω=det(M T ωM )=(det M )2det ω.Thus det M =±1,which means area conservation.η(i)=(q(i),p(i))be the place of the system in the phase space at time step i,i.e.t i=i∆t. At the next time step,the algorithm brings the system to pointη(i+1)=(q(i+1),p(i+1)). Let M be the Jacobian matrix that connectsη(i)withη(i+1).The method is symplectic if MωM T=ω.The good news is that the Verlet-family algorithms are symplectic.To show this is the case,let us consider the leap-frog algorithm.Let q(i)correspond to q(i·∆t)and p(i) correspond to p((i−12)·∆t).Then the leap-frog algorithm is,p(i+1)=p(i)−∂V∂q(i)∆t(17)q(i+1)=q(i)+p(i+1)m∆t=q(i)+1mp(i)−∂V∂q(i)∆t∆t(18)Thus,M=∂q(i+1)∂q(i)∂q(i+1)∂p(i)∂p(i+1)(i)∂q(i+1)(i)=1−1m∂2V∂q(i)∂q(i)∆t2∆tm−∂2V∂q(i)∂q(i)∆t1(19)Thus M is symplectic.Exercise:show that for M defined in Eq.(19),MωM T=ω(symplectic)and det M=1 (area preserving).Exercise:derive the Jacobian matrix M for velocity Verlet algorithm and show that it is symplectic and area preserving.3Reversibility and chaosBecause the Verlet-family algorithms(original Verlet,leap-frog and velocity Verlet)are sym-plectic,they are also time reversible.Nonetheless,it is instructive to show the time reversibil-ity explicitly.If we run the leap-frog algorithm in reverse time,then the iteration can be written as,˜q(i+1)=˜q(i)+˜p(i)m∆t(20)˜p(i+1)=˜p(i)−∂V∂q(i+1)∆t(21)(22)where i is the iteration number,˜q(i)=q(−i∆t),˜p(i)=p((−i−12)∆t).It is not difficult toshow that if˜q(i)=q(i+1)and˜p(i)=−p(i+1),then˜q(i+1)=q(i)and˜p(i+1)=−p(i).Exercise:show that the leap-frog algorithm is reversible.Exercise:show that the velocity Verlet algorithm is reversible.Exercise:show that the Forward Euler method:q(i+1)=q(i)+∆t·p(i)/m,p(i+1)=p(i)−∆t·∂V/∂q(i),is not reversible and not symplectic.Since every integration step using the Verlet algorithm is time reversible,then in principle an entire trajectory of Molecular Dynamics simulation using the Verlet algorithm should be time reversible.To be specific,let the initial condition of the simulation be q(0)and p(0). After N steps of integration,the system evolves to point q(N)and p(N).Time reversibility means that,if we run the reverse Verlet iteration starting from q(N)and−p(N),after N steps we should get exactly to point q(0)and−p(0).In practice,however,we never have perfect reversibility.This is due to the combined effect of numerical round offerror and the chaotic nature of trajectories of many Hamiltonian systems.The chaotic nature of many Hamiltonian system can be illustrated by considering the two trajectories with very close initial conditions.Let the distance between the two points in phase space at time0be d(0).At sufficiently large t,the two phase space points will diverge exponentially,i.e.d(t)∼d(0)eλt,whereλis called the Lyapunov exponent and this behavior is called Lyapunov instability.Given that we can only represent a real number on a computer withfinite precision,we are in effect introducing a small(random) perturbation at every step of the integration.Therefore,sooner or later,the numerical trajectory will deviate from the analytical trajectory(if the computer had infinite precision) significantly5.While the analytical trajectory is reversible,the one generated by a computer is not.This is why there exist a maximum number of iterations N,beyond which the original state cannot be recovered by running the simulation backwards.This by itself does not present too big a problem,since we seldom have the need to recover the initial condition of the simulation in this way.However,this behavior does mean that we will not be able to follow the“true”trajectory of the original system forever on a computer.Sooner or later, the“true”trajectory and the simulated one will diverge significantly.This is certainly a problem if we would like to rely entirely on simulations to tell us where is our satellite or space explorer.However,this usually does not present a problem in Molecular Simulations, in which we are usually not interested in the exact trajectories of every atom.We will return to this point in future lectures when we discuss Molecular Dynamics simulations involving many atoms.5Notice that the analytical trajectory here is not the same as the“true”trajectory of the original system. The analytical trajectory is obtained by advancing the system in discrete jumps along the time axis.。
《分子模拟设计》课件
目录
• 分子模拟设计概述 • 分子模拟设计的基本方法 • 分子模拟设计的应用领域 • 分子模拟设计的挑战与展望 • 分子模拟设计案例分析
01
CATALOGUE
分子模拟设计概述
定义与特点
定义
分子模拟设计是指利用计算机模 拟技术，对分子结构和性质进行 预测和设计的过程。
蒙特卡洛模拟
总结词
基于概率统计的模拟方法
详细描述
蒙特卡洛模拟是一种基于概率统计的模拟方法，通过随机抽样和统计计算来获 得系统的性质。该方法适用于模拟复杂系统，能够考虑系统的随机性和不确定 性。
分子力学模拟
总结词
基于势能函数的模拟方法
详细描述
分子力学模拟是一种基于势能函数的模拟方法，通过势能函数来描述分子间的相互作用和分子结构。该方法适用 于快速计算分子的结构和性质，常用于药物设计和材料科学等领域。
材料的界面行为等多个方面。
高分子材料的模拟设计有助于缩短新材料研发周期、 降低研发成本，提高新材料开发的成功率。
高分子材料的模拟设计是利用分子模拟技术对 高分子材料的结构和性质进行预测和优化的一 种方法。
通过模拟高分子材料的结构和性质，可以预测材 料的力学性能、热性能、电性能等，从而优化材 料的设计和制备工艺。
在生物大分子模拟中，研究人员可以使用分子模拟设计来研究蛋白质、 核酸和糖等生物大分子的结构和动力学性质。
这有助于理解这些大分子在细胞中的功能和相互作用的机制，以及与疾 病相关的生物大分子的异常行为。
04
CATALOGUE
分子模拟设计的挑战与展望
计算资源的限制
计算资源不足
高性能计算机和计算集群的资源有限，难以满足 大规模分子模拟的需求。
最新分子模拟技术导论教学讲义ppt课件
本章要求
Xi’an University of Science & Technology
教学目的
演示： 让学生了解分子模拟技术的最新进展与应用
教学要求
掌握 分子模拟技术中常用的方法原理； 精通 HyperChem 软件的操作； 探索 对一些反应过程进行分子动力学模拟尝试； 延伸 利用分子模拟技术设计防晒剂 （课外科技活动）
Dept. of Chemical Science & Engineering
Email: Jansweili@ Phone: 029—85583997
2. 分子模拟技术的计算方法
Xi’an University of Science & Technology
量子力学方法
量子力学方法是基于量子力学的分子模拟，它借助计算分子结 构中各微观参数，如电荷密度、键序、轨道、能级等与性质的 关系，设计出具有特定性能的新分子。它们的共同点是对电子 的相互作用采用量子化学的知识进行描述，而不是采用经验性 的势能函数来表示，这种方法有很强的理论基础。
著名的从头计算程序有系列Gaussian程序
Gaussian 98
Gaussian 2003等
ChemComp
Dept. of Chemical Science & Engineering
Email: Jansweili@ Phone: 029—85583997
2. 分子模拟技术的计算方法
Xi’an University of Science & Technology
半经验方法
是对从头计算中的许多积分采用经验参数替代的简化方法，所使用的 经验参数是通过对实验数据的拟合得到的。半经验方法采用了价电子 近似，假定分子中各原子的内层电子可以看作对分子不极化的原子实 的一部分，而只处理价电子，这样进一步减少了计算时间。主要用于 计算构象能与结构的Ｘ射线结果分析，以此分析平衡态性质
《分子模拟教程》课件
人工智能和机器学习技术将在分子模拟中发挥越 来越重要的作用，例如用于优化模拟参数、预测 性质等。
多尺度模拟
目前分子模拟主要集中在原子或分子级别，未来 将进一步发展多尺度模拟方法，将微观尺度和宏 观尺度相结合，以更全面地理解物质性质和行为 。
跨学科融合
分子模拟将与生物学、医学、材料科学等更多学 科领域进行交叉融合，为解决实际问题提供更多 可能性。
环境科学
在环境科学领域，分子模拟可用于研究污 染物在环境中的迁移转化机制，为环境保 护提供理论依据。
THANKS.
分子动力学模拟的常见算法
Verlet算法
一种基于离散时间步长的算法，用于计算分子位置和速度。
leapfrog算法
一种常用的分子动力学模拟算法，具有数值稳定性和计算效率高的特 点。
Parrinello-Rahman算法
一种基于分子力场的算法，可以用于模拟大尺度分子体系的运动。
Langevin动力学算法
材料科学
通过模拟材料中分子的运动和相互作 用，可以研究材料的力学、热学和电 学等性质，为材料设计和优化提供依 据。
03
Monte Carlo模拟
Monte Carlo模拟的基本概念
随机抽样
Monte Carlo模拟基于随 机抽样的方法，通过大量 随机样本的统计结果来逼 近真实结果。
概率模型
Monte Carlo模拟建立概 率模型，模拟系统的状态 变化和行为。
通过模拟药物分子与靶点分子的相互作用，预测 药物活性并优化药物设计。
材料科学
研究材料中分子的结构和性质，预测材料的物理 和化学性质。
生物大分子模拟
模拟生物大分子的结构和动力学行为，如蛋白质 、核酸等，有助于理解其功能和性质。
分子模拟PPT—第一章 概论
对1998 年诺贝尔化学奖 划时代的评价
瑞典皇家科学院的评价：
“ ··量子化学已发展成为广大化学家都能使用 · 的工具，将化学带入一个新时代 — 实验 与理论能携手协力揭示分子体系的性质。 化学不再是一门纯实验科学了”
I2 在 光 滑 球 模 型
不 Ar 同溶 半剂 径中 光的 滑振 球动 内弛 豫 时 间
Radius / nm 1.2 1.5 1.8 2.2 Shift / cm-1 1.20 1.79 1.76 1.75
受限于不同半径的光滑球内I2在Ar溶剂中的振动光谱位移
单 壁 碳 纳 米 管 模 型
Radius / nm Shift / cm-1
0.68 (10,10)
1.02 (15,15) 1.36 (20,20) 1.70 (25,25) 2.04 (30,30)
3.42
3.53 3.54 3.55 3.55
受限于不同半径的碳纳米管内I2在 Ar溶剂中的振动光谱位移密度为0.5 g/cm3
纳米反应器
自然界生命体系中的化学变化 绝佳的反应环境
R
+
R
R
R R
product shape selectivity
• “三十年前，如果说并非大多数化学家，那末至少 是有许多化学家嘲笑量子化学研究，认为这些工 作对化学用处不大，甚至几乎完全无用。现在的 情况却是完全两样了…。当90年代行将结束之际， 我们看到化学理论和计算研究的巨大进展，导致 整个化学正在经历一场革命性的变化。Kohn和 Pople是其中的两位最优秀代表”
1986：李远哲：“ 在十五年前，如果理论 结果与实验有矛盾，那么经常证明是理论结果错 了。但是最近十年则相反，常常是实验错了。… 量子力学有些结果是实验工作者事先未想到的，
结构生物学大分子模拟课件
现在学习的是第12页，共69页
到2007年3月，已知的蛋白质序列超过60万 条，而测定了三维结构的蛋白质仅为4万多个
蛋白质三维结构的测定已经成为生命科学 发展的“瓶颈”
现在学习的是第13页，共69页
蛋白质结构的基本概念
蛋白质的二级结构单元 – a螺旋 – b折叠股（b-strand）和b折叠片（b-sheet）
Turns Random coil (its neither!)
现在学习的是第18页，共69页
酰胺平面
肽键的键长介于单键和双键之间，具有部分双键的性质，不能自由旋转。
现在学习的是第19页，共69页
a螺旋
α－螺旋结构
最常见、含量最丰富的二级结构
现在学习的是第20页，共69页
现在学习的是第21页，共69页
现在学习的是第35页，共69页
蛋白质分子表面的环区域
现在学习的是第36页，共69页
蛋白质分子表面的环区域（LOOP）
大多数蛋白质分子是由a-螺旋和/或b回折构成的。这些a-螺旋和/
或b回折通常由不同长度和不规则形状的环区域相连接。 螺旋与回折的组合形成了分子的稳定的疏水内核，环区
域位于分子的表面。主链上的这环区域的C＝O和NH基团一
特点：不依赖于已知的结构模式，是一种普适的解 决方法，目前尚处于探索阶段
现在学习的是第42页，共69页
同源和比较模建方法
根据蛋白质结构的相似性，以已知蛋白质的结构为 模板构建未知蛋白质的三维结构
– 在进化过程中蛋白质三维结构的保守性远大 于序列的保守性，当两个蛋白质的序列同源 性/相似性高于35%时，一般情况下它们的三 维结构基本相同
针对局域能量极小，不是整个系统 能够计算含有大量原子的体系 简单有效，目前应用得最广泛
中国矿业大学分子模拟课件第一章
分 子 模 拟牛继南njn0516@2011.3第一章引言•分子模拟(molecular simulation)根据物理和化学的基本原理, 建立一种以计算数据(由计算机来执行)来代替实验测量的研究方法,获取相关的物理和化学信息.•分子模拟的作用①模拟材料的结构②计算材料的性质③预测材料行为④验证试验结果(重现试验过程)⑤从微观角度认识材料�总之,是为了更深层次理解结构,认识各种行为.•分子模拟的优势:可以降低试验成本具有较高的安全性实现通常条件下较难或无法进行的试验超低温(低于-100℃)超高压(大于100Mpa)研究过快或过慢的反应时间“让飞一会儿……”“时间,你妈喊回家吃饭……”从微观角度认识材料•分子模拟的重要性:"总之,我愿意强调我的信念:计算化学的年代已经到来,成千上百的化学家以计算机代替实验室,来获得众多的化学信息.唯一的障碍是你必须偿付机时费."-R.S.Mulliken的诺贝尔获奖感言(1966年)•分子模拟的应用情况1980年以来全球每年关于”Molecular Simulation”的文章(Engineering Village)1980年以来各国关于”Molecular Simulation”的文章(Engineering Village )Engineering Village Search query: (molecular simulation)Source title•Journal Of Chemical Physics 6019•Journal Of Physical Chemistry B3700•Physical Review Letters1544•Physical Review E - Statistical, Nonlinear, And Soft Matter Physics 1487•Journal Of Physical Chemistry A1238•Macromolecules1173•Journal Of The American Chemical Society1108•Journal Of Physics Condensed Matter876•Proceedings Of Spie - The International Society For Optical Engineering813•Materials Research Society Symposium – Proceedings 778•Surface Science756•Molecular Simulation755•Langmuir726•Nuclear Instruments And Methods In Physics Research, Section B: Beam Interactions With Materials And Atoms642•Journal Of Non-Crystalline Solids608•Journal Of Physical Chemistry C604•Journal Of Computational Chemistry570•Journal Of Applied Physics569•Molecular Physics563•Biochemistry505•Computational Materials Science499•Computer Physics Communications488•Fluid Phase Equilibria488•Chemical Physics Letters483•Polymer468•Journal Of Biological Chemistry435•Applied Physics Letters431•Materials Research Society Symposium Proceedings 388•Aiche Annual Meeting, Conference Proceedings 383•Biopolymers376•Applied Surface Science344Engineering Village Search query: (molecular simulation)Author affiliation•Department Of Chemistry, University Of California 223•Sandia National Laboratories179•Univ Of California138•Theoretical Division, Los Alamos National Laboratory 114•Department Of Chemical Engineering, University Of California 109•Lawrence Livermore National Laboratory108•Department Of Chemistry, University Of Cambridge 103•Los Alamos National Laboratory98•Pacific Northwest National Laboratory75•Department Of Chemistry, Stanford University74•Department Of Engineering Mechanics, Tsinghua University73•Department Of Chemical Engineering, Vanderbilt University72•Accelerator Laboratory, University Of Helsinki72•Department Of Chemical Engineering, Princeton University71•Kyoto Univ71•Department Of Chemistry, Northwestern University71•Pennsylvania State Univ70•Department Of Chemistry, New York University69•Aiaa69•Department Of Chemistry And Biochemistry, University Of California66•Cornell Univ65•Department Of Chemical And Biomolecular Engineering, National University Of Singapore64•Centre For Molecular Simulation, Swinburne University Of Technology63•Lawrence Livermore Natl. Laboratory62•Osaka Univ62•Department Of Physics And Astronomy, University College London 62•Department Of Chemistry, Princeton University62•Department Of Chemistry, Columbia University61•Department Of Chemical Engineering, University Of Massachusetts61•Massachusetts Inst Of Technology61•Institute For Materials Research, Tohoku University60•Department Of Chemical Engineering, University Of Michigan60•Materials Science Division, Argonne National Laboratory60•Tohoku Univ59•Department Of Chemistry, Indian Institute Of Technology59•Division Of Physical Chemistry, Arrhenius Laboratory, Stockholm University57•Department Of Chemical Engineering, Massachusetts Institute Of Technology57•Center For Molecular Modeling, Department Of Chemistry, University Of Pennsylvania56•Department Of Chemical Engineering, Tsinghua University56•Division Of Engineering, Brown University55•Department Of Chemistry, University College London55•Institute Of Fluid Science, Tohoku University55•Institute Of Physics, University Of Silesia55•School Of Materials Science And Engineering, Georgia Institute Of Technology54•Research School Of Chemistry, Australian National University54•Ieee53•Department Of Chemistry, University Of Michigan53•Institut Laue-Langevin53•Department Of Materials Science And Engineering, University Of Florida52•Oak Ridge National Laboratory52•Max-Planck-Institut Fur Polymerforschung51•Univ Of Tokyo50•Argonne National Laboratory50•Department Of Chemical Engineering, University Of Tennessee49•Department Of Chemistry, Zhejiang University49•Graduate School, Chinese Academy Of Sciences48•Department Of Mechanical Engineering, University Of Tokyo47•Department Of Biochemistry, University Of Oxford46•Center For Polymer Studies, Department Of Physics, Boston University44•Univ Of Michigan38•Physical Chemistry38725•Computer Applications29119•Atomic and Molecular Physics22202•Mechanics16094•Physical Properties of Gases, Liquids and Solids 15486•Organic Compounds13792•Chemical Reactions11443•Chemical Operations10703•Chemical Products Generally10607•Mathematics9635•Inorganic Compounds8927•Thermodynamics7951•Mathematical Statistics7896•Crystalline Solids7283•Numerical Methods6268•Electricity: Basic Concepts and Phenomena6154•Chemistry5962•Organic Polymers5795•Light/Optics5341•Crystal Lattice5075•Chemical Agents and Basic Industrial Chemicals 4968•Quantum Theory; Quantum Mechanics4435•Probability Theory4262•Semiconductor Devices and Integrated Circuits 4114•Biological Materials and Tissue Engineering 3979•Polymeric Materials3955•Materials Science3933•Strength of Building Materials; Mechanical Properties 3900•Nanotechnology3478•Fluid Flow, General3437•Solid State Physics3074•Semiconducting Materials2846•Biology2656•Optical Devices and Systems2642•Crystal Growth2613•Biochemistry2595•Computer Software, Data Handling and Applications 2433•Classical Physics; Quantum Theory; Relativity2365•High Energy Physics2307•Nonferrous Metals and Alloys excluding Alkali and Alkaline Earth Metals2198•Colloid Chemistry2139•Calculus1988•Heat Transfer1974•Metallography1958•Metallurgy and Metallography1933•Minerals1896•Computer Theory, Includes Formal Logic, Automata Theory, Switching Theory, Programming Theory1518•Electronic Structure of Solids1510•Polymer Products1447•Algebra1437•Control Systems1432•Precious Metals1352•Coating Materials1325•Compound Semiconducting Materials1308•Magnetism: Basic Concepts and Phenomena1295•Heat Treatment Processes1285•Polymerization1253•Data Processing and Image Processing1231•Medicine and Pharmacology1202•Metallurgy1172Engineering Village Search query: (molecular simulation)Classification code“给一个支点,我就能撬起地球.”-阿基米德“Give me an enough powerful computer, I can simulate the whole world.”•分子模拟框架原理基础模型方案软件工具方法目的Next分子模拟的原理量子力学经典牛顿力学统计物理返回分子模拟的模型•分子模拟中最关键的步骤,最能体现人的思想,即实验设计.•这里的模型指的是建模,不仅仅是在电脑上建立可视化图形!•是如何针对要解决的问题建立合适的物理模型,这个模型要结合所用的软件,合情合理,且高效省时.返回分子模拟软件•量子力学Gaussian, Gammes, VASP, Abinit, Castep......•经典力学Lammps, Gromacs, Namd, Dl_Poly, GULP, Tinker, Towhee, Moldy, Materials Explorer, Discovery Studio, Charmm, Amber......返回布朗动力学蒙地卡罗模拟分子动力学量子力学模拟分子力学分子模拟经典力学模拟DFT半经验方法从头算方法分子模拟方法返回•why not 量子力学, but经典力学?•量子力学是利用波函数来研究微观粒子的运动规律的一个物理学分支学科，它主要研究原子、分子、凝聚态物质，以及原子核和基本粒子的结构、性质的基础理论，它与相对论一起构成了现代物理学的理论基础.•正是量子力学的出现,许多现象才得以真正地被解释，新的、无法直觉想象出来的现象被预言.•至今基本上所有物理间基本相互作用都可以在量子力学的框架内描写.•量子力学以电子的非定域化为基础,一切电子的行为以波函数来表示,要得到波函数必须解薛定谔方程,但分子中的电子数目众多,使得这项工作十分困难.•利用从头算方法通常只能计算上百个原子.•为了提高运算效率,人们又提出半经验方法,但即便如此目前能解的体系仍然在上千个电子的级别.•量子力学or经典力学?•因此用量子力学方法模拟大体系,如上千个原子的聚合物,是不切实际的,即便这种模拟能够进行,大多数情况下得到的大部分信息都要被舍弃.•这是因为在模拟大体系的时候,常常是为了得到统计性质如模量,扩散系数等,而这些量主要依赖于原子核的位置,或者说是主要依赖于原子核的平均配置,此时计算量极大的电子运动信息并不需要.•玻恩-奥本海默近似（Born-Oppenheimer approximation)•根据量子力学的波恩-奥本海默近似，分子的能量可以近似看作构成分子的各个原子的空间坐标的函数，简单地讲就是分子的能量随分子构型的变化而变化，而描述这种分子能量和分子结构之间关系的就是分子力场(Force Field)函数。
生物大分子结构模拟实验PPT课件
（一）三级结构 • 定义
整条肽链中全部氨基酸残基的相对空间位置。 即肽链中所有原子或基团在三维空间的排布位置。
•主要化学键：
疏水作用、离子键、氢键和 范德华力等
24
•三级结构的主要特点：
A.含多种二级结构单元 B. 有明显的折叠层次 C. 是紧密的球状或椭球状实体 D.分子表面有一空穴(活性部位) E. 疏水侧链埋藏在分子内部
纤连蛋白分子的结构域
27
• 结构域(domain)在空间上相对独立
结构域
免疫球蛋白结构域
轻链 重链
28
蛋白质的四级结构 (quaternary structure)
29
（一）具有四级结构形式的蛋白质由多亚基 组成
• 亚基 (subunit)：有些蛋白质分子含有二条或多条 多肽链，每一条多肽链都有完整的三级结构，此多肽链 就是蛋白质分子的亚基 。一般来说亚基不具有生物活 性，只有当这些亚基聚合成一个完整的蛋白质分子后， 才具有生物活性。
（一）肽单元
参与肽键的6个原子C1、C、O、N、H、C2位于同一平面， C1和C2在平面上所处的位置为反式(trans)构型，此同一 平面上的6个原子构成了所谓的肽单元 (peptide unit) 。 12
（二）蛋白质二级结构的主要形式
1. -螺旋 ( -helix ) 2. -折叠 ( -pleated sheet ) 3. -转角 ( -turn ) 4. 无规卷曲 ( random coil )
单糖 五碳糖： 核 糖 脱氧核糖 是核酸的组成成分
寡糖：含3 ~ 十几个糖分子，是构成细胞膜的成分
多糖
糖原 淀粉
3
脂类
• 脂类包括：脂肪酸、类固醇、磷 脂、糖脂等。 • 特点：难溶于水，而易溶于有机溶剂。
分子模拟PPT—第五章 分子动力学模拟运用
1 −E Pj E = ∑ E 2 e j ∑ j Q J j
2 j
k BT
kB ∂ −E =− E je j ∑ Q ∂ (1 T ) j
k BT
=−
kB ∂ ∂E ∂ ln Q ( EQ) = − k B − kB E Q ∂ (1 T ) ∂ (1 T ) ∂ (1 T )
2
∂E = k BT + E2 ∂T
ˆ′ ˆ ˆ ˆ C A (ν ) = A∗ (ν ) A(ν ) = A(ν ) ′ C A (τ ) = 1 2τ run
2τ run −1
2
(v = 0,1,L 2τ run − 1)
∑ ν
=0
ˆ A(ν ) exp(i 2πντ / 2τ run )
2
自相关函数的计算
傅里叶变换计算相关函数的步骤：
第五章
分子动力学模拟 计算的应用
本章内容
运动轨迹分析 热力学特性的计算 径向分布函数 相关函数的计算
运动轨迹分析
结构图像 （可视化图形软件) 几何参量的时间关系曲线 （grace，origin，excel）
键长： rab
= ( xa − xb )2 + ( ya − yb )2 + ( za − zb )2
∞
时间相关函数
物理意义：物理量随时间改变后与其起始的相关性 自相关函数
C A (t ) = A(t ) ⋅ A(0) = A(T + t ) ⋅ A(T ) CB (t ) = B(t ) ⋅ B(0) = B (T + t ) ⋅ B (T )
A ( t ) ⋅ A (0) C A (t ) C A (t ) = = C A (0) A (0) ⋅ A (0)
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
2021/3/7
CHENLI
18
四、结论
1.温度越高扩散越快。 随着温度的升高， 水在聚乙烯醇中和尼龙6中的扩散系数都 增加。
素 6. 课题提出
2021/3/7
CHENLI
3
1.聚合物对水的阻隔性能表征
聚合物与水的相互作用的宏观表现：
1、阻水聚合物 没有明确的界限 2、吸水聚合物 对水的阻隔性或者吸附性的强弱主要是通过水的 透过量，水的吸附量或者水的扩散程度来衡量。
水的阻隔性---扩散系数D（水在聚合物中的扩散）
本质：水在聚合物中的运动（吸附，扩散）实际 上是一种传质现象，它服从fick第二定律。
293.15K
303.15K
323.15K
图2.7 不同温度下水分子在尼龙6中均方位移随时间的变化
2021/3/7
CHENLI
15
三、结果与讨论
1、实验结果与分析 2、活化能的计算
298K时： 水在聚乙烯醇中：D=a/6=2.61 (10-8cm2/s) 水在尼龙6中： D=a/6=0.08 (10-8cm2/s) 323K时： 水在聚乙烯醇中：D=a/6=2.92 (10-8cm2/s) 水在尼龙6中： D=a/6=0.18 (10-8cm2/s) 373K时： 水在聚乙烯醇中：D=a/6=3.11 (10-8cm2/s) 水在尼龙6中： D=a/6=0.66 (10-8cm2/s)
2021/3/7
CHENLI
8
二、实验部分
1、实验原理(扩散定律) 2、实验步骤
2021/3/7
图2.1实验模拟工艺图
CHENLI
9
1、实验原理
扩散定律（fick第二定律）：
单位体积内粒子浓度随时间的变化：
考虑扩散系数受浓度的影响：
扩散系数与均方位移的关系：
扩散系数与均方位移的微分有关：
模拟得到MSD对时间的图，a代表MSD曲线的斜率，所以
小分子物质结构的影响
W．H．Jiang和姜珍考察了小分子物质的尺寸对扩散系数的影响，发现扩散物 质分子尺寸越大，扩散阻力也就越大，相应地扩散系数就越小。
聚合物结构和物理形态的影响
聚合物的结晶度越大，小分子物质在其中的扩散系数就会越小。聚合物的结构 形态对小分子物质在其中的扩散也有很大的影响，如在极性聚合物中极性分子 的扩散要比非极性分子的扩散容易，而在非极性聚合物中极性分子的扩散要比 非极性分子的扩散难。
2021/3/7
CHENLI
10
2、实验步骤
1.建立初始结构 2.建立无定形晶胞 3.进行晶胞的驰豫 4.分子动力学的运行和分析 5.输出数据并计算扩散系数D
组成 调节 运行 计算
2021/3/7
CHENLI
11
建立初始结构
a
b
c
图2.2 建立的初始结构，a-水分子，b-尼龙6，c-聚乙烯醇
优化结构
a
b
图2.3优化的结构图，a-尼龙6，b-聚乙烯醇
2021/3/7
CHENLI
12
建立无定形晶胞
a
b
图2.4聚合物的晶胞图，a-尼龙6的晶胞，b-聚乙烯醇晶胞
晶胞的驰豫
当一个无规则晶胞生成时，分子可能不是等价地分布在晶胞中，这样就造成了 真空区。为了矫正这个，要进行能量最小化来优化晶胞。最小化过后，要进行
2021/3/7
CHENLI
4
2. 实验法测定的介绍（物理变化量）
本体平衡法
测定聚合物重量与时间的曲线，利用数学分析方法计算扩散系数
脉冲梯度场核磁共振法
某原子在两个脉冲时间间隔内在磁场中的回声振幅衰减而获得聚合 物体系中溶剂自扩散系数
激光全息技术
利用光的干涉原理将物体光波完全精确记录，并能使之再现还原的 二次成像技术
2021/3/7
CHENLI
5
3. 分子动力学模拟法的介绍
量子化学法 分子动力学法 分子力学法
4.分子模拟中扩散系数的测定方法
Einstein法 微分一区限变分法 聚类分析法
2021/3/7
CHENLI
6
5.影响小分子物质（水）在聚合物 中扩散系数的因素
温度的影响
W．H．Jiang和姜珍等报道了温度对扩散系数的影响，发现温度越高扩散系数 越 大，这是由于温度的升高加剧了分子运动。
表中数据我们通过图形来对比分析：
2021/3/7
CHENLI
17
3.2活化能的计算
扩散系数与温度的关系符合Arrhenius方 程：
D=D0exp（-Ed/RT）或lnD=lnD0- Ed/RT 以lnD对1/T作图，由斜率可求出活化能Ed。
a
b
图3.2聚合物中lnD与T的关系，a-聚乙烯醇，b-尼龙6
高分子树脂 对水的阻隔性能研究
2021/3/7
余勇 指导老师:王孟
化学化工学院 高分子材料与工程专业071班
CHENLI
1
大纲
一.前言 二.实1/3/7
CHENLI
2
一. 前言
1.聚合物对水的阻隔性能表征 2.实验法测定的介绍 3. 分子动力学模拟法的介绍 4.分子模拟中扩散系数的测定方法 5.影响小分子物质在聚合物中扩散系数的因
2021/3/7
CHENLI
7
6. 课题提出
实验法去测定小分子物质在聚合物体系中的扩散系 数有很多困难跟不足。
1.耗时长 2.误差大3.操作条件苛刻4.设备费用高 分子模拟法提供了材料的结构和动力学的细节
1.时间短2.精确度高3. 操作方便4.适用范围广
本实验就通过分子动力学（MD）模拟水分子在聚 乙烯醇和在聚酰胺（尼龙6）中的扩散行为，根据软 件的动力学分析来得到水分子在聚合物中扩散的均方 位移MSD，由得到的均方位移与时间的关系作图， 来得到水分子的扩散系数D。
分子动力学模拟来平衡晶胞。
分子动力学的运行和分析
分析水分子在晶胞中的移动情形，得到水分子在均方位移MSD.
2021/3/7
CHENLI
13
输出数据并计算扩散系数
293.15K
303.15K
323.15K
图2.6不同温度下水分子在聚乙烯醇中均方位移随时间的变化
2021/3/7
CHENLI
14
2021/3/7
CHENLI
16
3.1实验结果与分析
水在不同温度下在两种聚合物中的扩散系数：
温度T
293.15K 303.15K 323.15K
水在聚乙烯醇中的扩散系数 D(10-8cm2/s) 2.61 2.92 3.11
水在尼龙6中的扩散系数 D(10-8cm2/s) 0.08 0.18 0.66