分子模拟一般性步骤.doc

Below is presented a generalised procedure for performing a simulation. The
exact steps and processes involved will vary depending on exactly what is being
attempted. Use as a general guide only!

1> Clearly identify the property / phenomena of interest to be studied by performing the simulation.

2>Select the appropriate tools to be able to perform the simulation and observe the property / phenomena of interest. It is important to read and familiarise
yourself with publications by other researchers on similar systems. Tools include:
- software to perform the simulation with, consideration of force field may
influence this decision.
- force field which describes how the atoms / particles within the system interact with each other. Select one that is appropriate for the system being studied and the property / phenomena of interest. Very important and non-trivial step!

3>Obtain / generate the initial coordinate file for each molecule to be placed
within the system.

4>Generate the raw starting structure for the system by placing the molecules
within the coordinate file as appropriate. Molecules may be specifically placed or arranged randomly.

5>Obtain / generate the topology file for the system, using (for example)
pdb2gmx, PRODRG or your favourite text editor in concert with chapter 5 of the GROMACS Manual.

6>Describe a simulation box (e.g. using editconf) whose size is appropriate for
the eventual density you would like, fill it with solvent (e.g. using genbox), and
add any counter-ions needed to neutralize the system (e.g. using grompp and
genion). In these steps you may need to edit your topology file to stay current
with your coordinate file.

7>Run an energy minimisation simulation on the system (using grompp and
mdrun). This is required to sort out any bad starting structures caused during
generation of the system, which may cause the production simulation to crash.

8>Select the appropriate simulation parameters for the equilibration simulation (defined in .mdp file). You need to be consistent with how force field was
derived. You may need to simulate at NVT with position restraints on your
solvent and/or solute to get the T almost right, then relax to NPT to fix the
density, then move further (if needed) to reach your production simulation
ensemble (e.g. NVT, NVE).

9>Run the equilibration simulation for sufficient time so that the system relaxes sufficiently to allow the production run to be commenced (using grompp and
mdrun, then g_energy and trajectory visualisation tools).

10>Select the appropriate simulation parameters for the production simulation (defined in .mdp file) in particular be careful not to re-generate the velocities.
You still need to be consistent with how the force field was derived and how to
measure the property / phenomena of interest.

11>Run the production simulation for sufficient time so that property / pheno-
mena of interest can be observed in required detail (using grompp/tpbconv and mdrun).

12>Analyse / visualise the resulting trajectory and data files to obtain information on the property / phenomena of interest.

分子模拟一般性步骤 [转载]
以下是做模拟的一般性步骤，具体的步骤和过程依赖于确定的系统或者是软件，但这不影响我们把它当成一个入门指南：

1）首先我们需要对我们所要模拟的系统做一个简单的评估, 三个问题是我们必须要明确的：
做什么（what to do）为什么做（why to do）怎么做（how to do）

2）选择合适的模拟工具，大前提是它能够实现你所感兴趣的目标，这需要你非常谨慎的查阅文献，看看别人用这个工具都做了些什么，有没有和你相关的，千万不要做到一半才发现原来这个工具根本就不能实现你所感兴趣的idea，切记！

考虑1：软件的选择，这通常和软件主流使用的力场有关，而软件本身就具体一定的偏向性，比如说，做蛋白体系，Gromacs，Amber，Namd均可；做DNA, RNA体系，首选肯定是Amber；做界面体系，Dl_POLY比较强大，另外做材料体系，Lammps会是一个不错的选择

考虑2：力场的选择。力场是来描述体系中最小单元间的相互作用的，是用量化等方法计算拟合后生成的经验式，有人会嫌它粗糙，但是它确确实实给我们模拟大系统提供了可能，只能说关注的切入点不同罢了。常见的有三类力场：全原子力场，联合力场，粗粒化力场；当然还有所谓第一代，第二代，第三代力场的说法，这里就不一一列举了。

  再次提醒注意：必须选择适合于我们所关注体系和我们所感兴趣的性质及现象的力场。

3）通过实验数据或者是某些工具得到体系内的每一个分子的初始结构坐标文件，之后，我们需要按我们的想法把这些分子按照一定的规则或是随机的排列在一起，从而得到整个系统的初始结构，这也是我们模拟的输入文件。

4）结构输入文件得到了，我们还需要力场参数输入文件，也就是针对我们系统的力场文件，这通常由所选用的力场决定，比如键参数和非键参数等势能函数的输入参数。

5）体系的大小通常由你所选用的box大小决定，我们必须对可行性与合理性做出评估，从而确定体系的大小，这依赖于具体的体系，这里不细说了。

6）由于初始构象可能会存在两个原子挨的太近的情况（称之为bad contact），所以需要在正式模拟开始的第一步进行体系能量最小化，比较常用的能量最小化有两种，最速下降法和共轭梯度法，最速下降法是快速移除体系内应力的好方法，但是接近能量极小点时收敛比较慢，而共轭梯度法在能量极小点附近收敛相对效率高一些，所有我们一般做能量最小化都是在最速下降法优化完之后再用共轭梯度法优化，这样做能有效的保证后续模拟的进行。

7）以平衡态模拟为例，你需要设置适当的模拟参数，并且保证这些参数设置和力场的产生相一致，举个简单的例子，gromos力场是用的范德华势双截断来定范德华参数的，若你也用gromos力场的话也应该用双截断来处理范德华相互作用。常见的模拟思路是，先在NVT下约束住你的溶质（剂）做限制性模拟，这是一个升温的过程，当温度达到你的设定后, 接着做NPT模拟，此过程将调整体系的压强进而使体系密度收敛。
经过一段时间的平衡模拟，在确定系统弛豫已经完全消除之后，就可以开始取数据了。如何判断体系达到平衡，这个问题是比较技术性的问题，简单的讲可以通过以下几种方式，一，看能量（势能，动能和总能）是否收敛；二，看系统的压强，密度等等是否收敛；三看系统的RMSD是否达到你能接受的范围，等等。

8）运行足够长时间的模拟以确定我们所感兴趣的现象或是性质能够被观测到，并且务必确保此现象出现的可重复性。

9）数据拿到手后，很容易通过一些可视化软件得到轨迹动画，但这并不能拿来发文章。真正的工作才刚刚开始——分析数据，你所感兴趣的现象或性质只是表面，隐含在它们之中的机理才是文章中的主题。

参考文献: 1. 陈正隆《分子模拟的理论与实践》讲习班教材
                     2. Steps to Perform a Simulation (GMX) http://www.mdbbs.org/thread-54-1-1.html

附录：分子模拟的简要介绍

Molecular modelling is a collective term that refers to theoretical methods and computational techniques to model or mimic the behaviour of molecules. The techniques are used in the fields of computational chemistry, computational biology and materials science for studying molecular systems ranging from small chemical systems to large biological molecules and material assemblies. The simplest calculations can be performed by hand, but inevitably computers are required to perform molecular modelling of any reasonably sized system. The common feature of molecular modelling techniques is the atomistic level description of the molecular systems; the lowest level of information is individual atoms (or a small group of atoms). This is in contrast to quantum chemistry (also known as electronic structure calculations) where electrons are considered explicitly. The benefit of molecular modelling is that it reduces the complexity of the system, allowing many more particles (atoms) to be considered during simulations.

Molecular mechanics is one aspect of molecular modelling, as it is refers to the use of classical mechanics/Newtonian mechanics to describe the physical basis behind the models. Molecular models typically describe atoms (nucleus and electrons collectively) as point charges with an associated mass. The interactions between neighbouring atoms are described by spring-like interactions (representing chemical bonds) and van der Waals forces. The Lennard-Jones potential is commonly used to describe van der Waals forces. The electrostatic interactions are computed based on Coulomb's law. Atoms are assigned coordinates in Cartesian space or in internal coordinates, and can also be assigned velocities in dynamical simulations. The atomic velocities are related to the temperature of the system, a macroscopic quantity. The collective mathematical expression is known as a potential function and is related to the system internal energy (U), a thermodynamic quantity equal to the sum of potential and kinetic energies. Methods which minimize the potential energy are known as energy minimization techniques (e.g., steepest descent and conjugate gradient), while methods that model the behaviour of the system with propagation of time are known as molecular dynamics.

E = E_bonds + E_angle + E_dihedral + E_non − bonded

E_non − bonded = E_electrostatic + E_vanderWaals

This function, referred to as a potential function, computes the molecular potential energy as a sum of energy terms that describe the deviation of bond lengths, bond angles and torsion angles away from equilibrium values, plus terms for non-bonded pairs of atoms describing van der Waals and electrostatic interactions. The set of parameters consisting of equilibrium bond lengths, bond angles, partial charge values, force constants and van der Waals parameters are collectively known as a force field. Different implementations of molecular mechanics use slightly different mathematical expressions, and therefore, different constants for the potential function. The common force fields in use today have been developed by using high level quantum calculations and/or fitting to experimental data. The technique known as energy minimization is used to find positions of zero gradient for all atoms, in other words, a local energy minimum. Lower energy states are more stable and are commonly investigated because of their role in chemical and biological processes. A molecular dynamics simulation, on the other hand, computes the behaviour of a system as a function of time. It involves solving Newton's laws of motion, principally the second law, F = ma. Integration of Newton's laws of motion, using different integration algorithms, leads to atomic trajectories in space and time. The force on an atom is defined as the negative gradient of the potential energy function. The energy minimization technique is useful for obtaining a static picture for comparing between states of similar systems, while molecular dynamics provides information about the dynamic processes with the intrinsic inclusion of temperature effects.

Molecules can be modelled either in vacuum or in the presence of a solvent such as water. Simulations of systems in vacuum are referred to as gas-phase simulations, while those that include the presence of solvent molecules are referred to as explicit solvent simulations. In another type of simulation, the effect of solvent is estimated using an empirical mathematical expression; these are known as implicit solvation simulations.

Molecular modelling methods are now routinely used to investigate the structure, dynamics and thermodynamics of inorganic, biological, and polymeric systems. The types of biological activity that have been investigated using molecular modelling include protein folding, enzyme catalysis, protein stability, conformational changes associated with biomolecular function, and molecular recognition of proteins, DNA, and membrane complexes.