LAMMPS中的系综（NPT/NVT）命令

语法：

fix ID group-ID style_name keyword value ...

• ID 是fix名，group-ID是fix所作用的群组名。这两个是所有fix命令都需要的。

• style_name可以选择NVT、NPT或者NPH

• 关键词keyword和对应的数值value可以添加

keyword = temp or iso or aniso or tri or x or y or z or xy or yz or xz or couple
          or tchain or pchain or mtk or tloop or ploop or nreset or drag or dilate 
          or scalexy or scaleyz or scalexz or flip or fixedpoint or update
  temp values = Tstart Tstop Tdamp
    Tstart,Tstop = external temperature at start/end of run
    Tdamp = temperature damping parameter (time units)
  iso or aniso or tri values = Pstart Pstop Pdamp
    Pstart,Pstop = scalar external pressure at start/end of run (pressure units)
    Pdamp = pressure damping parameter (time units)
  x or y or z or xy or yz or xz values = Pstart Pstop Pdamp
    Pstart,Pstop = external stress tensor component at start/end of run (pressure units)
    Pdamp = stress damping parameter (time units)
  couple = none or xyz or xy or yz or xz
  tchain value = N
    N = length of thermostat chain (1 = single thermostat)
  pchain values = N
    N length of thermostat chain on barostat (0 = no thermostat)
  mtk value = yes or no = add in MTK adjustment term or not
  tloop value = M
    M = number of sub-cycles to perform on thermostat
  ploop value = M
    M = number of sub-cycles to perform on barostat thermostat
  nreset value = reset reference cell every this many timesteps
  drag value = Df
    Df = drag factor added to barostat/thermostat (0.0 = no drag)
  dilate value = dilate-group-ID
    dilate-group-ID = only dilate atoms in this group due to barostat volume changes
  scalexy value = yes or no = scale xy with ly
  scaleyz value = yes or no = scale yz with lz
  scalexz value = yes or no = scale xz with lz
  flip value = yes or no = allow or disallow box flips when it becomes highly skewed
  fixedpoint values = x y z
    x,y,z = perform barostat dilation/contraction around this point (distance units)
  update value = dipole or dipole/dlm
    dipole = update dipole orientation (only for sphere variants)
    dipole/dlm = use DLM integrator to update dipole orientation (only for sphere variants)

例子：

fix 1 all nvt temp 300.0 300.0 100.0
fix 1 water npt temp 300.0 300.0 100.0 iso 0.0 0.0 1000.0
fix 2 jello npt temp 300.0 300.0 100.0 tri 5.0 5.0 1000.0
fix 2 ice nph x 1.0 1.0 0.5 y 2.0 2.0 0.5 z 3.0 3.0 0.5 yz 0.1 0.1 0.5 xz 0.2 0.2 0.5 xy 0.3 0.3 0.5 nreset 1000

详细说明：

这些命令用来对Nose-Hoover型的non-Hamiltonian运动方程进行时间积分。该运动方程用以从正则系综、等温等压系综、等焓系综采样生成原子的位置和速度，从而对所定义的群组内的原子在每一个时间步进行位置和速度的更新。（添注：事实上，系综最简单的理解是通过命令的字母，譬如NVT中N代表系统的原子数、V代表系统的体积、T代表系统的温度，因此NVT就代表恒原子数恒温恒体积的系综，具体的积分迭代形式可以查阅相关的分子动力学模拟专著。）

控温器（thermostatting）和控压器(barostatting)通过添加一些动态变量来实现，这些动态变量与粒子速度（控温）和模拟区域的尺寸（控压）相关联。在基本单向的控温和控压之外，这些命令还可以与粒子控温器耦合创建控温器链，或者与控压变量耦合创建另一条控温器链。控压器与盒子的总体积、或者包括xy、xz及yz的盒子夹角相互耦合。来自于控压器的外部压强可以定义为一个标量（等压系综）也可以定义为对称应力张量（固定应力系综）。如果命令被正确使用，时间平均的温度和应力张量会与通过“Tstart/Tstop”和“Pstart/Pstop”所指定的目标值相符。（添注：简单的理解来说，控温改变粒子速度；控压改变系统盒子的尺寸。注意比较控压器定义为标量或者是对称应力张量的区别，在后者盒子会有剪切自由度，而前者情况下，盒子在控压过程中一直保持只是长方体的block。）

所使用的运动方程是Shinoda等人的方程，将Martyna、Tobias和Klein的流体静力学方程与Parrinello中的Parrinello和Rahman提出的应变能结合在一起。 时间积分框架采用Tuckerman等人推导的时间可逆测度保持 （time-reversible measure-preserving）Verlet和rRESPA积分器。

——————————————————————————————————————

fix命令下通过nvt和npt的控温参数是用temp关键词来进行设定。其它控温相关的关键词包括tchain、tloop和drag，下面对此进行讨论。

控温器仅仅应用在拥有平移自由度的粒子。平移自由度可以在进行控温前移除偏速度（bias velocity）,参见下面进行的讨论。期望的温度在运行的每一个时间步都随着Tstart到Tstop不断上升。Tdamp参数的单位是时间单位，用来决定温度的释放速度的快慢。举例来说，Tdamp的数值如果设定为10.0，这意味着释放温度在大概10个时间步的范围内进行（譬如tau，fa及ps，这根据units命令中的设定时间单位来确定）。只有在fix命令中所指定的group，他们的速度和位置才会通过指定系综下进行积分更新。

NOTE
Nose-Hoover 控温器不会在Tdamp取任意值得情况下都能很好地工作。如果Tdamp取值太小了，温度会剧烈的抖动；如果取值太大了，温度需要很长的时间才能取得平衡。对于大多数模型，一个很好地选择是Tdamp取值大约为100个时间步。值得注意的是，对于绝大部分时间设定，这不等同于100个时间单位。一个简单的办法来确定这项，是通过直接变量定义的方法来确定，例子如下：（添注：这很重要，注意区别时间单位和时间步的区别，这两者只有dt=1的时候才是一致的，举例来说，如果在metal单位制下，时间单位是ps，Tdamp的单位也是ps。此时如果时间步设定为0.001，也就是1fs，这个时候，Tdamp取100个时间步，那就是0.1ps，那么Tdamp=0.1，而不是100。）

fix 1 all nvt temp 300.0 300.0 $(100.0*dt)

——————————————————————————————————————

在该fix命令下npt和nph的控制器参数可以使用一个或一个以上，分别为iso, aniso, tri, x, y, z, xy, xz, yz 和 couple 等关键词。 这些关键词可以用来特定控制外部应力张量的六个分量，或者耦合这些分量一起因此他们（应力分量）所代表的尺寸可以在一个恒压模拟中一起变化。

其它相关的关键词pchain, mtk, ploop, nreset, drag, 和dilate 会在接下来进行讨论。

正交模拟盒子（长方体）有三个可以调整的尺寸（x, y, z）。三斜（非正交）模拟盒子可以有六个可以调整的尺寸 （x, y, z, xy, xz, yz)。通过命令creat_box, read data 和 read_restart 命令可以指定模拟盒子是正交还是非正交的，并且明白xy，xz，yz等倾角参数的内涵。（添注：相关内容可以参见我的知乎专栏 泛柏舟：lammps中怎么建立斜方晶并进行计算）

六种应力张量分量下的目标压强可以独立的通过x, y, z, xy, xz, yz 关键词来控制，分别对应于模拟盒子的六个尺寸。对于每一种分量，外部压强或张量分量在每一个时间步都由Pstart向Pstop演化。如果目标压强是指定为一个分量，那么该分量对应的盒子尺寸也会在模拟中发生变化。例如，如果使用关键词y，那么盒子在y方向的长度就会发生变化。如果xy被使用，那么盒子xy夹角参数也会发生改变。在未指定应力分量的情况下，对应的盒子尺寸就不会发生改变，但是你仍然有机会可以通过另外的fix deform命令来改变尺寸。（添注：fix deform在NVT系综中同样可以改变盒子尺寸，很多人困惑的是，如果是NVT，那么体积不变尺寸不变，为什么fix deform可以改变尺寸呢，可以理解为每一次deform后的体积是不变的，此时在NVT中进行驰豫。就相当于拉伸过程中，每次拉伸后都是一个独立的NVT系综。）

值得注意的是，如果使用xy, xz, 或者yz 关键词，模拟盒子必须是非正交三斜的，尽管初始参数可以设置为0.0 。（添注：简单理解，必须允许盒子的夹角有自由度，如果你定义成正交的长方体，相当于就不允许夹角有自由度，初始参数只是结构的初始构型。）

在所有控压器关键词中，Pdamp参数操作起来和Tdamp基本是一致的，用来决定压强释放的时间尺度。例如，取值10.0就意味着在大约10个时间单位下释放（进行控制）一次压强。

NOTE
Nose-Hoover 控压器不会在Pdamp取任意值得情况下都能很好地工作。如果Pdamp取值太小了，压强和体积会剧烈抖动；如果取值太大了，压强需要很长的时间才能取得平衡。对于大多数模型，一个很好地选择是Tdamp取值大约为1000个时间步。然而，Pdamp的单位是时间单位，对于绝大部分模拟的单位设定中，这不等同于时间步。

不管在fix命令中所指定的group群组包括哪些原子（只有这些原子需要进行时间积分），压强或者应力张量的作用是全局性质的，是作用在整体所有原子上的。类似的，当模拟盒子的尺寸发生改变时，所有的原子会缩放到新的位置。除非，对于dilate-group-ID下的原子使用了关键词dilate。这种方法在某些情况下是十分必要的，例如，对于固体基质中原子进行固定而控制周边流体的压强。而使用该选项需要十分谨慎，由于它会按照不符合物理规律的方式膨胀一些原子而对另一些原子不采取这一作用，这一会使得在一对原理膨胀起源的原子（一个已膨胀，一个未膨胀）之间引入不合理的大的瞬时位移。同样需要说明的是，对于没有在fix nvt中所定义群组的原子，可以采用其他的时间积分，包括fix nve 和fix nvt，不管他们是否受到膨胀作用。

——————————————————————————————————————

关键词couple允许将两个或三个压强张量的对角线分量耦合作用。用关键词指定的值来确定耦合的对象。 例如，xz表示应力张量的Pxx和Pzz分量耦合。 Xyz表示所有3个对角线分量都已耦合。 耦合意味着两件事：瞬时应力将作为相应对角线分量的平均值进行计算，并且耦合盒子尺寸将以锁定步长一起更改，这意味着耦合尺寸将在每个时间步上以相同的百分比进行扩张或收缩。 任何耦合尺寸的Pstart，Pstop，Pdamp参数必须相同。 耦合xyz可以用于2d仿真； z维度只是被忽略。

——————————————————————————————————————

iso，aniso和tri关键词只是等同于一起指定其他几个关键字的快捷方式。

关键词iso意味着三个对角线分量一起耦合来计算静水压强，并且膨胀/收缩作用在三个维度上同时作用，当使用“iso Pstart Pstop Pdamp”时等同于设定以下四个关键词(添注：更加简单的理解，此时压强在三个方向是各向同性的，意味着Pstart和Pstop在三个方向都是一致的。)：

x Pstart Pstop Pdamp
y Pstart Pstop Pdamp
z Pstart Pstop Pdamp
couple xyz

关键词aniso意味着x，y和z尺寸是相互独立的进行控制的，分别采用Pxx， Pyy和Pzz。当使用“aniso Pstart Pstop Pdamp”时等同于设定以下四个关键词(添注：更加简单的理解，此时压强在三个方向是各向异性的，意味着Pstart和Pstop在三个方向都是相互独立的，可以大小不一样。):

x Pstart Pstop Pdamp
y Pstart Pstop Pdamp
z Pstart Pstop Pdamp
couple none

关键词tri意味着x，y，z，xy，xz和yz尺寸是分别独立的控制的，采用“tri Pstart Pstop Pdamp”默认采用以下七个关键词进行设定（添注：值得注意的是，tri适用于三斜非正交盒子，此时允许盒子的所有自由度。）：

x Pstart Pstop Pdamp
y Pstart Pstop Pdamp
z Pstart Pstop Pdamp
xy 0.0 0.0 Pdamp
yz 0.0 0.0 Pdamp
xz 0.0 0.0 Pdamp
couple none

——————————————————————————————————————

在某些情况下（譬如固体），当向系统施加Nose/Hoover控压和控温器的时候，系统的压强（体积）和温度会意料之外的剧烈振荡。此时，采用关键词drag可以阻缓这些振荡的发生，但是这将修改Nose/Hoover方程。当drag值设定为0.0时，意味着Nose/Hoover方程不发生改变。当采用非零值时，此时会向方程添加一个drag项；drag值设定的越大，阻尼效应越明显。短期执行并监视压强和温度是确定阻尼项是否起作用的最佳方法。 通常，介于0.2到2.0之间的值足以抑制几个周期后的振荡。 请注意，使用drag关键词将干扰能量守恒，并且还将更改位置和速度的分布，以使其与名义上的NVT，NPT或NPH系综不完全相符。

控制初始振荡的另一种方法是使用链式控温器。 关键词tchain用以确定粒子控温器中控温器的数量。 数值1对应于原始的Nose/Hoover控温器。 关键词pchain用以指定控温器链中控温器的自由度和控温器数量。 数值0对应于控压器变量中没有进行控温。

关键词mtk控制运动方程式中是否包括由Martyna，Tuckerman和Klein提到的校正项。 如果不指定，则将复制原始的Hoover控压器，其体积概率分布函数与真实的NPT和NPH系综相差1 / V。 因此，使用yes更为正确，但是在许多情况下，差异可以忽略不计。