关于分子动力学模拟中邻区列表算法的优化理论 侯吉旋 Laboratoire de Physique, UMR 5182 CNRS, Ecole Normale Suprieure de Lyon, 46, Alle dItalie, F-69364 Lyon Cedex 07, France 在过去的半个世纪里，分子动力学方法已经成功地应用到许多科学领域并取得了众多成果。但由于计算机的计算能力有限，大尺度的分子动力学模拟一直是一个难题。对于一个含有 N 个粒子的可加系统，每一步运算都需要计算 N(N-1)/2 个粒子相互作用。然而对于短程相互作用体系，如 Lennard-Jones 系统，每个粒子只与距离小于截断半径 R cut 范围内的粒子相互作用，因此在实际运算过程中只需要计算大约 (4 p R cut 3 r /3)N/2 个粒子相互作用即可， r 为系统的密度。可见大部分计算时间都浪费在对结果没有贡献的粒子间相互作用上。 为此在 1967 年 Verlet 采用了一种邻区列表算法，大大缩短了短程相互作用系统的计算机模拟的计算时间。在这个算法中，引入了一个比截断半径 R cut 稍大的列表半径 R list ，两者之差叫做皮肤半径 D R list -R cut ，见图 1 。在模拟的第一步，每个粒子的半径为 R list 的邻区内的粒子编号都储存在一个列表里，在接下来的运算中，我们只需要考虑该粒子与之相对应的列表中的粒子的相互作用，这样每步的运算量正比于 N 。直到有一个粒子的位移大于皮肤半径的一半，即 D /2 ，则列表需要更新，以免在列表半径以外的粒子进入到相互作用区域内，那么这一步的运算量正比于 N 2 。由于仅在需要更新列表的时候运算量与没有采用邻区列表算法时候的运算量相当，而其他步数都节省了很多运算时间，因此邻区列表算法大大加速了分子动力学模拟。 邻区列表算法示意图 皮肤半径 D 大小的选取直接影响了计算时间的长短。如果 D 太小，则列表需要经常更新，那么就无法节省计算时间。如果 D 太大，以至于列表里面几乎涵盖了整个体系大部分粒子，尽管列表不需要更新，但是每一步的计算量和不使用邻区列表的时候一样，也无法节约时间。因此需要选择一个最优的 D ，使得计算时间最小。 尽管邻区列表算法被广泛应用，但是很少有文章系统地研究过邻区列表算法的优化问题。我们提供了一种选择最优化参数的计算方法。 通过分别使用自由粒子近似和扩散近似对所需模拟的时间进行计算 , 再对两种近似计算进行比较 . 我们研究了更新间隔和皮肤半径 D 的关系。当 D 较小的时候，更新间隔是 D 的一次函数，对应于自由粒子近似；当 D 较大的时候，更新间隔是 D 的二次函数，对应了扩散近似。 更新间隔与皮肤半径的关系 the solid line is given by free particle approximation; the dashed lineis given bydiffusion approximation. 现在来看看我们最关心的计算机所消耗的时间。 下图显示了不同浓度下计算时间和皮肤半径大小的关系。正如所预期的，在皮肤半径 D 小的时候，自由粒子描述与模拟数据符合得很好，而在皮肤半径 D 大的时候扩散描述与模拟数据符合得很好。同时，在密度小的时候，自由粒子描述与大部分数据都很接近，而扩散描述只与小部分数据符合。而当密度升高，自由粒子描述与模拟数据的符合程度逐渐降低，而扩散描述与模拟数据的符合程度逐渐升高。从图中我们可以看到，能让计算时间最小的皮肤半径的值介于自由粒子描述的最优化点和扩散描述的最优化点之间。密度越小，实际模拟的最优化点和自由粒子描述的最优化点越接近。而密度越大，实际模拟的最优化点和扩散描述的最优化点越接近。 The CPU time as a function of the skin radius for different density. 有了我们的理论之后，做分子动力学模拟的你就不需要盲目的尝试了。如果你的系统处在低密度状态，使用自由粒子近似就知道怎么选择参数让计算时间最小。如果你的系统是高密度状态，那么用扩散近似就知道怎么选择参数了。一般情况下，最优的参数选择都会落在这两个近似给出的最优值之间。 References 侯吉旋，司黎明 , Optimization Theory for Neighbor List Algorithmin Fluid System Simulation . 物理化学学报 , 2009, 25 (03): 430-434.
Melipal 发表于2008-11-14 星期五 5:00 流体运动的多样性和不稳定性长期以来都是流体力学重要的话题。由于描述方程组往往是非线性的，对相当一部分情况，解析计算比较困难，因此通过实验或数值模拟来直接观察流体的行为也是重要的研究手段。下面这些图片都选自美国物理学会流体力学分会举办的年度流体运动图片展，有的是实验结果，也有计算机模拟图象。其意义不仅仅在于科学，更有着独特的欣赏性。 左图：风洞中的振荡流。图片提供：C. Stern 1 , S. P. R. Czitrom 2 , and R. Godoy 2 , (? 1 Facultad de Ciencias y, 2 Instituto de Ciencias del Mar y Limnologa, Universidad Nacional Autnoma de Mxico)? 1 、 水流流过球体后产生的旋涡 ：管中稳定的水流流过直径1厘米的球体时产生的周期性旋涡结构。激光照射出的流体形态显示，旋涡呈反向纤维状。实验时雷诺数的选取要将将使流体尾迹有周期性行为。本实验根据球体直径将雷诺数选为320，远小于湍流的临界雷诺数2000左右，但球体后方的尾流仍相当复杂。 图片提供：T. Leweke, M. Provansal, D. Ormires, and R. Lebescond (IRPHE, CNRS/Universits Aix-Marseille, France) 2 、 共轴射流附近的剪切不稳定性 ：共轴慢速圆射流和快速窄射流界面上产生的纵向和横向剪切不稳定性。下图的实验中两道射流的速度比选为3，雷诺数约为20000。由横向不稳定性产生的纵向旋涡可以有效地混合射流附近的流体。如进一步加大速度比，则出现类似不稳定尾流的振荡，不再表现为射流行为。 图片提供：E. Villermaux, H. Rehab, and E. J. Hopfinger (LEGI-CNRS, Institut de Mcanique de Grenoble, BP 53X, 38041 Grenoble Cedex, France) 3 、 低重力下的水膜球 ：在DC-9飞机模拟的低重力环境下观察到的水膜气球行为。水膜被针尖刺破后，表面向某一方向喷出飞沫，随后剩余部分会出现长久的振动。 图片提供：M. M. Weislogel 1 and S. Lichter 2 (NASA Lewis Research Center 1 , Northwestern University 2 ) 4 、 蓟花冠冕 ：将一滴染色的水滴滴到甘油层上得到的结果。水滴先是分为外环和内层两部分，随后表面张力使水层产生涡旋，波动又加强了张力梯度，产生图中所示的叶状外观（各图时间间隔约5秒）。 图片提供：E. Tan and S. T. Thoroddsen (University of Illinois at Urbana-Champaign) 5 、 湍流的色彩 ：该实验所采用的流体在不同温度下会呈现不同的色彩，因而颜色的变化就表现了湍流的热量传输过程，图样则与瞬时传输系数相关。图中展现的是在涡轮驱动的下行水流抵达固壁时的温度分布。 图片提供：D. R. Sabatino and T. J. Praisner (Lehigh University) 6 、 射流扰动形成的空穴 ：连续的水流射入水池中，导致池中出现充有空气的空穴。控制喷嘴处的阀门可以增强入射水流，这导致了鼓包的出现。照片中展现了鼓包在张力和引力作用下随时间演化的情况。实验中的雷诺数约为12300。 图片提供：Y. Zhu, H. N. Ogiz, and A. Prosperetti (The Johns Hopkins University) 7 、 涡旋的产生 ：涡旋可以通过抽动水流产生。下图两个涡旋分别由两端的盘状物旋转生成。染色的水流可以帮助人们理解湍流中涡丝的动力学行为。 图片提供：Philippe Petitjeans (Laboratoire de Physique et Mcanique des Milieus Htrognes, Ecole Suprieure de Physique te de Chimie Industrielles, Paris, France) 8 、 三维流体与固体的作用 ：使用微元法和欧拉-拉格朗日标记进行的三维流体模拟。流体流出后在引力作用下流向三条短堤围起的障碍物。模拟中表现了障碍物后方和周边区域的水波结构。 图片提供：Peter E. Raad and Razvan Bidoae (Southern Methodist University, Dallas, Texas) 9 、 微观流体的湍动 ：46微米厚的液晶层的流动。液晶夹在两块玻璃平板之间，两块玻璃板的内表面连以电极。当电场强度平缓增加时，液晶先是分成6束流（右上），再变为弱的湍流状态（左下）或是方格状对流元（右下），下侧两图中液晶的物理参数均处在混沌与湍流的过渡范围内。 图片提供：T. Peacock and T. Mullin (University of Manchester) 10 、 受迫固壁射流的双螺旋不稳定性 ：对平面固壁射流（Wall Jet）传输过程的分析。可视化采用了米（Mie）散射的方法。其中a、c两图为实验结果，d为实验装置介绍，b、e、f、g为计算机模拟结果。 图片提供：M. Visbal 1 , D. Gaitonde 1 , and S. Gogineni 2 ( 1 Air Force Research Laboratory, WrightPatterson AFB, 2 Innovative Scientific Solutions, Inc., Dayton, Ohio) 11 、 剪切层不稳定波与斜激波作用产生的声波 ：剪切层不稳定性与激波元的作用可以使超音速流产生噪音。下图是为了了解作用过程而构造的二维模型，由斜激波和超音速剪切层组成，剪切层中有不稳定的波动。图中黄色表示强压缩区域，红色为压缩区域，蓝色为膨胀区，灰色为音速区，绿色为涡度等高线。 图片提供：Ted A. Manning and Sanjiva K. Lele (Stanford University Department of Aeronautics and Astronautics, Stanford, California 94305-4035) 12 、 圆形振膜上水滴的雾化 ：1厘米直径的水滴在400毫秒的时间内被粉碎迅速。水滴表面的不稳定性引起了表面波，波峰喷射出更小的次级水滴，瓦解了原水滴，而原水滴与振膜的耦合会影响次级水滴的演化行为。 图片提供：Bojan Vukasinovic, Ari Glezer, and Marc K. Smith (Georgia Institute of Technology) 13 、 涡旋的不稳定并合 ：初始条件为两个共转的层流涡旋。当雷诺数增加至2000后，三维不稳定性使涡旋瓦解。下图上侧为并合前的情况，下侧为并合后，左侧为侧视图，右侧为附视图。 图片提供：P. Meunier and T. Leweke (CNRS/Universits Aix-Marseille, France) 14 、 负电雾化 ：室温下静电作用瓦解油滴的照片。充电由插入液体内的电极完成。停止通电后，静电作用导致了液体的雾化。 图片提供：Dimitris E. Nikitopoulos (Louisiana State University) and Arnold J. Kelly (CIC, Inc.) 15 、 涡旋射流旋转对称性的破坏 ：涡旋参数S是流体绕轴向与径向运动速度相对大小的度量。以下左图S取0.38，右图取0.49，可见轴向生成的不稳定结构。对于有涡旋存在的射流，不稳定性的成因与非涡旋流类似，但涡旋的存在放大了涡旋反方向的不稳定旋涡，却减小了涡旋方向的不稳定，而且还增大了轴向的不稳定形变。 图片提供：Thomas Loiseleux and Jean-Marc Chomaz, CNRS-Ecole polytechnique 16 、 粒子负载流的模拟 ：使用粒子模拟方法再现含有固体粒子的球形不可压缩流体在重力与粘滞力共同左右下的行为。图中红色代表流体，白色代表固体颗粒。其中上面三图考虑了初始的涡度，下面三图不计初始涡度。 图片提供：J. H. Walther (ETH Zrich), S.-S. Lee (Stanford) and P. Koumoutsakos (ETH Zrich and NASA Ames) 17 、 流体中的奇点 ：下图是甘油与水的混合液在垂直于表面且振幅接近临界值的正弦激发作用下，从表面波最低点的消失一直到向上运动的射流形成全过程的多次曝光照片。右侧的插图则是表面波瓦解的过程。由于瑞利不稳定性的存在，射流末端破裂成了液滴。 图片提供：by B. W. Zeff 1 , J. Fineberg 1,2 , and Daniel P. Lathrop 1 ( 1 University of Maryland, College Park, 2 Hebrew University of Jerusalem) 年度流体运动图片展的官方网站是： http://pof.aip.org/pof/gallery/ ，这项活动自从1985年来已经举办了20余届。 标签： 流体 , 艺术 转载原创文章请注明，转载自： 科学松鼠会 本文链接: http://songshuhui.net/archives/4090.html
标签: 流体
