空间站广场

论文

Notebooks

比赛

课程

Apps

我的主页

我的Notebooks

我的论文库

我的足迹

我的工作空间

任务

节点

文件

数据集

镜像

项目

数据库

公开

ABACUS DFT+U计算模拟实例：TiO2极化子（TiO2−）能态密度的计算模拟

notebook

bohr47b42a

更新于 2024-06-21

推荐镜像 :ABACUS:3.6.4

推荐机型 :c2_m4_cpu

ABACUS DFT+U计算模拟实例：TiO2极化子（TiO2−）能态密度的计算模拟

目标

背景

1. 使用VESTA手动构建TiO2单胞晶体结构（参照《VESTA软件的晶体学可视化实验教学设计与应用》）

1.1 查询TiO2晶体结构参数

1.2 VESTA生成TiO2晶体

2. TiO2结构优化(relax)

2.1 STRU文件中的输入参数设置

2.2 INPUT文件中的输入参数设置

2.3 KPT文件中的输入参数设置

2.4 运行ABACUS

3. TiO2自洽计算(scf)

3.1 STRU文件中的输入参数设置

3.2 INPUT文件中的输入参数设置

3.3 KPT文件中的输入参数设置

3.4 运行ABACUS

4. TiO2非自洽计算（输出能态密度）

4.1 INPUT文件中的输入参数设置

4.2 KPT文件中的输入参数设置

4.3 运行ABACUS

5. 画能态密度TDOS/PDOS

5.1 安装abacus-plot

5.2 画总能态密度TDOS

5.2.1 不同U值的TDOS

5.2.2 对比VASP结果

5.3 画投影能态密度PDOS

6.参考文献

ABACUS DFT+U计算模拟实例：TiO2极化子（TiO2 $^{-}$ ）能态密度的计算模拟

©️ Copyright 2024 @ Authors
作者：陈俊含
邮箱：chenjh@aisi.ac.cn
创建日期：2024-06-20
最近修改：2023-06-21
ABACUS版本：v.3.6.5
共享协议：本作品采用知识共享署名-非商业性使用-相同方式共享 4.0 国际许可协议进行许可。

代码

文本

🎯 本教程旨在快速熟悉怎样使用ABACUS中的DFT+U功能实现在半导体材料TiO2极化子的能态密度计算。

面向对使用ABACUS中DFT+U功能进行半导体材料极化子研究感兴趣的读者。
一键运行，可以快速在实践中检验您的想法。
丰富完善的注释，对于入门者友好。

代码

文本

目标

上手熟悉ABACUS中的DFT+U运算功能

在学习本教程后，您将能够：

学会在ABACUS中运行DFT+U功能所需的输入参数设置
在TiO2 $^{-}$ 体系上对ABACUS中的DFT+U功能进行标定

阅读该教程【最多】约需 30 分钟，让我们开始吧！

代码

文本

背景

向过渡金属氧化物半导体中掺入额外的一个电子（或空穴），由于过渡金属d轨道上的电子有强关联相互作用，很容易使得额外的电子局域在某一个或者某几个相邻的过渡金属原子上（与此相关的话题还有Anderson局域化、Mott绝缘体)。过渡金属氧化物中极化子的形成和迁移是一个非常重要的研究课题,对过渡金属氧化物材料性质有直接的影响。

代码

文本

1. 使用VESTA手动构建TiO2单胞晶体结构（参照《VESTA软件的晶体学可视化实验教学设计与应用》）

代码

文本

1.1 查询TiO2晶体结构参数

代码

文本

金红石二氧化钛，四方晶系，空间群P42/mnm，晶格参数a=b=0.4653nm, c=0.2969nm, α=β=γ=90°。每个晶胞中含有两个Ti原子和4个O原子，Ti1(0,0,0),Ti2(0.5,0.5,0.5) (我们这里将两个Ti原子区分为Ti1和Ti2是为了表明在极化子情形下其不同的电荷环境)，O(0.8,0.8,0.5),O(0.2,0.2,0.5),O(0.3,0.7,0),O(0.7,0.3,0)；

代码

文本

1.2 VESTA生成TiO2晶体

代码

文本

VESTA -> File -> New Structure:

a) 在"Phase"-"Title"栏中输入TiO2；

b) 在"Unit cell"栏中选择"System"="Tetragonal","Space Group"="No. 136 P42/mnm"，在“Lattice parameters"输入"a=4.653 Angstrom，c=2.969 Angstrom“；

c) 在"Structure parameters"栏中输入上述Ti原子和O原子的坐标即可。

代码

文本

通过VESTA可视化，可以初步了解金红石二氧化钛中Ti原子和O原子的大致位置

alt

代码

文本

2. TiO2结构优化(relax)

代码

文本

了解完TiO2的基本结构之后，我们接下来使用ABACUS DFT+U来对晶体结构优化（这里我们使用的优化方法是relax，固定晶轴，仅优化原子位置）。

代码

文本

[ ]

mkdir /opt/abacus-develop/examples/dos/TiO2_dos //创建一个文件夹，用于实施DOS的计算（注意，"/opt/abacus-develop"取决于你自己安装ABACUS的位置）

代码

文本

[ ]

cd /opt/abacus-develop/examples/dos/TiO2_dos //进入到相应目录下

代码

文本

2.1 `STRU`文件中的输入参数设置

代码

文本

STRU文件如下：

代码

文本

[ ]

ATOMIC_SPECIES

Ti1 47.867 Ti_ONCV_PBE-1.0.upf//赝势（PBE泛函）

Ti2 47.867 Ti_ONCV_PBE-1.0.upf

O 15.999 O_ONCV_PBE-1.0.upf

NUMERICAL_ORBITAL

Ti_gga_8au_100Ry_4s2p2d1f.orb//数值原子轨道基组，DZP（轨道和赝势都匹配地存放在"../../../test/PP_ORB“相对路径下，见INPUT文件pseudo_dir, orbital_dir）

Ti_gga_8au_100Ry_4s2p2d1f.orb

O_gga_7au_100Ry_2s2p1d.orb

LATTICE_CONSTANT

8.662//晶格常数，单位bohr

LATTICE_VECTORS

1.00 0.00 0.00//晶轴a,b,c

0.00 1.00 0.00

0.00 0.00 0.64

ATOMIC_POSITIONS

Direct

Ti1//原子类型名称

1.0//初猜原子磁矩1.0，对应Ti3+ (原子构型[Ar]3d1)

1//原子数1

0.00//x 坐标 0.00//y 坐标 0.00//z 坐标 1//移动x方向 1//移动y方向 1//移动z方向

Ti2//原子类型名称

0.0//初猜原子磁矩0.0，对应Ti4+ (原子构型[Ar])

0.50 0.50 0.50 1 1 1

0.0

0.79//打破高对称性 0.80 0.50 1 1 1

0.19//打破高对称性 0.20 0.50 1 1 1

0.29//打破高对称性 0.70 0.00 1 1 1

0.69//打破高对称性 0.30 0.00 1 1 1

代码

文本

注意到：与VESTA中设定的数值不同，我们故意把O原子的x坐标统一减去0.01,是为了打破晶体的高对称性，防止计算结果收敛到高对称性的局部极小值。

代码

文本

2.2 `INPUT`文件中的输入参数设置

代码

文本

INPUT文件如下：

代码

文本

[ ]

INPUT_PARAMETERS

suffix relax//输出文件的后缀（计算结果储存文件夹名称为：OUT.$suffix）

nelec_delta 1//额外加入一个电子， TiO2(-)

nbands 40//能带数

//注意：这是赝势和轨道存放的相对路径

pseudo_dir ../../../test/PP_ORB

orbital_dir ../../../test/PP_ORB

calculation relax//结构优化，固定晶轴长度

relax_nmax 500//最大离子步数500

force_thr_ev 0.001//原子受力收敛阈值，默认值是0.0257（lcao基本建议0.01即可）

ecutwfc 100//波函数截断能量

scf_thr 1.0e-8//scf收敛阈值（减小电子步收敛阈值是为了配合离子步的收敛阈值

scf_nmax 200//scf最大步数200

smearing_method gaussian//高斯smearing

smearing_sigma 0.0037//smearing程度，单位是Rydberg，默认值0.01

mixing_type broyden//broyden混合

mixing_beta 0.2//混合参数0.2

mixing_restart 1e-3//drho小于1e-3时清除密度矩阵历史

mixing_dmr 1//混合密度矩阵

ks_solver genelpa//kohm-sham方程求解器

basis_type lcao//数值原子轨道基组

gamma_only 0

symmetry 0

nspin 2//两个自旋通道(对应于quantum chemistry里的unrestricted calculation)

#Parameter DFT+U

dft_plus_u 1

orbital_corr 2 2 -1//按照STRU文件中的原子顺序，2 2 对应Ti1和Ti2原子上有U， -1对应O原子上没有U

hubbard_u 4.0 4.0 0.0//Ti原子上的U值大小为4.0 ev

out_bandgap 1//输出能带间隙

out_chg 1//输出收敛电荷

out_mul 1//输出mulliken电荷

代码

文本

2.3 `KPT`文件中的输入参数设置

代码

文本

KPT文件如下：

代码

文本

[ ]

K_POINTS

0//自动生成k点采样

Gamma//以gamma点为中心

4 4 4 0.0 0.0 0.0//每个倒格矢上有4个k点采样，中心无平移

代码

文本

2.4 运行ABACUS

代码

文本

[ ]

OMP_NUM_THREADS=1 mpirun -np 16 abacus | tee relax.output //运行abacus，开始结构优化计算

代码

文本

3. TiO2自洽计算(scf)

代码

文本

在得到优化的TiO2结构后，我们作一步自洽计算

代码

文本

3.1 `STRU`文件中的输入参数设置

代码

文本

读取结构优化的原子坐标，STRU文件如下：

代码

文本

[ ]

ATOMIC_SPECIES

Ti1 47.867 Ti_ONCV_PBE-1.0.upf//赝势

Ti2 47.867 Ti_ONCV_PBE-1.0.upf

O 15.999 O_ONCV_PBE-1.0.upf

NUMERICAL_ORBITAL

Ti_gga_8au_100Ry_4s2p2d1f.orb//数值原子轨道基组，DZP（轨道和赝势都匹配地存放在"../../../test/PP_ORB“相对路径下，见INPUT文件pseudo_dir, orbital_dir）

Ti_gga_8au_100Ry_4s2p2d1f.orb

O_gga_7au_100Ry_2s2p1d.orb

LATTICE_CONSTANT

8.662//晶格常数，单位bohr

LATTICE_VECTORS

1.00 0.00 0.00//晶轴a,b,c

0.00 1.00 0.00

0.00 0.00 0.64

ATOMIC_POSITIONS

Direct

Ti1 #原子类型名称

1.0000 #初始磁矩

1 #原子数

1.0000000000 0.0000000000 0.0000000000 m 1 1 1 mag 0.8153//覆盖初始磁矩

Ti2 #原子类型名称

0.0000 #初始磁矩

1 #原子数

0.5000000000 0.5000000000 0.5000000000 m 1 1 1 mag 0.3655//覆盖初始磁矩

O #原子类型名称

0.0000 #初始磁矩

4 #原子数

0.7942045448 0.7942045448 0.5000000000 m 1 1 1 mag -0.0431//覆盖初始磁矩

0.2057954552 0.2057954552 0.5000000000 m 1 1 1 mag -0.0431//覆盖初始磁矩

0.3008771922 0.6991228078 0.0000000000 m 1 1 1 mag -0.0474//覆盖初始磁矩

0.6991228078 0.3008771922 0.0000000000 m 1 1 1 mag -0.0474//覆盖初始磁矩

代码

文本

3.2 `INPUT`文件中的输入参数设置

代码

文本

INPUT文件如下：

代码

文本

[ ]

INPUT_PARAMETERS

suffix dos

nelec_delta 1//额外加入一个电子， TiO2(-)

nbands 40//能带数

//注意：这是赝势和轨道存放的相对路径

pseudo_dir ../../../test/PP_ORB

orbital_dir ../../../test/PP_ORB

calculation scf//自洽计算

ecutwfc 100//波函数截断能量阈值

scf_thr 1.0e-6//scf收敛阈值

scf_nmax 200//scf最大步数200

smearing_method gaussian//高斯smearing

smearing_sigma 0.0037//smearing程度，单位是Rydberg，默认值0.01

mixing_type broyden//broyden混合

mixing_beta 0.2//混合参数0.2

mixing_restart 1e-3//drho小于1e-3时清除密度矩阵历史

mixing_dmr 1//混合密度矩阵

ks_solver genelpa//kohm-sham方程求解器

basis_type lcao//数值原子轨道基组

gamma_only 0

symmetry 0

nspin 2//两个自旋通道(对应于quantum chemistry里的unrestricted calculation)

#Parameter DFT+U

dft_plus_u 1

orbital_corr 2 2 -1//按照STRU文件中的原子顺序，2 2 对应Ti1和Ti2原子上有U， -1对应O原子上没有U

hubbard_u 4.0 4.0 0.0//Ti原子上的U值大小为4.0 ev

out_bandgap 1//输出能带间隙

out_chg 1//输出收敛电荷

out_mul 1//输出mulliken电荷

代码

文本

3.3 `KPT`文件中的输入参数设置

代码

文本

KPT文件如下：

代码

文本

[ ]

K_POINTS

0//自动生成k点采样

Gamma//以gamma点为中心

4 4 4 0.0 0.0 0.0//每个倒格矢上有4个k点采样，中心无平移

代码

文本

3.4 运行ABACUS

代码

文本

[ ]

OMP_NUM_THREADS=1 mpirun -np 16 abacus | tee scf.output //运行abacus，开始自洽计算

代码

文本

4. TiO2非自洽计算（输出能态密度）

代码

文本

使用和自洽计算相同的STRU文件

代码

文本

4.1 `INPUT`文件中的输入参数设置

代码

文本

INPUT文件如下：

代码

文本

[ ]

INPUT_PARAMETERS

suffix dos

nelec_delta 1//额外加入一个电子， TiO2(-)

nbands 40//能带数

//注意：这是赝势和轨道存放的相对路径

pseudo_dir ../../../test/PP_ORB

orbital_dir ../../../test/PP_ORB

calculation nscf//非自洽计算

ecutwfc 100//波函数截断能量阈值

scf_thr 1.0e-6//scf收敛阈值

scf_nmax 200//scf最大步数200

smearing_method gaussian//高斯smearing

smearing_sigma 0.0037//smearing程度，单位是Rydberg，默认值0.01

mixing_type broyden//broyden混合

mixing_beta 0.2//混合参数0.2

mixing_restart 1e-3//drho小于1e-3时清除密度矩阵历史

mixing_dmr 1//混合密度矩阵

ks_solver genelpa//kohm-sham方程求解器

basis_type lcao//数值原子轨道基组

gamma_only 0

symmetry 0

nspin 2//两个自旋通道(对应于quantum chemistry里的unrestricted calculation)

#Parameter DFT+U

dft_plus_u 1

orbital_corr 2 2 -1//按照STRU文件中的原子顺序，2 2 对应Ti1和Ti2原子上有U， -1对应O原子上没有U

hubbard_u 4.0 4.0 0.0//Ti原子上的U值大小为4.0 ev

init_chg file//读取“OUT.dos”文件中的scf计算结果作为初猜

out_dos 1//输出能态密度

out_bandgap 1//输出能带间隙

out_mul 1//输出mulliken电荷

代码

文本

4.2 `KPT`文件中的输入参数设置

代码

文本

KPT文件如下：

代码

文本

[ ]

K_POINTS

0//自动生成k点采样

Gamma//以gamma点为中心

6 6 6 0.0 0.0 0.0//每个倒格矢上有6个k点采样，中心无平移（适当增加非自洽计算的k点采样，得到质量更高的能态密度结果）

代码

文本

4.3 运行ABACUS

代码

文本

[ ]

OMP_NUM_THREADS=1 mpirun -np 16 abacus | tee nscf.output //运行abacus，开始非自洽计算

代码

文本

5. 画能态密度TDOS/PDOS

代码

文本

[ ]

使用abacus-plot工具画TDOS，PDOS

(具体使用和安装也可详见abacus（本地使用版本为v3.6.5）目录下的"tools/plot-tools/README.md")

代码

文本

5.1 安装abacus-plot

代码

文本

[ ]

cd /opt/abacus-develop/tools/plot-tools//首先进入到abacus（本地使用版本为v3.6.5）目录下的"tools/plot-tools/

代码

文本

[ ]

python3 setup.py install

代码

文本

5.2 画总能态密度TDOS

代码

文本

[ ]

cd /opt/abacus-develop/examples/dos/TiO2_dos/OUT.dos //进入到计算scf和nscf的输出文件夹下

代码

文本

创建TDOS的config.json文件

代码

文本

[ ]

cat config.json

{

"tdosfile": "TDOS",

"efermi": 0,//费米能级，在running_scf.log中搜索“EFERMI”

"energy_range": [-10,10],//x轴范围

"dos_range": [0,5],//y轴范围

"tdosfig": "tdos.png"//保存图片格式和名称

}

代码

文本

[ ]

abacus-plot -d //画TDOS

代码

文本

5.2.1 不同U值的TDOS

代码

文本

U=1ev

代码

文本

alt

代码

文本

U=4ev

代码

文本

alt

代码

文本

U=6ev

代码

文本

alt

代码

文本

U=8ev

代码

文本

alt

代码

文本

5.2.2 对比VASP结果

代码

文本

alt

代码

文本

VASP输入参数：PAW_PBE Ti_sv, PAW_PBE O, LDAU=2, U=4ev, J=0ev 对比结论：ABACUS U=1ev结果与VASP U-J=4ev结果最相似

代码

文本

5.3 画投影能态密度PDOS

代码

文本

创建PDOS的config.json文件

代码

文本

[ ]

cat config.json

{

"pdosfile": "PDOS",

"efermi": 0,//费米能级，在running_scf.log中搜索“EFERMI”

"energy_range": [-10,10],

"dos_range": [-5,5],

"species": {"Ti1": [2], "Ti2": [2],"O": [1]},//画Ti1原子上的d轨道态密度，Ti2原子上的d轨道态密度和O原子上的p轨道态密度(0=s, 1=p, 2=d)

"pdosfig": "pdos.png"

}

代码

文本

[ ]

abacus-plot -d -p//画PDOS

代码

文本

alt

代码

文本

6.参考文献

代码

文本

[1] Wang Y C , Jiang H .Constrained Density Functional Theory Plus the Hubbard $U$ Correction Approach for the Electronic Polaron Mobility: A Case Study of TiO2[J]. Chinese J. Chem. Phys., 2021.https://doi.org/10.1007/s44214-023-00041-1

[2]司聪慧,张忠华,卢启芳,等.VESTA软件的晶体学可视化实验教学设计与应用[J].实验室研究与探索, 2023, 42(7):182-187.

[3] 赵天琦. 快速开始 ABACUS｜自洽能带态密度结构优化
[4] 郑大也. ABACUS磁性材料计算使用教程

代码

文本

notebook

点个赞吧