【LAMMPS翻译系列】minimize命令

minimize命令通过不断迭代调整原子坐标的方式对体系的能量进行最小化。

使用语法

minimize etol ftol maxiter maxeval
  • etol = 能量的停止容差(无单位)
  • ftol = 力的停止容差(力的单位)
  • maxiter = 能量最小化器minimizer的最大迭代次数
  • maxeval = 计算力或能量的最大次数

使用举例

minimize 1.0e-4 1.0e-6 100 1000
minimize 0.0 1.0e-8 1000 100000

使用介绍

该命令通过不断迭代调整原子坐标的方式对体系的能量进行最小化。当满足任何一个停止判据(能量或力)时,迭代过程就会终止。迭代结束时刻的构型在很大程度上说就是局域势能最小的构型。更精确地说,这个时候的构型接近目标函数(下面会介绍)的临界点,它有可能是,也有可能不是局域最小值。

优化过程所使用的优化算法是由min_style命令进行设置的。其他的选项可以使用min_modify命令进行设置。minimize命令可以与run命令交替使用,从而可以交替着进行弛豫和动力学。能量最小化的过程限制了原子在一个迭代步中移动的距离,使用动力学你就可以将高度重叠的原子(具有很大的能量和力)相互推开。

还有一种对体系进行弛豫的方法是使用较小的或受限的时间步,对体系运行动力学。或者在运行动力学的时候使用fix viscous命令施加阻尼力,慢慢地拖动体系的总动能。在运行动力学的时候,可以使用 pair_style soft 势函数分开重叠的原子。

能量最小化算法cg、sd、hftn先在外层迭代中确定坐标改变的搜索方向。然后基于线搜索(line search)算法进行内层迭代。线搜索算法通过不断的计算力和能量来确定新的坐标。目前的程序使用的是一种回溯算法,就计算力的次数来说可能不是最优的算法,但经过测试确是最稳定的算法。Nocedal and Wright's Numerical Optimization 对回溯算法进行了介绍(Procedure 3.1 on p 41)。

能量最小化算法quickmin和fire是使用欧拉积分步运行带阻尼的分子动力学。因此,需要定义时间步。虽然使用更大的时间步执行效率会更高一些,但一般使用与running dynamics相同的时间步。

能量最小化的过程实际上是对下面的目标函数进行最小化的过程。该目标函数是体系的总势能,它是体系中N个原子坐标的函数。



其中第一项是所有无键作用的对势相互作用的总和,也包括长程库伦相互作用。第二项到第五项分别代表键(bond)、角(angle)、二面角(dihedral)和improper的相互作用。最后一项考虑了某些fix命令通过约束或在原子上施加力而带来的作用,比如通过与壁面进行作用。参考fix命令,了解它是如何影响能量最小化的。

能量最小化的起始点是开始进行能量最小化之时的原子构型。




如果下面任何一个判据满足了,能量最小化过程就会停止。

  • 外层迭代的能量改变值小于etol
  • 全局力矢量的长度小于ftol
  • the 2-norm (length) of the global force vector is less than the ftol
  • 由于步长返回0.0而导致线性搜索失败
  • 外层迭代次数或时间步超过了maxiter
  • 计算力的次数超过了maxeval

对于第一个判据,指定的能量容差etol是没有单位的。当两个相邻的迭代步的能量差除以总的能量值小于或等于该容差时,这个判据就满足了。举个例子,如果将etol设置为1.0e-4,就是说能量容差等于总能的1/10000。使用带阻尼的动力学进行能量最小化的时候,这个判据只有在速度被重置为0之后的一些步才会被检查,否则会使能量最小化过程较早的收敛。

对于第二个判据,指定的力容差ftol是力的单位,因为它是所有原子的全局力矢量的长度。对于N个原子而言,全局力矢量是一个维度为3N的矢量。由于在能量最小化之后,力矢量的很多元素都接近0,你可以认为是ftol是任何原子任何方向上受力的上边界(译注:上边界即最高要求,所有要求都达到才可以;下边界即最低要求,只要达到一个就行)举个例子,如果将ftol设置为1.0e-4,就是说通过能量最小化之后,没有任何原子的在x,y,z方向上受力大约1.0e-4(力的单位)。

etol和ftol都可以同时或单独设置为0.0。在这种情况下,其他的判据会用来结束能量最小化过程。

在能量最小化的过程中,最外层的迭代次数被作为时间步处理。实际的输出也是以它为参考的,比如热力学输出或dump输出或重启动文件输出。

使用thermo_style custom命令中的关键字fmax或fnorm,可以有效的监控能量最小化的过程。Note that these outputs will be calculated only from forces on the atoms, and will not include any extra degrees of freedom, such as from the fix box/relax command.

在能量最小化结束后,程序会打印一段统计摘要信息,介绍满足了何种收敛判据,以及能量、受力、最终的线性搜索和跌打次数等。下面是一个例子:

Minimization stats:
  Stopping criterion = max iterations
  Energy initial, next-to-last, final = 
       -0.626828169302     -2.82642039062     -2.82643549739
  Force two-norm initial, final = 2052.1 91.9642
  Force max component initial, final = 346.048 9.78056
  Final line search alpha, max atom move = 2.23899e-06 2.18986e-05
  Iterations, force evaluations = 2000 12724
  • 其中的三个能量值分别是进行能量最小化之前的能量、能量最小化结束之前一步的能量、能量最小化结束时的能量。这就是判据参数etol所要检查比较的。
  • 其中的“Force max component”是指能量最小化之前和之后的全局力矢量的长度。这是判据参数ftol说要检查比较的。
  • 其中的“Force max component”是指全局力矢量中最大元素的绝对值。
  • 其中的“Final line search alpha”,即线性搜索参数alpha,乘以最后迭代步时的最大力元素(即力矢量中的最大值),就是最后迭代步中原子最大移动距离(max atom move)。如果线搜索不能使能量降低,alpha就会是0.0。即使alpha不是0.0,但如果原子最大移动距离相比于原子坐标太小,就会导致最后一步原子几乎没有移动,进而能量改变很微小,而使得能量最小化过程结束。换句话说,即使受力不为0,但如果这些力太小了,就可能使原子发生移动的距离小于机器精度,从而就不能再使能量更小了。

判据参数maxitermaxeval会分别检查比较迭代次数和力的计算次数。




注意:在LAMMPS中有些力场(势函数)会因为不连续或其他近似,使得在进行能量最小化的时候不能使用精度较高的容差。举个例子,在进行能量最小化的是会,你需要使用那些在截断距离过零点的势类型,即使你后面又切换成运行分子动力学。如果你不这样做,但任何一对原子之间的距离从+epsilon截断变成-epsilon时,系统的总能会出现不连续,从而导致能量最小化过程表现不稳定[原文: the total energy of the system will have discontinuities when the relative distance between any pair of atoms changes from cutoff+epsilon to cutoff-epsilon and the minimizer may behave poorly]。 Some of the manybody potentials use splines and other internal cutoffs that inherently have this problem. The long-range Coulombic styles (PPPM, Ewald) are approximate to within the user-specified tolerance, which means their energy and forces may not agree to a higher precision than the Kspace-specified tolerance. In all these cases, the minimizer may give up and stop before finding a minimum to the specified energy or force tolerance.

Note that a cutoff Lennard-Jones potential (and others) can be shifted so that its energy is 0.0 at the cutoff via the pair_modify command. See the doc pages for inidividual pair styles for details. Note that Coulombic potentials always have a cutoff, unless versions with a long-range component are used (e.g. pair_style lj/cut/coul/long). The CHARMM potentials go to 0.0 at the cutoff (e.g.pair_style lj/charmm/coul/charmm), as do the GROMACS potentials (e.g. pair_style lj/gromacs).

If a soft potential (pair_style soft) is used the Astop value is used for the prefactor (no time dependence).

The fix box/relax command can be used to apply an external pressure to the simulation box and allow it to shrink/expand during the minimization.

只有一小部分fix命令(一般是用来施加力约束的)会在能量最小化的过程中被激活。具体类型的fix命令的页面会告诉你它是否会属于这类命令。

注意:有些fix命令会与能量最小化过程中被激活,计算一部分势能。如果要把这部分能量包括到体系的总势能中,必须先把fix_modify命令中的energy选项打开。细节可以参考具体类型的fix页面。

使用限制

这里列出了还没有在程序中实现的功能,看看你是否知道如何实现它们:

在能量最小化的时候使用fix shake命令会出现错误,因为该命令关闭键的作用,而这些作用却被包括在了体系的势能中。 The effect of a fix shake can be approximated during a minimization by using stiff spring constants for the bonds and/or angles that would normally be constrained by the SHAKE algorithm.

能量最小化的过程也不支持Fix rigid命令。虽然定义了这个命令并不会出现错误,但能量最小化的时候却不能使定义为刚体的部分保持为刚体。需要注意,如果刚体内部的键、角被关闭了(比如通过命令 neigh_modify exclude),它们就不会再对总势能有贡献了,这可能跟想象的不一样。

Pair potentials that produce torque on a particle (e.g. granular potentials or the GayBerne potential for ellipsoidal particles) are not relaxed by a minimization. More specifically, radial relaxations are induced, but no rotations are induced by a minimization, so such a system will not fully relax.

相关命令

min_modifymin_stylerun_style

标签: lammps, lammps翻译

相关文章推荐

添加新评论 (无需注册,可直接评论)

已有 9 条评论

  1. xxx007980632547

    minimize和用其他的nvt nve跑动力学有什么区别呢??

    1. 准确说,minimize是分子静力学,而后者是分子动力学。

  2. 辉辉豆

    请问下如何根据自己建立的工件选择合适的参数呢?有时候会遇到运行程序后提醒我WARNING: Resetting reneighboring criteria during minimization (../min.cpp:168) 参数的选择会不会影响我弛豫后的结果?

    1. 这个命令的参数更多是影响精度,设置得小一些就好。这里一般不会出问题。

      1. 辉辉豆

        恩恩 好的,谢谢博主,希望博主有时间能贴一些有关系综和弛豫方面的总结~

  3. 请问minimize 1.0e-4 1.0e-6 100 1000
    这里的 100 和1000是什么

    1. admin

      文中已经解释了,分别是能量和力的最大迭代步数。

      1. 走走停停边边

        请问下老师请问minimize 1.0e-4 1.0e-6 100 1000中前两个参数的指数-4和-6一般是如何确定的,按手册的理解是不是指能量和力的偏差范围分别是小于0.0001和0.000001?我的理解没错吧?

        1. 设置的足够小就可以了。一般1e-6就已经很小了。