ASE使用笔记(1)

本文主要是介绍ASE使用笔记(1)，希望对大家解决编程问题提供一定的参考价值，需要的开发者们随着小编来一起学习吧！

ASE使用笔记(1

- 写在前面
- 关于安装ASE
- ASE的简介
- 官网的第一个demo
- 关于ASE调用其他第一原理计算软件如VASP
- 对于材料结构上的操控(感觉还是和吸附有点像)
- 旋转操作
- CIF和XYZ,POSCAR的转换
- 总结

ASE(Atomic Simulation Environment) is a set of tools and Python modules for setting up, manipulating, running, visualizing and analyzing atomistic simulations. The code is freely available under the GNU LGPL license.

写在前面

作为一个刚刚由机器学习转入第一原理的人，对一些物理图像还是不怎么清晰的，但是，工欲善其事必先利其器，既然要做DFT相关的事情了，把和DFT相关的这些Python lib琢磨一下，应该会提高后期工作的效率。

关于安装ASE

pip install 即可，中间遇到的报错如VC++14之类的问题可以网上搜索解决办法进行规避or下载提示需要的VC++14,很多与和第一原理做接口的Python库都需要这个（如Pymatgen也遇到了这样的问题）
对于VC++14的下载：复制这段内容后打开百度网盘手机App，操作更方便哦
链接：https://pan.baidu.com/s/1I0OTvHfeNhaqrnY_nCqfyQ 提取码：8n85
官网链接：https://wiki.fysik.dtu.dk/ase/index.html

ASE的简介

官网解释的很清楚了，不再赘述。

官网的第一个demo

计算吸附能，氮气分子在铜板板上的
本笔记是在Jupyter Notebook的IPython环境下运行的。

#### 如何得到一个氮气分子
from ase import Atoms
d = 1.10
molecule = Atoms('2N', positions=[(0., 0., 0.), (0., 0., d)])
#### 如何得到铜板板
from ase.build import fcc111
slab = fcc111('Cu', size=(4,4,2), vacuum=10.0)#### 使用effective medium theory (EMT) calculator（有效介质理论）
from ase.calculators.emt import EMT
slab.calc = EMT()
molecule.calc = EMT()
#### 得到各自体系的孤立能量
e_slab = slab.get_potential_energy()
e_N2 = molecule.get_potential_energy()
#### 结构优化，这个例子是在on-top位置的
h = 1.85
add_adsorbate(slab, molecule, h, 'ontop')#### 固定（fix）铜板板是不能动的，只有N可以调整
from ase.constraints import FixAtoms
constraint = FixAtoms(mask=[a.symbol != 'N' for a in slab])
slab.set_constraint(constraint)
#### 设置受力收敛标准为0.005
from ase.optimize import QuasiNewton
dyn = QuasiNewton(slab, trajectory='N2Cu.traj')
dyn.run(fmax=0.00005)

                Step[ FC]     Time          Energy          fmax
*Force-consistent energies used in optimization.
BFGSLineSearch:    0[  0] 14:33:13       11.689927*       1.0797
BFGSLineSearch:    1[  2] 14:33:13       11.670814*       0.4090
BFGSLineSearch:    2[  4] 14:33:13       11.625880*       0.0409
BFGSLineSearch:    3[  6] 14:33:13       11.625228*       0.0309
BFGSLineSearch:    4[  7] 14:33:13       11.625212*       0.0039
BFGSLineSearch:    5[  8] 14:33:13       11.625212*       0.0000

## 把结构信息存储和读取
from ase.io import write
write('slab.xyz', slab)
from ase.io import read
slab_from_file = read('slab.xyz')
## 结构的可视化
from ase.visualize import view
view(slab)

结构的可视化

## 分子动力学，但是他的温度设置是如何的没有仔细看
from ase.md.verlet import VelocityVerlet
from ase import units
dyn = VelocityVerlet(molecule, dt=1.0 * units.fs)
for i in range(100):pot = molecule.get_potential_energy()kin = molecule.get_kinetic_energy()print('%2d: %.5f eV, %.5f eV, %.5f eV' % (i, pot + kin, pot, kin))dyn.run(steps=20)
0: 0.45841 eV, 0.32004 eV, 0.13837 eV1: 0.45945 eV, 0.41613 eV, 0.04332 eV2: 0.45790 eV, 0.31736 eV, 0.14054 eV3: 0.45831 eV, 0.34529 eV, 0.11302 eV4: 0.45924 eV, 0.39462 eV, 0.06462 eV5: 0.45870 eV, 0.34393 eV, 0.11477 eV6: 0.45933 eV, 0.40648 eV, 0.05285 eV7: 0.45741 eV, 0.27477 eV, 0.18264 eV8: 0.45975 eV, 0.44898 eV, 0.01077 eV9: 0.45782 eV, 0.27953 eV, 0.17828 eV
10: 0.46031 eV, 0.45943 eV, 0.00088 eV
11: 0.45762 eV, 0.26915 eV, 0.18847 eV
12: 0.45981 eV, 0.45633 eV, 0.00348 eV
13: 0.45737 eV, 0.26596 eV, 0.19141 eV
14: 0.45971 eV, 0.42801 eV, 0.03171 eV
15: 0.45845 eV, 0.32330 eV, 0.13515 eV
16: 0.45943 eV, 0.41328 eV, 0.04615 eV
17: 0.45795 eV, 0.32109 eV, 0.13686 eV
18: 0.45825 eV, 0.34118 eV, 0.11707 eV
19: 0.45928 eV, 0.39782 eV, 0.06146 eV
20: 0.45866 eV, 0.34048 eV, 0.11818 eV
21: 0.45939 eV, 0.41026 eV, 0.04913 eV
22: 0.45740 eV, 0.27294 eV, 0.18446 eV
23: 0.45976 eV, 0.45045 eV, 0.00931 eV
24: 0.45778 eV, 0.27753 eV, 0.18025 eV
25: 0.46032 eV, 0.45992 eV, 0.00040 eV
26: 0.45765 eV, 0.27054 eV, 0.18711 eV
27: 0.45981 eV, 0.45540 eV, 0.00441 eV
28: 0.45737 eV, 0.26702 eV, 0.19035 eV
29: 0.45965 eV, 0.42475 eV, 0.03491 eV
30: 0.45849 eV, 0.32661 eV, 0.13188 eV
31: 0.45940 eV, 0.41036 eV, 0.04903 eV
32: 0.45801 eV, 0.32490 eV, 0.13311 eV
33: 0.45819 eV, 0.33711 eV, 0.12107 eV
34: 0.45931 eV, 0.40097 eV, 0.05834 eV
35: 0.45862 eV, 0.33705 eV, 0.12157 eV
36: 0.45946 eV, 0.41395 eV, 0.04551 eV
37: 0.45739 eV, 0.27126 eV, 0.18612 eV
38: 0.45978 eV, 0.45182 eV, 0.00795 eV
39: 0.45775 eV, 0.27564 eV, 0.18211 eV
40: 0.46032 eV, 0.46022 eV, 0.00010 eV
41: 0.45768 eV, 0.27206 eV, 0.18562 eV
42: 0.45980 eV, 0.45435 eV, 0.00544 eV
43: 0.45737 eV, 0.26823 eV, 0.18914 eV
44: 0.45959 eV, 0.42137 eV, 0.03822 eV
45: 0.45853 eV, 0.32995 eV, 0.12858 eV
46: 0.45937 eV, 0.40739 eV, 0.05198 eV
47: 0.45806 eV, 0.32878 eV, 0.12928 eV
48: 0.45813 eV, 0.33310 eV, 0.12503 eV
49: 0.45934 eV, 0.40406 eV, 0.05527 eV
50: 0.45858 eV, 0.33364 eV, 0.12494 eV
51: 0.45952 eV, 0.41755 eV, 0.04197 eV
52: 0.45738 eV, 0.26973 eV, 0.18765 eV
53: 0.45979 eV, 0.45309 eV, 0.00670 eV
54: 0.45772 eV, 0.27386 eV, 0.18385 eV
55: 0.46032 eV, 0.46032 eV, 0.00000 eV
56: 0.45771 eV, 0.27370 eV, 0.18401 eV
57: 0.45979 eV, 0.45320 eV, 0.00658 eV
58: 0.45738 eV, 0.26960 eV, 0.18778 eV
59: 0.45953 eV, 0.41788 eV, 0.04164 eV
60: 0.45857 eV, 0.33332 eV, 0.12525 eV
61: 0.45934 eV, 0.40435 eV, 0.05499 eV
62: 0.45812 eV, 0.33273 eV, 0.12539 eV
63: 0.45807 eV, 0.32915 eV, 0.12892 eV
64: 0.45937 eV, 0.40710 eV, 0.05226 eV
65: 0.45854 eV, 0.33026 eV, 0.12827 eV
66: 0.45959 eV, 0.42105 eV, 0.03854 eV
67: 0.45737 eV, 0.26835 eV, 0.18902 eV
68: 0.45980 eV, 0.45425 eV, 0.00555 eV
69: 0.45768 eV, 0.27221 eV, 0.18548 eV
70: 0.46032 eV, 0.46024 eV, 0.00009 eV
71: 0.45775 eV, 0.27547 eV, 0.18228 eV
72: 0.45978 eV, 0.45195 eV, 0.00783 eV
73: 0.45739 eV, 0.27111 eV, 0.18627 eV
74: 0.45946 eV, 0.41429 eV, 0.04517 eV
75: 0.45861 eV, 0.33673 eV, 0.12189 eV
76: 0.45931 eV, 0.40126 eV, 0.05805 eV
77: 0.45818 eV, 0.33673 eV, 0.12145 eV
78: 0.45801 eV, 0.32526 eV, 0.13275 eV
79: 0.45939 eV, 0.41009 eV, 0.04931 eV
80: 0.45849 eV, 0.32692 eV, 0.13158 eV
81: 0.45965 eV, 0.42444 eV, 0.03521 eV
82: 0.45737 eV, 0.26713 eV, 0.19024 eV
83: 0.45981 eV, 0.45530 eV, 0.00450 eV
84: 0.45765 eV, 0.27068 eV, 0.18698 eV
85: 0.46032 eV, 0.45996 eV, 0.00036 eV
86: 0.45778 eV, 0.27735 eV, 0.18043 eV
87: 0.45976 eV, 0.45058 eV, 0.00918 eV
88: 0.45740 eV, 0.27278 eV, 0.18462 eV
89: 0.45940 eV, 0.41061 eV, 0.04879 eV
90: 0.45865 eV, 0.34016 eV, 0.11850 eV
91: 0.45928 eV, 0.39811 eV, 0.06116 eV
92: 0.45824 eV, 0.34079 eV, 0.11745 eV
93: 0.45796 eV, 0.32145 eV, 0.13651 eV
94: 0.45942 eV, 0.41301 eV, 0.04641 eV
95: 0.45845 eV, 0.32361 eV, 0.13484 eV
96: 0.45971 eV, 0.42771 eV, 0.03200 eV
97: 0.45737 eV, 0.26605 eV, 0.19132 eV
98: 0.45981 eV, 0.45625 eV, 0.00356 eV
99: 0.45762 eV, 0.26927 eV, 0.18835 eV

关于分子的一些优化和计算，本人不研究非周期体系，有时间了可以看一下

关于ASE调用其他第一原理计算软件如VASP

$ export ASE_VASP_COMMAND="mpiexec $HOME/vasp/bin/vasp_std"   ##可执行路径
$ export VASP_PP_PATH=$HOME/vasp/mypps   ##yanshiPP
$ export ASE_VASP_VDW=$HOME/<path-to-vdw_kernel.bindat-folder>  ##vdw文件路径

from ase.build import molecule
atoms = molecule('N2')    ## 这里还是拿个氮气作为例子
atoms.center(vacuum=5)
from ase.calculators.vasp import Vasp2calc = Vasp2(xc='pbe',  # Select exchange-correlation functionalencut=520, # Plane-wave cutoffkpts=(1, 1, 1)) # k-pointsatoms.calc = calc
en = atoms.get_potential_energy()  # This call will start the calculation
print('Potential energy: {:.2f} eV'.format(en))Potential energy: -16.59 eV

对于材料结构上的操控(感觉还是和吸附有点像)

from math import sqrt
from ase import Atoms
a = 3.55
atoms = Atoms('Ni4',cell=[sqrt(2) * a, sqrt(2) * a, 1.0, 90, 90, 120],pbc=(1, 1, 0),scaled_positions=[(0, 0, 0),(0.5, 0, 0),(0, 0.5, 0),(0.5, 0.5, 0)])
atoms.center(vacuum=5.0, axis=2)

atoms.cell

Cell([[5.020458146424487, 0.0, 0.0], [-2.5102290732122423, 4.347844293440141, 0.0], [0.0, 0.0, 10.0]])

atoms.positions

array([[ 0.        ,  0.        ,  5.        ],[ 2.51022907,  0.        ,  5.        ],[-1.25511454,  2.17392215,  5.        ],[ 1.25511454,  2.17392215,  5.        ]])

行吧，和POSCAR基本吻合

from ase.visualize import view
atoms.write('slab.xyz')
view(atoms)

在这里插入图片描述

还可以看看扩胞后的，在命令行上使用的话需要输入

ase gui -r 3,3,2 slab.xyz

当然，我这个环境需要在命令行前加入一个感叹号

!ase gui -r 3,3,2 slab.xyz   ## 扩胞查看

## 如果加点东西进去（先加个原子进去）
h = 1.9
relative = (1 / 6, 1 / 6, 0.5)
absolute = np.dot(relative, atoms.cell) + (0, 0, h)
atoms.append('Ag')
atoms.positions[-1] = absolute
view(atoms)

在这里插入图片描述

## 创建了一个水分子进来
import numpy as np
from ase import Atoms
p = np.array([[0.27802511, -0.07732213, 13.46649107],[0.91833251, -1.02565868, 13.41456626],[0.91865997, 0.87076761, 13.41228287],[1.85572027, 2.37336781, 13.56440907],[3.13987926, 2.3633134, 13.4327577],[1.77566079, 2.37150862, 14.66528237],[4.52240322, 2.35264513, 13.37435864],[5.16892729, 1.40357034, 13.42661052],[5.15567324, 3.30068395, 13.4305779],[6.10183518, -0.0738656, 13.27945071],[7.3856151, -0.07438536, 13.40814585],[6.01881192, -0.08627583, 12.1789428]])
c = np.array([[8.490373, 0., 0.],[0., 4.901919, 0.],[0., 0., 26.93236]])
W = Atoms('4(OH2)', positions=p, cell=c, pbc=[1, 1, 0])
W.write('WL.traj')

在这里插入图片描述

旋转操作

W.rotate(90, 'z', center=(0, 0, 0))
view(W)
W.wrap()
###########The wrap() method only works if periodic boundary conditions are enabled. We have a 2 percent lattice mismatch between Ni(111) and the water, so we scale the water in the plane to match the cell of the slab. The argument scale_atoms=True indicates that the atomic positions should be scaled with the unit cell. The default is scale_atoms=False indicating that the cartesian coordinates remain the same when the cell is changed
W.set_cell(slab.cell, scale_atoms=True)
zmin = W.positions[:, 2].min()
zmax = slab.positions[:, 2].max()
W.positions += (0, 0, zmax - zmin + 1.5)
interface = slab + W
interface.center(vacuum=6, axis=2)
interface.write('NiH2O.traj')
view(interface)

在这里插入图片描述

CIF和XYZ,POSCAR的转换

ASE的这个功能已经被Pymatgen整合了，两者用起来都很方便。
拿出你的一个cif文件，我这个从MP里直接下载了个钙钛矿。

f = open('TiPbO3.cif', 'r')
print(f.read())

# generated using pymatgen
data_TiPbO3
_symmetry_space_group_name_H-M   'P 1'
_cell_length_a   3.86121500
_cell_length_b   3.86121500
_cell_length_c   4.62706100
_cell_angle_alpha   90.00000000
_cell_angle_beta   90.00000000
_cell_angle_gamma   90.00000000
_symmetry_Int_Tables_number   1
_chemical_formula_structural   TiPbO3
_chemical_formula_sum   'Ti1 Pb1 O3'
_cell_volume   68.98476581
_cell_formula_units_Z   1
loop__symmetry_equiv_pos_site_id_symmetry_equiv_pos_as_xyz1  'x, y, z'
loop__atom_site_type_symbol_atom_site_label_atom_site_symmetry_multiplicity_atom_site_fract_x_atom_site_fract_y_atom_site_fract_z_atom_site_occupancyTi  Ti0  1  0.50000000  0.50000000  0.52092300  1Pb  Pb1  1  0.00000000  0.00000000  0.96921200  1O  O2  1  0.50000000  0.50000000  0.14244900  1O  O3  1  0.00000000  0.50000000  0.62665800  1O  O4  1  0.50000000  0.00000000  0.62665800  1

import warnings
warnings.filterwarnings("ignore")
!ase gui TiPbO3.cif -o POSCAR  ##这个参数 -o表示out 即输出为什么
f = open('POSCAR', 'r')
print(f.read())

Ti Pb  O 1.00000000000000003.8612150000000001    0.0000000000000000    0.00000000000000000.0000000000000000    3.8612150000000001    0.00000000000000000.0000000000000000    0.0000000000000000    4.62706100000000031   1   3
Cartesian1.9306075000000000  1.9306075000000000  2.41034249730300010.0000000000000000  0.0000000000000000  4.48460304593200031.9306075000000000  1.9306075000000000  0.65912021238900000.0000000000000000  1.9306075000000000  2.89958479213800051.9306075000000000  0.0000000000000000  2.8995847921380005

发现这个POSCAR和我要的版本不一样，我的vasp是5.4.4，后来发现大师兄的一个博客提到了这个东西是需要去修改ASE的源代码的，大师兄的博客连接https://www.bigbrosci.com/2020/08/31/A17/
1）找到ASE的安装目录，打开 ase/io 目录下的vasp5.py 文件；

2）将vasp.py中vasp5 = False 全部改为vasp5 = True；

3）保存vasp.py 即可。

大师兄牛皮，解决问题。