Lattice Lattice NCOUL OFF
Lattice Lattice COUL OFF
Probe Lattice NCOUL OFF
Probe Lattice COUL SUM FAST FIXED EWALD SFACTOR KMAX@15 KAPPA@6.7 LOCUT@1e-10
Probe Probe NCOUL OFF
Probe Probe COUL OFF
首先需要指出的是, 这里除了 Probe 与 Lattice的COUL作用之外, 其他的相互作用全部关闭.
来仔细看看这句:
Probe Lattice COUL SUM FAST FIXED EWALD SFACTOR KMAX@15 KAPPA@6.7 LOCUT@1e10
(1) “Probe” 和 “Lattice”, 这个就像“路人甲”一样, 你想怎么称呼都可以;
(2) “COUL”是从 ff.F90过来的, 在那个源程序里面定义了体系的能量从COUL和NONCOUL计算得到;
(3) “SUM”, “FAST”, “FIXED”,“EWALD”来自于ssdrive.F90, 我们来一一解读一下.
(a) SUM
Line 94, 定于了 “SUM”: “Type(SSSumParams), Pointer :: sum”. 很容易理解, 在任意位置上, Probe与Lattice的COUL作用, 需要累加起来.
(b) FAST
Line 495,是注释行, 关于“FAST”是这样说的, “fast -- True (False) => evaluate "Fast" (Slow) interactions”
字面上的理解就是, 用了“FAST”那么计算的时候, 关于相互作用的计算就能“加速”,至于如何加速? 不知道, 也不需要知道.
(c) FIXED
没有从源程序找到很好的对应, 但是这个是和Control文件中对于Lattice的处理对应的, 我们将Lattice文件Fixed处理, 那么这里也是需要告诉mapmake文件这个信息. 除此之外, 没啥用处.
(d) EWALD
这个指定了如果生成EMAP, 前面有提到过, GCMC的EMAP只能通过EWALD得到, 而on-the-fly COUL计算(i.e. 不用EMAP, 在GCMC中直接计算COUL)只能WFCOUL而不能使用EWALD.
Line542再次证明了这一点认识: “on the fly ewald sums are implemented for md only”
也就是说, MD时, on the fly COUL计算, 是可以采用EWALD SUM的, 但是MC(包括GCMC)不行....
(4) SFACTOR, KMAX, KAPPA, LOCUT是EWALD方法中的参数. 具体的原理,参见Frenkel&Smit 的第十二章, 或者MuSiC的相关文档. 你也可以从ewald.F90 的subroutine (Line 152 to Line 269) ewald_init看到程序中是怎样指定这些参数的. 另外, sssum.F90也很值得一看, 如果你想了解更多EWALD是怎样在MuSiC中实现的话.
下面一一分析一下这四个参数.
(a) SFACTOR
Line 205 to Line 207:
Case ('sfactor')
!** We want to use the structure factor.
eparams%useSF = .TRUE.
它是一个逻辑变量, 默认值是.FALSE., 所以在这里需要加入这个关键词来激活它, 变为.TRUE.
激活这个逻辑变量以后, 程序做了一些什么呢? 参加 Line 284 那边的注释.
!** Call the structure factor routine. This calculates the k space
!** vectors we will use later when computing the Ewald summation. If
!** the coordinates of the ions are fixed, this needs to be calculated
!** only once. If they change, it must be recalculated.
(b) KMAX
Line 184: “Number of k "boxes" in each direction (x,y,z)”. 你可能注意到在程序里面, KAMX is hard coded to be 5000, 你可能回想增大这个数值, 来处理更大的体系. 但是, 不幸的是, 你的内存很有可能比KMAX参数更快地达到阀值, 所以, 修改这个参数在实际使用上, 没有意义, 我没有测试过, 但是我觉得这里 KMAX@15, 其实15这个数值, 是没什么意义的, 只是需要一个input而已, 你用15, 10还是8, 没有区别.
(c) KAPPA 是和 EWALD SUM里面的Gaussian function有关的一个参数, 它定义了 Gaussian Function的带宽.
如果你看Frenkel&Smit 的第十二章, 这里程序中的KAPPA参数值等于Frenkel&Smit Chap.12中的alpha^0.5.
(d) LOCUT
这个很好理解, 请看ewald.F90文件中的这一行:
Line 200, “ !** Low end cutoff for the reciprocal space k vectors.”
接下来的问题就是, 我们应该如何选择合适的KAPPA和KMAX参数? 这两个参数的选择是基于EWALD SUM 计算时对于实空间和倒空间的“利用”和平衡. 一般而言, KMAX必须让SUMMATION在倒空间的计算满足随着KAPPA参数收敛, 同时尽量让实空间的计算也能快速收敛. 比方说, 减小KAPPA参数可以使得倒空间的收敛加快, 但是此时实空间的收敛就会减慢. 在非常小的KAPPA参数下, 实空间的EWALD SUM 和WFCOUL方法一样. 在其他MAP MAKE 参数相同的情况下, 测试不同的KAPPA数值, 看生成MAP需要的不同时间, 来找到最优值, 但是一般来讲, 这个对于MC(包括GCMC, 和MC里面相应的EMAP的生成)没有实际意义. 但是对于MD来讲, 由于on the fly COUL计算可以采用EWALD SUM 方法, 这样的测试还是能有明显提速功效的.
好了, 对于生成EMAP时EWALD 方法的介绍差不多了解这些就足够了. 还有值得提醒的就是, MuSiC中支持的EWALD 方法只是对正交盒子有效 (orthorhombic unit cell). 这个对于GCMC问题不大, zeolite, MOF等材料, 大多数还是orthombic的, 但是如果是自己生成的初始结构, 想要算吸附, 然后又想要使用EMAP的时候, 就需要注意这一点!