预存
Document
当前位置:文库百科文章详情
ASE原子环境模拟之原子构建和计算单元
来源: 时间:2022-12-16 16:33:54 浏览:3725次

ASE允许指定不同的计算代码作为原子尺度的计算单元。在这个入门练习中,我们通过ASE的基本概念和工作流程,最终计算N2的结合曲线

这些教程经常使用电子结构代码GPAW。使用其他受支持的代码也可以很好地完成它们,但需要稍作调整。




Python


在ASE中,通过编写和运行Python脚本来执行计算。在ASE文档中可以找到一个非常简短的Python入门教程。如果你是Python新手,最好了解一下基本语法、数据类型和导入之类的东西。




Atoms


我们建立一个分子,进行DFT计算。我们可以手动输入原子的化学符号,并且猜测原子的坐标(以Ångström为单位),来建立简单分子。比如N2:


from ase import Atoms

atoms = Atoms('N2', positions=[[0, 0, -1], [0, 0, 1]])

为了防止错误,我们使用ASE GUI可视化我们的分子:


from ase.visualize import view

view(atoms)

同样地,我们以某种格式保存原子,通常是ASE自己的trajectory格式:


from ase.visualize import view

view(atoms)

然后在终端上运行GUI:


ase gui myatoms.traj

ASE支持相当多的不同格式。对于完整的列表,运行:

$ ase info --formats



Calculators


接下来让我们进行一个电子结构计算。ASE使用calculators模块进行计算。calculators模块为进行实际计算的不同后端提供抽象接口。通常,calculators通过调用外部电子结构代码或力场代码来工作。要运行计算,我们必须首先创建一个calculator,然后将其附加到Atoms对象。这里我们使用GPAW并设置一些计算参数:

from gpaw import GPAW

calc = GPAW(mode='lcao', basis='dzp', txt='gpaw.txt', xc='LDA')

atoms.calc = calc

不同的电子结构软件具有不同的输入参数。GPAW可以使用实空间网格(mode='fd')、平面波(mode='pw')或局部原子轨道(mode='lcao')来表示波函数。在这里,我们调用更快但不太精确的LCAO模式,以及标准双zeta极化基组('dzp')。GPAW和许多其他代码也需要一个基本原胞(或模拟盒子)。因此,我们把原子放在一个盒子的正中间,在每个原子周围留下3  Å的空间:

atoms.center(vacuum=3.0)

print(atoms)

输出将显示模拟框(或cell)的坐标,该框也可以在GUI中看到。

一旦Atoms有了一个具有适当参数的计算器,我们就可以计算能量和力:


e = atoms.get_potential_energy()

print('Energy', e)

f = atoms.get_forces()

print('Forces')

print(f)

这就给出了以eV为单位的能量eV和以eV/Å为单位的力。(ASE还提供了atoms.get_kinetic_energy()来表示原子运动时的动能。在DFT计算中,我们通常只需要Kohn-Sham基态能量,也就是计算器提供的“势能”。)

调用get_potential_energy()get_forces()会触发一个自洽计算,并输出大量的输出文本。检查gpaw.txt文件,可以查看选择了哪些计算参数。请注意get_forces()调用实际上并没有触发新的计算—力是根据已经计算过的基态得到的,因此我们只运行了一次计算。




结合曲线


ASE的优点是可以编写脚本。atoms.positions是表示原子位置的numpy数组:

print(atoms.positions)

我们可以通过在数组元素中添加或赋值来移动原子。然后,我们可以再次调用atoms.get_potential_energy()atoms.get_forces()来触发新的计算。

这样我们就可以实现任何一系列的计算。当运行多个计算时,我们通常希望将它们写入一个文件。我们可以使用标准的轨迹格式来编写多种计算(原子和能量),如下所示:


from ase.io.trajectory import Trajectory

traj = Trajectory('mytrajectory.traj', 'w')

...

traj.write(atoms)

练习:编写一个循环,以小的步长移动N2分子的其中一个原子,得到体系能量对距离的曲线。将每一步的构型存盘并且可视化。平衡距离是多少呢?

注意原子不要太靠近模拟盒子的边缘(否则电子会挤压盒子,增加能量和/或破坏计算)。

如果我们忘记编写轨迹,我们还可以在gpaw.txt文件上运行ASE GUI,尽管它的输出精度是有限的。

尽管GUI将为我们绘制能量曲线,但发布高质量的曲线通常需要一些手工修改。ASE提供两个函数来读取轨迹或其他文件:

  • ase.io.read()读取并返回最后一个构型,如果指定了index关键字,则可能是一组构型。 

  • ase.io.iread()读取多个构型,一次一个。

使用ase.io.iread()读取构型,例如:


for atoms in iread('mytrajectory.traj'):

   print(atoms)

练习:使用matplotlib绘制结合曲线(能量对距离的函数)。你将需要搜集随着轨迹循环的计算中每一步能量和距离。如果刚才已经计算过能量,那么atoms就已经包含能量,所以调用atoms_get_potential_energy()可以简单地获得能量而不需要再做计算。

选做练习:为了计算更加正确的结合能,设置一个孤立的N原子并计算它的能量。然后计算由原子到分子的分子结合能:E_结合 = E_N2 - 2E_N。你可以在计算开始之前设置atoms.set_initial_magnetic_moments([3])来告诉GPAW你的原子是具有自旋极化的。




练习答案


from ase import Atoms

from ase.io import Trajectory

from gpaw import GPAW

atoms = Atoms('N2', positions=[[0, 0, -1], [0, 0, 1]])

atoms.center(vacuum=3.0)

calc = GPAW(mode='lcao', basis='dzp', txt='gpaw.txt')

atoms.calc = calc

traj = Trajectory('binding_curve.traj', 'w')

step = 0.05

nsteps = int(3 / step)

for i in range(nsteps):

   d = 0.5 + i * step

   atoms.positions[1, 2] = atoms.positions[0, 2] + d

   atoms.center(vacuum=3.0)

   e = atoms.get_potential_energy()

   f = atoms.get_forces()

   print('distance, energy', d, e)

   print('force', f)

   traj.write(atoms)



import matplotlib.pyplot as plt

from ase.io import iread

energies = []

distances = []

for atoms in iread('binding_curve.traj'):

   energies.append(atoms.get_potential_energy())

   distances.append(atoms.positions[1, 2] - atoms.positions[0, 2])

ax = plt.gca()

ax.plot(distances, energies)

ax.set_xlabel('Distance [Å]')

ax.set_ylabel('Total energy [eV]')

plt.show()



from ase import Atoms

from gpaw import GPAW

atoms = Atoms('N')

atoms.center(vacuum=3.0)

atoms.set_initial_magnetic_moments([3])

calc = GPAW(mode='lcao', basis='dzp')

atoms.calc = calc

atoms.get_potential_energy()


计算狗·模拟计算

评论 / 文明上网理性发言
12条评论
全部评论 / 我的评论
最热 /  最新
全部 3小时前 四川
文字是人类用符号记录表达信息以传之久远的方式和工具。现代文字大多是记录语言的工具。人类往往先有口头的语言后产生书面文字,很多小语种,有语言但没有文字。文字的不同体现了国家和民族的书面表达的方式和思维不同。文字使人类进入有历史记录的文明社会。
点赞12
回复
全部
查看更多评论
相关文章

基础理论丨一文了解XPS(概念、定性定量分析、分析方法、谱线结构)

2020-05-03

晶体结构可视化软件 VESTA使用教程(下篇)

2021-01-22

手把手教你用ChemDraw 画化学结构式:基础篇

2021-06-19

【科研干货】电化学表征:循环伏安法详解(上)

2019-10-25

【科研干货】电化学表征:循环伏安法详解(下)

2019-10-25

Zeta电位的基本理论、测试方法和应用

2020-08-24

项目推荐/Project
第一性原理-电荷密度

第一性原理-电荷密度

第一性原理-扩散能垒

第一性原理-扩散能垒

热门文章/popular

基础理论丨一文了解XPS(概念、定性定量分析、分析方法、谱线结构)

晶体结构可视化软件 VESTA使用教程(下篇)

手把手教你用ChemDraw 画化学结构式:基础篇

【科研干货】电化学表征:循环伏安法详解(上)

电化学实验基础之电化学工作站篇 (二)三电极和两电极体系的搭建 和测试

【科研干货】电化学表征:循环伏安法详解(下)

微信扫码分享文章