快速寻找近邻关系的计算方法

  • strict warning: Non-static method view::load() should not be called statically in /home/vasp/wwwroot/drupal-6.33/sites/all/modules/views/views.module on line 906.
  • strict warning: Declaration of views_handler_argument::init() should be compatible with views_handler::init(&$view, $options) in /home/vasp/wwwroot/drupal-6.33/sites/all/modules/views/handlers/views_handler_argument.inc on line 744.
  • strict warning: Declaration of views_plugin_row::options_validate() should be compatible with views_plugin::options_validate(&$form, &$form_state) in /home/vasp/wwwroot/drupal-6.33/sites/all/modules/views/plugins/views_plugin_row.inc on line 134.
  • strict warning: Declaration of views_plugin_row::options_submit() should be compatible with views_plugin::options_submit(&$form, &$form_state) in /home/vasp/wwwroot/drupal-6.33/sites/all/modules/views/plugins/views_plugin_row.inc on line 134.

 在绝大多数的计算方法中,寻找近邻关系都是一个必要的环节。虽然不需要复杂的数学知识,但是如何做到高效寻找仍然是一个有趣的问题。其寻找效率对于分子动力学以及紧束缚方法而言均有比较大的影响。本文仅就这一问题作一个粗浅的探讨。

对于周期性保持得比较完好的晶体,特别是正交单胞的情况,可以利用MIC(minimum image convention)方法比较快地解决。考虑任意两个原子k和l,其间距约化为:
R(kl) = R(kl)-L*INT(R(kl)/L)。之后再和截断半径比较即得结果。如果基矢不正交,那么MIC应用起来会有失效的危险。一种可能的扩展方法是找出每个基矢在其他两个基矢所决定的平面法线上的投影,用它除截断半径,从而近似得出需要考虑的格点(边界处需要特别对待)。

但是对于分子动力学或者紧束缚所经常处理的原子数目极大的缺陷体系而言,上述方法并不适用,因为MIC无法处理超单胞的情况。最直接的方法当然是遍历法,利用两层循环比较每一对原子的距离,考虑小于截断半径的原子对。但是这种方法的所需时间是O(n^2)。虽然实现起来极为简单,但是实际测试表明当体系包含有5000个原子时,就需要将近10分钟去完成检验(SGI Origin3000,单CPU)。因此,需要寻找更好的处理方法。与高效排序算法类似,如果可以将这个超大的体系分而治之(Divide & Conquer),每个原子只需要在包括它自己所在的区域在内的最近邻的27个区域内寻找最近邻原子。这样,寻找最近邻就优化为O(n)复杂度的问题了。在最好的情况下,70000个原子的体系可以在4秒左右的时间内完成。要获得比较好的计算效率,首要的问题在于如何确定这些分区的大小,或者更直接一点,如何保证这些分区包含有合适的原子个数。如果只需要最近邻,那么分区的边长可以定为1.2*(Volume/Natom)^(1/3)。其中Volume是体系的体积,Natom是原子总数。如果工作需要寻找给定截断半径Rcut内的近邻原子,分区的边长即可设为Rcut。除去某些基矢间夹角过小的情况外,这种设定对于大多数体系都是有效的,也是安全的。

除此而外,一些编程上的细节也需要强调一下。

1、这种算法需要构建三个映射关系:给定原子存在于哪个分区、给定分区包含哪些原子,以及给定分区的最近邻26个分区的编号。前两个关系可以同时给出。非正交的情况下,只要计算原子相对于体系左下角的距离在各基矢上的投影,再除以分区边长即可。

2、对于周期性边界条件,需要另外指定数组,用于标定MIC中基矢的平移次数。多数情况下为+1或-1。

在我完善测试程序之后,会将原程序贴到网上。希望会给大家一些启示。