130 lines
10 KiB
Markdown
130 lines
10 KiB
Markdown
---
|
||
title: '使用 VASP 计算声子'
|
||
date: 2025-02-26T10:46:07+08:00
|
||
draft: true
|
||
summary: some summary
|
||
---
|
||
|
||
通常来说，计算声子的原理是这样的[^2]：
|
||
1. 计算力矩阵。
|
||
“力矩阵”是一个实数矩阵，长宽都是三倍的原子数目（因为每个原子有三个方向可以活动）。
|
||
它的每一个元素都描述这样一个事情：一个原子移动一小段距离（简谐近似），另外一个原子会感受到多大的力？
|
||
或者也可以说，每个元素都是系统能量对原子位置的二阶导数。
|
||
显然它应该是一个对称矩阵。
|
||
2. 计算它的特征值和特征向量。
|
||
得到的特征值就是声子的频率（可能要乘以一个与质量有关的系数，记不清了），
|
||
特征向量就是这个声子模式中各个原子的振动幅度（同样要乘以一个与质量有关的系数）。
|
||
3. 好了，算完了。就这样。之后也可以再去计算声子态密度或者拉曼张量等，但这些都是后面的事情了。
|
||
|
||
首先，为什么是这个矩阵的特征值？我是这样理解的：
|
||
* 要精确描述系统的振动，最好是假定这个振动和系统本身都满足某种周期性，这样一个周期内原子个数是有限个，就可以用有限个量来描述这个振动。
|
||
使用声子的语言描述振动时就是这样做的，我们取的“周期”就是一个晶胞，原子排列满足这个周期性；
|
||
振动状态（原子在振动过程中的位移）不一定满足这个周期（想象一下相邻晶格振动方向相反的情况，这种情况我们也需要研究），
|
||
但经过“调整”后也满足这个周期（类似于 bloch 函数那样，乘上某一个波矢的平面波）。
|
||
* 在一个周期内的原子的运动状态组成的线性空间里，我们可以取不同的基来描述它。
|
||
力矩阵的特征向量就是满足条件的一组基，因此将原子的振动状态总是理解为力矩阵的特征向量的线性组合是可行的。
|
||
* 之所以使用这组特征向量（而不是别的某个矩阵的特征向量），是因为这组特征向量有这样一个特殊的性质：
|
||
在简谐近似的前提下考虑，其中任意一个向量所代表的振动状态是可以稳定存在的。
|
||
也就是说，如果系统现的振动状态对应其中的某个向量，那之后系统会一直是这个状态（一直对应于这个向量）。
|
||
而其它矩阵的特征向量（或者说，这一组特征向量的线性组合，只要不是恰好在同一个特征子空间里），都没有这个性质。
|
||
再换句话说，就是：这组特征向量描述的振动是比较“稳定”的（指对它的描述与时间无关），
|
||
只有需要考虑破坏了对称性或者非简谐效应的时候才会有变化（而这些东西的影响都相对来说比较小）。
|
||
|
||
一切看起来都十分美好：只要想办法算出来这个力矩阵，剩下的就都是一些数学处理了
|
||
（并且这个处理对计算机来说并不困难，随便一个家用计算机都可以在几秒或者几分钟内处理完成）。
|
||
唯一的问题就是如何计算这个力矩阵。
|
||
|
||
一个直观的方法是这样的：
|
||
把每个原子都向各个方向移动一小段距离，然后计算系统的能量或者各个原子的受力，然后就得到了矩阵的一行或者几行（如果考虑对称性的话）。
|
||
事实上也的确可以这样做，这就是所谓的“有限位移法”，
|
||
洋文叫“[finite differences](https://www.vasp.at/wiki/index.php/Phonons_from_finite_differences)”[^1]。
|
||
有一些小问题需要考虑，比如晶胞往往比较小、相邻晶胞中的同一个原子会相互影响，解决办法是在一个比较大的晶胞中计算；
|
||
比如施加的位移往往比较小（为了简谐近似），所以需要比较高精度的计算（不然这点受力就会被计算中的误差抹掉）；
|
||
比如这个方法只能直接得到Γ点的声子，解决办法还是扩胞（得到几个点的声子，剩下的插值）。
|
||
总之都是一些可以解决的小问题。
|
||
然后就可以计算得到想要的性质，和实验或者文献比对一下，绝大多数都符合，简直完美。
|
||
然后就可以去研究剩下那一点点对不上的东西，然后就发现：不对劲。
|
||
|
||
原理上来说，是这样的：扩胞的时候我们需要扩到足够大（一般边长为十几个埃），使得相邻晶胞中的同一个原子不会互相影响。
|
||
但有时并不能做到这样：有些原子可能带电，位移后可能导致了长程的电场（库伦相互作用）。
|
||
更具体地说，在非极性材料中，这种效应不会出现（我们只考虑简谐作用，也就是假定电场（受力）与位移成正比）；
|
||
在理想晶体中，虽然在局域（移动的原子周围）会产生电场，但很快就被电子云屏蔽掉。
|
||
当这两个条件都不满足时，例如在极性比较强的半导体中，这种效应就会显现出来，
|
||
这时，十几个埃的晶胞并不足以使得相邻晶胞中的原子不会互相影响。
|
||
|
||
以 Gamma 点的声子为例。
|
||
计算这时的 dynamical matrix 时，我们只需要把计算出来的原始数据中的某一行加起来，
|
||
也就是，将超胞中，某一个原子朝着某个方向移动后，处在各个晶胞中的另外一个原子在某个方向上受到的力，把这一列数据加起来，
|
||
就得到了 dynamical matrix 的一个元素。
|
||
如果不是 Gamma 点，需要做一个相位的变换，这个不表，
|
||
我们只需要考虑 Gamma 点的情况，就足以说明为什么没有消除这种远程效应后，会导致计算出来的声子频率不对了。
|
||
相加的过程中，可以预见，如果没有这样的长程效应的话，较远处的受力应该会很小；如果是手算的话，完全可以忽略这些数据。
|
||
倘若你发现这些数据并不小，那么就说明超胞取得太小了；你应该再扩大超胞，才能得到正确的结果，这时会有更多的数据加入计算。
|
||
但一旦出现这种长程效应，十几个埃的晶胞就不够用了，可能得扩大到几十个埃或者上百个埃，才能套用这个方法；
|
||
或者说，如果一定要用十几个埃的超胞来算的话，相加的时候就少加了几个数据[^3]。这就是为什么会出现计算出来的声子频率不对的情况。
|
||
|
||
解决办法是这样的：
|
||
把原子之间的受力拆成两个部分，一个是长程的库伦相互作用，受到电子云的屏蔽作用，大概是指数衰减[^4]；
|
||
另外一个是量子力学的作用，就像对交换关联能的估计那样，它是相当局域的，只在一个化学键那么长的范围内有作用。
|
||
我们只要考量算出来的原子受力随着第几个晶格的变化，就可以估计出这两个部分分别有多大，然后用不同的方案计算 dynamical matrix 就可以了。
|
||
|
||
VASP 已经包含了这个功能，只需要设置 `LPHON_POLAR` `PHON_DIELECTRIC` `PHON_BORN_CHARGES` 这几个参数就可以了，
|
||
具体怎么设置在 vasp wiki 上有介绍。
|
||
此外，vasp wiki 上还有一段[关于这个方法的原理的讨论](https://www.vasp.at/wiki/index.php/Phonons:_Theory#Long-range_interatomic_force_constants_.28LO-TO_splitting.29)，但我没看懂。
|
||
但我即使没看懂，给它捉了一个虫，说明这个页面可能真的是没什么人看。
|
||
|
||
其实有一个偷懒的解决办法：如果只要 Gamma 点的声子的话，可以只用单胞来计算，得到的结果反而是准确的。
|
||
因为算出来直接就是 Gamma 的 dynamical matrix，什么效应都考虑进去了。
|
||
当然这样就真的只能算 Gamma 点的声子[^5]。
|
||
|
||
想象这样一个情景：振动本身以某个周期在进行。
|
||
你可以把这个周期想象成 VASP 计算时的超胞，但它的边长可以大一些，比如有几百几千个晶胞。
|
||
每个周期中，不仅原子的初始位置相同，振动导致的位移也是相同的。
|
||
我们假定这个超胞中有 $N$ 个单胞，每个单胞中有 $n$ 个原子。
|
||
那么，就可以将各个原子在各个方向上的位移写成一个 $3nN$ 维的列向量，记为 $u$。
|
||
进一步，把它对时间的二阶导数记为 $\ddot{u}$。
|
||
|
||
在我们考虑中（简谐近似），哈密顿量长这个样子：
|
||
|
||
$$
|
||
H = \frac{1}{2} M \ddot{u}^2 + \frac{1}{2}\Phi u^2
|
||
$$
|
||
|
||
$M$ 是各个原子的质量，$\Phi$ 是力矩阵，定义为能量对位移的二阶导数：
|
||
|
||
$$
|
||
\Phi_{ij} = \frac{\partial^2 H}{\partial u_i \partial u_j}
|
||
$$
|
||
|
||
套入哈密顿力学的某个方程（忘记那个叫什么了，反正就是那个东西），可以得到：
|
||
|
||
$$
|
||
M \ddot{u} + \Phi u = 0
|
||
$$
|
||
|
||
如果你不习惯哈密顿力学的说法（其实在经典力学里，我也不习惯），那也可以把上面这个式子理解成 $F = ma$。
|
||
注意到原子受力是 $-\Phi u$，就也可以得到这个式子。
|
||
总之这就是我们要解的方程。
|
||
|
||
如上所说：要研究振动，最好是假定有一定的“周期性”。
|
||
因此我们取：
|
||
|
||
$$
|
||
u = \Re{u_0 e^{i \vec{q}\cdot\vec{r}}}
|
||
$$
|
||
|
||
高温
|
||
异质结构 位错 三个方面 n p 杂质 O 空位 对能带的影响
|
||
gan 射频器件 失配 应力 位错
|
||
意义 技术路线 存在问题 （加电磁波）
|
||
应力 温度 电磁波
|
||
|
||
[^1]: 写到这里我才意识到应该翻译成“有限差分”而不是“有限位移（displacement）”。
|
||
但是我感觉我也在哪里看到过“有限位移”的说法（也许是我记错了），但反正就是那个意思，就这样吧。
|
||
[^2]: 这并不是一个非常严格的描述，甚至可以说是错误的描述：我们甚至没有区分力矩阵和 dynamical matrix。
|
||
[^3]: 其实不仅是少算了几个距离较远的数据，同时也多算了几个距离较近的数据（隔壁晶胞中移动的那个原子导致的效应也被算进来了）。
|
||
也正是由于这些多算了的数据，导致之后偷懒的办法可以使用。
|
||
[^4]: 单个原子的作用是指数衰减，但如果考虑到，距离越远，施加作用的原子数量越多（正比于球面面积）的话，
|
||
这个作用大概是 $\frac{r^2}{e^r}$ 的衰减。
|
||
[^5]: 如果用的“单胞”其实是超胞的话，同时需要不能让 VASP 根据平移对称性自动寻找到更小的单胞，我没找到怎么做。
|