计算化学公社

标题: 如何基于从头算结合计算几何计算元素在两相间的分配比 [打印本页]

作者
Author: nusiew 时间: 2021-8-24 19:58
标题: 如何基于从头算结合计算几何计算元素在两相间的分配比
本帖最后由 nusiew 于 2021-8-26 08:59 编辑

请问大家有没有听说过 hull 和 inhedron 两个程序？

各位老师、同学、前辈，大家晚上好！最近看文献，发现一种利用轨迹输出文件（由从头算分子模拟得到），结合计算几何，计算元素在两相间分配比的方法。对于其利用计算几何计算分配比的部分非常感兴趣，在此向大家求教。论文作者关注的领域是地核与岩浆中元素的分配，采用类似的方法一共发表了3篇文章，最近的一篇是2020年，在此举其2012年发表在PNAS上的一篇对其研究方法进行说明。

作者利用vasp进行第一性原理分子动力学模拟，想要计算碳元素在铁水和硅酸盐熔体两相间的分配比。模拟和计算的过程是：先将260个原子（包含125个O、51个Mg、44个Fe 、36个Si和4个C原子）随机均匀分布在模拟盒子中，跑若干步后盒子中的原子会合并成铁水和硅酸熔体两相（如下图所示），之后便是利用计算几何的方法，结合 hull 和 inhedron 两个程序对Fe相成分和硅酸盐熔体组成进行了计算。不知道大家有没有听说过 hull 和 inhedron 这两个程序？由于本人对计算几何不太了解，对论文中计算几何部分的操作步骤也不熟悉，不知道有没有相关方面的教程或书籍？希望大家多多指导。论文中具体的成分分析过程大致如下。

在平衡之后的每个 FPMD 步骤中，首先为 Fe 簇构建多面体，然后利用计算几何学中的概念和方法判断原子是否被包围在多面体中，从而获得包围在 Fe 相中的原子数。基本上，需要三个步骤。首先，计算 Fe 簇的质心并将其移动到模拟单元的中心。通过使用模拟单元的周期性边界条件，所有原子都相应地移动。这是使用计算几何的准备步骤。其次，通过使用程序“hull”，它采用随机增量算法。第三，通过使用随机射线交叉算法的程序“inhedron”确定原子是否在第二步中构造的多面体内部。在每个 FPMD 步骤中，Fe 相中的原子数在模拟的整个积累周期内相加并求平均值，以计算 Fe 相的组成。由于整个系统的整体组成是已知的，因此可以直接推导出硅酸盐熔体的组成。

模拟得到的两相，Fe相，初始成分和原子数，以及平衡两相成分和原子数如下：

(, 下载次数 Times of downloads: 29)
具体的计算方法大家可以参考这三篇文献
2012：Carbon and other light element contents in the Earth’s core based on first-principles molecular dynamics
https://www.pnas.org/content/109/48/19579

2009：Partitioning of Si and O between liquid iron and silicate melt: A two-phase ab-initio molecular dynamics study（分析过程最详细）
https://agupubs.onlinelibrary.wi ... 0.1029/2009GL039751

2020：Magnesium partitioning between silicate melt and liquid iron using first-principles molecular dynamics: Implications for the early thermal history of the Earth's core
https://www.sciencedirect.com/sc ... 821X19306260#br0620

作者
Author: Charlie-Ma 时间: 2021-8-26 10:39
本帖最后由 Charlie-Ma 于 2021-8-27 12:52 编辑

#关于inhedron程序
inhedron程序可以用来查询点与多面体的相对位置，即在多面体内部、外部还是表面。在模拟计算中可以用于判断原子与相或者团簇的相对位置。
inhedron程序的运行方式如下：
1. 首先，它会询问您多面体有多少个顶点。
2. 然后它让你输入每个顶点的 XYZ 坐标。
3. 然后它会询问您多面体有多少个三角形面。
4. 然后它让你为每个面输入三个索引（顶点从 0 开始编号，而不是 1）。
5. 然后它会要求您输入任意数量点的 XYZ 坐标，inhedron确定每个点是在多面体的内部、外部还是正好位于多面体上。

输出结果中
'V'：查询点a与多面体P的顶点重合。
'E'：查询点a在多面体P的边的相对内部。
'F'：查询点a位于多面体P的面的相对内部。
'i'：查询点a严格位于多面体P 的内部。
'o'：查询点a严格位于多面体P 的外部（或外部）。
该链接中对其如何使用说明的很详细，望参考
http://ptinpoly.pbworks.com/w/page/19749262/O'Rourke%20C%20Code

#关于hull程序
程序可以从此处下载：jhsa26/hull (github.com)
http://www.qhull.org/

Matlab中的alphashape函数也可以实现改需求
https://ww2.mathworks.cn/help/matlab/ref/alphashape.html

欢迎光临计算化学公社 (http://bbs.keinsci.com/)