二维六角晶格的声子谱

上篇文章中介绍的是通过第一性原理使用VASP和phonopy计算声子谱，而为了更加深刻理解声子，在此篇我将会以二维简单六角晶格为例，考虑最近邻的力常数矩阵，从而推导出动力学矩阵的解析表达式，进而通过求解本征值本征态而得到六角晶格的声子谱以及振动模式。

理论基础

在简谐近似下，晶体中原子的原子的运动方程可以写成

\[\begin{equation}\label{eq1} D(\boldsymbol{k}) \epsilon(\boldsymbol{k},\sigma) = \omega^{2}_{\boldsymbol{k},\sigma} \epsilon(\boldsymbol{k},\sigma) \end{equation}\]

其中 \(D(\boldsymbol{k})\) 是动力学矩阵，矩阵中的元素是

\[\begin{equation}\label{eq2} D_{\alpha,\alpha^{'}}(\boldsymbol{k}) = \sum_{\boldsymbol{R_{t'}-\boldsymbol{R_{t}}}} \frac{\Phi_{t \alpha,t' \alpha'}}{\sqrt{m_{\alpha} m_{\alpha'}}} e^{i (\boldsymbol{R_{t'}-\boldsymbol{R_{t}}}) \cdot \boldsymbol{k}} \end{equation}\]

\(D_{\alpha,\alpha^{'}}(\boldsymbol{k})\) 一个是 \(d \times d\) 的矩阵 (对于二维体系 \(d=2\)，对于三维体系 \(d=3\))。 \(\Phi_{t \alpha,t' \alpha'}\) 是力常数矩阵，且在简谐近似下和弹性系数矩阵 \(K_{t \alpha,t' \alpha'}\) 的关系是

\[\begin{equation}\label{eq3} \Phi_{t \alpha,t' \alpha'} = K_{t \alpha,t' \alpha'} \end{equation}\]

它们是第 \(t\) 个原胞中的第 \(\alpha\) 个原子相对于第 $t' $ 个原胞中的第 \(\alpha’\) 个原子。

沿着一个原子和它的 \(n\) 次近邻的原子的弹性系数矩阵 有着如下的形式

\[\begin{equation}\label{eq4} K^{(n)} = \begin{pmatrix} K_{L}^{(n)} & 0 & 0 \\ 0 & K_{T}^{(n)} & 0 \\ 0 & 0 & K_{Z}^{(n)} \\ \end{pmatrix} \end{equation}\]

通过旋转操作，我们可以得到第 \(n\) 次近邻之间所有的弹性系数矩阵。而为了得到动力学矩阵，我们需要得到所有原胞之间的相互作用之和。

在 \(x-y\) 平面内沿着 以\(z\) 正方向为旋转轴的转动矩阵是

\[\begin{equation}\label{eq5} R(\theta) = \begin{pmatrix} \cos\theta & -\sin\theta & 0 \\ \sin\theta & \cos\theta & 0 \\ 0 & 0 & 1 \end{pmatrix} \end{equation}\]

在旋转操作下 \(R(\theta) \boldsymbol{r} = \boldsymbol{r'}\) ，那么可以得到旋转之后的弹性系数矩阵

\[\begin{equation}\label{eq6} \begin{aligned} K^{(n)}(\boldsymbol{r'}) &= R(\theta) K^{(n)} (\boldsymbol{r}) R^{-1}(\theta) = R(\theta) K^{(n)}(\boldsymbol{r}) R^{\dagger}(\theta) \\ &= \begin{pmatrix} \cos\theta & -\sin\theta & 0 \\ \sin\theta & \cos\theta & 0 \\ 0 & 0 & 1 \end{pmatrix} \begin{pmatrix} K_{L}^{(n)} & 0 & 0 \\ 0 & K_{T}^{(n)} & 0 \\ 0 & 0 & K_{Z}^{(n)} \\ \end{pmatrix} \begin{pmatrix} \cos\theta & \sin\theta & 0 \\ -\sin\theta & \cos\theta & 0 \\ 0 & 0 & 1 \end{pmatrix} \\ &= \begin{pmatrix} K_{L}^{(n)}\cos^{2}\theta+K_{T}^{(n)}\sin^{2}\theta & (K_{L}^{(n)}-K_{T}^{(n)})\sin\theta \cos\theta & 0 \\ (K_{L}^{(n)}-K_{T}^{(n)})\sin\theta \cos\theta & K_{L}^{(n)}\sin^{2}\theta+K_{T}^{(n)}\cos^{2}\theta & 0 \\ 0 & 0 & K_{Z}^{(n)} \end{pmatrix} \end{aligned} \end{equation}\]

六角晶格结构

六角晶格的结构和第一布里渊区如下图所示

其中晶格基矢为

\[\begin{equation}\label{eq7} \begin{aligned} \boldsymbol{a}_{1} &= (a,0) \\ \boldsymbol{a}_{2} &= (\frac{a}{2},\frac{\sqrt{3}a}{2}) \end{aligned} \end{equation}\]

和最近邻的矢量是

\[\begin{equation}\label{eq8} \begin{aligned} \boldsymbol{R}_{1} &= \frac{a}{\sqrt{3}}(0,-1) \\ \boldsymbol{R}_{2} &= \frac{a}{2\sqrt{3}}(-\sqrt{3},1) \\ \boldsymbol{R}_{3} &= \frac{a}{2\sqrt{3}}(\sqrt{3},1) \end{aligned} \end{equation}\]

倒空间基矢为

\[\begin{equation}\label{eq9} \begin{aligned} \boldsymbol{b}_{1} &= \frac{2\pi}{\sqrt{3}a}(\sqrt{3},-1) \\ \boldsymbol{b}_{2} &= \frac{4\pi}{\sqrt{3}a}(0,1) \end{aligned} \end{equation}\]

第一布里渊区的高对称点可以定义为

\[\begin{equation}\label{eq10} \begin{aligned} \boldsymbol{\mathrm{K}} &= \frac{4\pi}{3a}(1,0) \\ \boldsymbol{\mathrm{K'}} &= \frac{2\pi}{3a}(1,\sqrt{3}) \\ \boldsymbol{\mathrm{M}} &= \frac{\pi}{\sqrt{3}a}(\sqrt{3},1) \end{aligned} \end{equation}\]

动力学矩阵的构造

在二维的情况下，也就是 \(d=2\)时，沿着键方向（可以理解为未发生转动之前的）最近邻的弹性系数矩阵

\[\begin{equation}\label{eq11} R(\theta_i) = \begin{pmatrix} K_{L}^{(1)} & 0 \\ 0 & K_{T}^{(1)} ​ \end{pmatrix} \end{equation}\]

转动矩阵可以写成

\[\begin{equation}\label{eq12} R(\theta_i) = \begin{pmatrix} \cos\theta_i & -\sin\theta_i \\ \sin\theta_i & \cos\theta_i ​ \end{pmatrix} \end{equation}\]

那么可以转动得到转动后的弹性系数矩阵

\[\begin{equation}\label{eq13} \begin{aligned} K^{(1)}(\boldsymbol{R}_i) &= R(\theta_i) K^{(1)} R^{-1}(\theta_i) \\ &= \begin{pmatrix} K_{L}^{(1)}\cos^{2}\theta_i+K_{T}^{(1)}\sin^{2}\theta_i & (K_{L}^{(1)}-K_{T}^{(1)})\sin\theta_i \cos\theta_i \\ (K_{L}^{(1)}-K_{T}^{(1)})\sin\theta_i \cos\theta_i & K_{L}^{(1)}\sin^{2}\theta_i+K_{T}^{(1)}\cos^{2}\theta_i \end{pmatrix} \end{aligned} \end{equation}\]

那么对于三个最近邻的矢量可以写成

\[\begin{equation}\label{eq14} \begin{aligned} K^{(1)}(\boldsymbol{R}_1) &= K^{(1)} R^{-1}(\frac{3\pi}{2}) \\ &= \begin{pmatrix} K_{T}^{(1)} & 0 \\ 0 & K_{L}^{(1)} \end{pmatrix} \\ K^{(1)}(\boldsymbol{R}_2) &= K^{(1)} R^{-1}(\frac{5\pi}{6}) \\ &= \begin{pmatrix} \frac{3}{4} K_{L}^{(1)}+\frac{1}{4} K_{T}^{(1)} & -\frac{\sqrt{3}}{4} (K_{L}^{(1)}-K_{T}^{(1)}) \\ -\frac{\sqrt{3}}{4} (K_{L}^{(1)}-K_{T}^{(1)}) & \frac{1}{4} K_{L}^{(1)}+\frac{3}{4} K_{T}^{(1)} \end{pmatrix} \\ K^{(1)}(\boldsymbol{R}_3) &= K^{(1)} R^{-1}(\frac{\pi}{6}) \\ &= \begin{pmatrix} \frac{3}{4} K_{L}^{(1)}+\frac{1}{4} K_{T}^{(1)} & \frac{\sqrt{3}}{4} (K_{L}^{(1)}-K_{T}^{(1)}) \\ \frac{\sqrt{3}}{4} (K_{L}^{(1)}-K_{T}^{(1)}) & \frac{1}{4} K_{L}^{(1)}+\frac{3}{4} K_{T}^{(1)} \end{pmatrix} \end{aligned} \end{equation}\]

非对角块

对于非对角元部分我们可以写成

\[\begin{equation}\label{eq15} \begin{aligned} D_{A,B} &= -\sum_{i=1}^{3} \frac{K^{(1)}(\boldsymbol{R}_i)}{\sqrt{m_{A} m_{B}}} e^{i \boldsymbol{R}_{i} \cdot \boldsymbol{k}} \\ &= -\frac{1}{\sqrt{m_{A} m_{B}}} \left[ K^{(1)}(\boldsymbol{R}_1) e^{i \boldsymbol{R}_{1} \cdot \boldsymbol{k}} + K^{(1)}(\boldsymbol{R}_2) e^{i \boldsymbol{R}_{2} \cdot \boldsymbol{k}}+ K^{(1)}(\boldsymbol{R}_3) e^{i \boldsymbol{R}_{3} \cdot \boldsymbol{k}} \right] \\ &= -\frac{1}{\sqrt{m_{A} m_{B}}} \begin{pmatrix} d_{xx} & d_{xy} \\ d_{yx} & d_{yy} \\ \end{pmatrix} \end{aligned} \end{equation}\]

其中

\[\begin{equation}\label{eq16} \begin{aligned} d_{xx} &= K_{T}^{(1)} e^{i \boldsymbol{R}_{1} \cdot \boldsymbol{k}} + (\frac{3}{4} K_{L}^{(1)}+\frac{1}{4} K_{T}^{(1)}) (e^{i \boldsymbol{R}_{2} \cdot \boldsymbol{k}} +e^{i \boldsymbol{R}_{3} \cdot \boldsymbol{k}} ) \\ &= -\frac{1}{4} (K_{L}^{(1)} - K_{T}^{(1)} )(2 e^{i \boldsymbol{R}_{1} \cdot \boldsymbol{k}} - e^{i \boldsymbol{R}_{2} \cdot \boldsymbol{k}} - e^{i \boldsymbol{R}_{3} \cdot \boldsymbol{k}} ) + \frac{1}{2} (K_{L}^{(1)}+ K_{T}^{(1)}) (e^{i \boldsymbol{R}_{1} \cdot \boldsymbol{k}} + e^{i \boldsymbol{R}_{2} \cdot \boldsymbol{k}} +e^{i \boldsymbol{R}_{3} \cdot \boldsymbol{k}} ) \\ &= -\frac{1}{4} h K_{-}^{(1)} + \frac{1}{2} f K_{+}^{(1)} \\ d_{xy} &= -\frac{\sqrt{3}}{4} (K_{L}^{(1)}-K_{T}^{(1)}) e^{i \boldsymbol{R}_{2} \cdot \boldsymbol{k}} + \frac{\sqrt{3}}{4} (K_{L}^{(1)}-K_{T}^{(1)}) e^{i \boldsymbol{R}_{3} \cdot \boldsymbol{k}} \\ &= -\frac{\sqrt{3}}{4} (K_{L}^{(1)}-K_{T}^{(1)}) ( e^{i \boldsymbol{R}_{2} \cdot \boldsymbol{k}} - e^{i \boldsymbol{R}_{3} \cdot \boldsymbol{k}} ) \\ &= -\frac{\sqrt{3}}{4} g K_{-}^{(1)} \\ d_{yx} &= d_{xy}\\ d_{yy} &= K_{L}^{(1)} e^{i \boldsymbol{R}_{1} \cdot \boldsymbol{k}} + (\frac{1}{4} K_{L}^{(1)}+\frac{3}{4} K_{T}^{(1)}) (e^{i \boldsymbol{R}_{2} \cdot \boldsymbol{k}} +e^{i \boldsymbol{R}_{3} \cdot \boldsymbol{k}} ) \\ &= \frac{1}{4} (K_{L}^{(1)} - K_{T}^{(1)} )(2 e^{i \boldsymbol{R}_{1} \cdot \boldsymbol{k}} - e^{i \boldsymbol{R}_{2} \cdot \boldsymbol{k}} - e^{i \boldsymbol{R}_{3} \cdot \boldsymbol{k}} ) + \frac{1}{2} (K_{T}^{(1)}+ K_{L}^{(1)}) (e^{i \boldsymbol{R}_{1} \cdot \boldsymbol{k}} + e^{i \boldsymbol{R}_{2} \cdot \boldsymbol{k}} +e^{i \boldsymbol{R}_{3} \cdot \boldsymbol{k}} ) \\ &= \frac{1}{4} h K_{-}^{(1)} + \frac{1}{2} f K_{+}^{(1)} \end{aligned} \end{equation}\]

其中相因子 \(f\), \(g\), \(h\) 是

\[\begin{equation}\label{eq17} \begin{aligned} f &= e^{i \boldsymbol{R}_{1} \cdot \boldsymbol{k}} + e^{i \boldsymbol{R}_{2} \cdot \boldsymbol{k}} + e^{i \boldsymbol{R}_{3} \cdot \boldsymbol{k}} \\ g &= e^{i \boldsymbol{R}_{2} \cdot \boldsymbol{k}} - e^{i \boldsymbol{R}_{3} \cdot \boldsymbol{k}} \\ h &= 2 e^{i \boldsymbol{R}_{1} \cdot \boldsymbol{k}} - e^{i \boldsymbol{R}_{2} \cdot \boldsymbol{k}} - e^{i \boldsymbol{R}_{3} \cdot \boldsymbol{k}} \end{aligned} \end{equation}\]

和

\[\begin{equation}\label{eq18} \begin{aligned} K_{+}^{(1)} &= K_{L}^{(1)} + K_{T}^{(1)} \\ K_{-}^{(1)} &= K_{L}^{(1)} - K_{T}^{(1)} \end{aligned} \end{equation}\]

这样我们就可以得到非对角块的矩阵

\[\begin{equation}\label{eq19} \begin{aligned} D_{A,B} &= -\frac{1}{\sqrt{m_{A} m_{B}}} \begin{pmatrix} -\frac{1}{4} h K_{-}^{(1)} + \frac{1}{2} f K_{+}^{(1)} & -\frac{\sqrt{3}}{4} g K_{-}^{(1)} \\ -\frac{\sqrt{3}}{4} g K_{-}^{(1)} & \frac{1}{4} h K_{-}^{(1)} + \frac{1}{2} f K_{+}^{(1)} \\ \end{pmatrix} \end{aligned} \end{equation}\]

而对于 \(D_{B,A}\) 则满足

\[\begin{equation}\label{eq20} D_{B,A} = D_{A,B}^{*} \end{equation}\]

在这里简单说明一下，\(D_{B,A}\) 的矩阵元中最近邻基矢变了 \(\boldsymbol{R'}_{i} = -\boldsymbol{R}_{i}\)，而弹性系数矩阵保持不变，这是因为

\[\begin{equation}\label{eq21} \begin{aligned} K^{(1)}(\boldsymbol{R'}_i) &= R(\theta_i + \pi) K^{(1)} R^{-1}(\theta_i + \pi) \\ &= \begin{pmatrix} \cos(\theta_i+\pi) & -\sin(\theta_i+\pi) \\ \sin(\theta_i+\pi) & \cos(\theta_i+\pi) \end{pmatrix} \begin{pmatrix} K_{L}^{(1)} & 0 \\ 0 & K_{T}^{(1)} \\ \end{pmatrix}\begin{pmatrix} \cos(\theta_i+\pi) & \sin(\theta_i+\pi) \\ -\sin(\theta_i+\pi) & \cos(\theta_i+\pi) \end{pmatrix} \notag\\ &= \begin{pmatrix} \cos\theta_i & -\sin\theta_i \\ \sin\theta_i & \cos\theta_i \end{pmatrix}\begin{pmatrix} K_{L}^{(1)} & 0 \\ 0 & K_{T}^{(1)} \\ \end{pmatrix}\begin{pmatrix} \cos\theta_i & \sin\theta_i \\ -\sin\theta_i & \cos\theta_i \end{pmatrix} \notag\\ &= K^{(1)}(\boldsymbol{R}_i) \end{aligned} \end{equation}\]

然后我们可以得到

\[\begin{equation}\label{eq22} \begin{aligned} D_{B,A} &= -\sum_{i=1}^{3} \frac{K^{(1)}(\boldsymbol{R'}_i)}{\sqrt{m_{A} m_{B}}} e^{i \boldsymbol{R'}_{i} \cdot \boldsymbol{k}} \\ &= -\frac{1}{\sqrt{m_{A} m_{B}}} \left[ K^{(1)}(\boldsymbol{R'}_1) e^{i \boldsymbol{R'}_{1} \cdot \boldsymbol{k}} + K^{(1)}(\boldsymbol{R'}_2) e^{i \boldsymbol{R'}_{2} \cdot \boldsymbol{k}}+ K^{(1)}(\boldsymbol{R'}_3) e^{i \boldsymbol{R'}_{3} \cdot \boldsymbol{k}} \right] \\ &= -\frac{1}{\sqrt{m_{A} m_{B}}} \left[ K^{(1)}(\boldsymbol{R}_1) e^{-i \boldsymbol{R}_{1} \cdot \boldsymbol{k}} + K^{(1)}(\boldsymbol{R}_2) e^{-i \boldsymbol{R}_{2} \cdot \boldsymbol{k}} + K^{(1)}(\boldsymbol{R}_3) e^{-i \boldsymbol{R}_{3} \cdot \boldsymbol{k}} \right] \\ &= D_{A,B}^{*} \end{aligned} \end{equation}\]

对角块

由于力常数满足下面的等式

\[\begin{equation}\label{eq23} \sum_{t' \alpha'} \Phi_{t\alpha,t' \alpha'} = 0 \Longrightarrow \sum_{t' \alpha'} K_{t\alpha,t' \alpha'} = 0 \end{equation}\]

所以对于 \(\ce{A}\) 原子可以得到

\[\begin{equation}\label{eq24} \begin{aligned} K^{(0)}(A) + \sum_{i=1}^{3} K^{(1)}(\boldsymbol{R}_i) &= 0 \\ \Longrightarrow K^{(0)}(A) + \begin{pmatrix} \frac{3}{2} K_{L}^{(1)}+\frac{3}{2} K_{T}^{(1)} & 0 \\ 0 & \frac{3}{2} K_{L}^{(1)}+\frac{3}{2} K_{T}^{(1)} \end{pmatrix} &= 0 \end{aligned} \end{equation}\]

可以得到在 \(\ce{A}\) 原子点位上

\[\begin{equation}\label{eq25} \begin{aligned} K^{(0)}(A) &= -\begin{pmatrix} \frac{3}{2} K_{L}^{(1)}+\frac{3}{2} K_{T}^{(1)} & 0 \\ 0 & \frac{3}{2} K_{L}^{(1)}+\frac{3}{2} K_{T}^{(1)} \end{pmatrix} \notag\\ &= - \frac{3}{2} \begin{pmatrix} K_{+}^{(1)} & 0 \\ 0 & K_{+}^{(1)} \end{pmatrix} \end{aligned} \end{equation}\]

所以动力学矩阵对角块 \(D_{A,A}\) 为

\[\begin{equation}\label{eq26} D_{A,A} = -\frac{K^{(0)}(A)}{m_{A}} = \frac{3}{2 m_{A}} \begin{pmatrix} K_{+}^{(1)} & 0 \\ 0 & K_{+}^{(1)} \end{pmatrix} \end{equation}\]

而对于 \(\ce{B}\) 原子而言，由于弹性系数矩阵未发生变化，同样可以得到

\[\begin{equation}\label{eq27} K^{(0)}(B) = K^{(0)}(A) = - \frac{3}{2} \begin{pmatrix} K_{+}^{(1)} & 0 \\ 0 & K_{+}^{(1)} \end{pmatrix} \end{equation}\]

所以动力学矩阵对角块 \(D_{B,B}\) 为

\[\begin{equation}\label{eq28} D_{B,B} = -\frac{K^{(0)}(B)}{m_{B}} = \frac{3}{2 m_{B}} \begin{pmatrix} K_{+}^{(1)} & 0 \\ 0 & K_{+}^{(1)} \end{pmatrix} \end{equation}\]

动力学矩阵

经过上面分别对于非对角块和对角块的分析，我们可以得到动力学矩阵的总表达式

\[\begin{equation}\label{eq29} \frac{1}{\sqrt{m_{A} m_{B}}} \begin{pmatrix} \sqrt{\frac{m_B}{m_A}} \frac{3}{2} K_{+}^{(1)} & 0 & \frac{1}{4} h K_{-}^{(1)} - \frac{1}{2} f K_{+}^{(1)} & \frac{\sqrt{3}}{4} g K_{-}^{(1)} \\ 0 & \sqrt{\frac{m_B}{m_A}}\frac{3}{2} K_{+}^{(1)} & \frac{\sqrt{3}}{4} g K_{-}^{(1)} & -\frac{1}{4} h K_{-}^{(1)} - \frac{1}{2} f K_{+}^{(1)} \\ \frac{1}{4} h^{*} K_{-}^{(1)} - \frac{1}{2} f^{*} K_{+}^{(1)} & \frac{\sqrt{3}}{4} g^{*} K_{-}^{(1)} & \sqrt{\frac{m_A}{m_B}} \frac{3}{2} K_{+}^{(1)} & 0 \\ \frac{\sqrt{3}}{4} g^{*} K_{-}^{(1)} & -\frac{1}{4} h^{*} K_{-}^{(1)} - \frac{1}{2} f^{*} K_{+}^{(1)} &0 & \sqrt{\frac{m_A}{m_B}} \frac{3}{2} K_{+}^{(1)} \\ \end{pmatrix} \end{equation}\]

其中相因子 \(f\), \(g\), \(h\) 的表达式在式子 \(\eqref{eq17}\)，而\(f^{*}\), \(g^{*}\), \(h^{*}\) 分别是它们的复共轭。\(K_{+}^{(1)}\) 和 \(K_{-}^{(1)}\) 的表达式在式子 \(\eqref{eq18}\)。

声子谱

在等到动力学矩阵之后，通过运动方程 \(\eqref{eq1}\) 得到久期方程求解，便可以知道每个k点的本征值和本征矢，也就得到了二维六角晶格的声子谱。

下面是是给定确定的参数后得到的二维六角晶格的声子谱，四个图变化的参数是 \(m_B\) 的质量（或者说是 \(m_A\) 和 \(m_B\) 的质量比）。