本文是以下两篇文章的阅读笔记

M. Shao, F. H. Da Jornada, C. Yang, J. Deslippe, and S. G. Louie, “Structure preserving parallel algorithms for solving the Bethe–Salpeter eigenvalue problem,” Linear Algebra and its Applications 488, 148–167 (2016)

C. Penke, A. Marek, C. Vorwerk, C. Draxl, and P. Benner, “High performance solution of skew-symmetric eigenvalue problems with applications in solving the Bethe-Salpeter eigenvalue problem,” Parallel Computing 96, 102639 (2020).

在Full BSE中

$$\left[\begin{array}{cc}A& B\\ B^{\ast }& A^{\ast }\end{array}\right]\left[\begin{array}{c}X\\ Y\end{array}\right]=\mathrm{\Omega }\left[\begin{array}{cc}I& 0\\ 0& -I\end{array}\right]\left[\begin{array}{c}X\\ Y\end{array}\right]$$

其中

$$\begin{array}{rl}{A}_{ai{\mathit{k}}_1,bj{\mathit{k}}_2}^{\mathrm{BSE}}=& (E_{a{\mathit{k}}_1}^{\mathrm{GW}}-E_{i{\mathit{k}}_1}^{\mathrm{GW}}){\delta }_{ab}{\delta }_{ij}{\delta }_{{\mathit{k}}_1{\mathit{k}}_2}+\alpha (ai{\mathit{k}}_1|V|jb{\mathit{k}}_2)-(j{\mathit{k}}_2,i{\mathit{k}}_1|W|a{\mathit{k}}_1,b{\mathit{k}}_2)\\ B_{ai{\mathit{k}}_1,bj{\mathit{k}}_2}^{\mathrm{BSE}}=& (ai{\mathit{k}}_1|V|bj{\mathit{k}}_2)-(b{\mathit{k}}_2,i{\mathit{k}}_1|W|a{\mathit{k}}_1,j{\mathit{k}}_2)\end{array}$$

可以验证，A是厄米矩阵，B是对称矩阵。对于复数Full BSE矩阵，文章提出了以下对角化步骤：

第一步，构造M

$$M=\left[\begin{array}{cc}\mathrm{\Re }(A+B)& \mathrm{\Im }(A-B)\\ -\mathrm{\Im }(A+B)& \mathrm{\Re }(A-B)\end{array}\right]$$

第二步，构造J

$$J=\left[\begin{array}{cc}0& I\\ -I& 0\end{array}\right]$$

第三步，Cholesky分解M

$$M=U^{T}U$$

第四步，构造反对称矩阵$C=UJ{U}^T$ 第五步，调用skew_eigenvectors对角化反对称矩阵，得到本征值和本征矢

$$Cz=\mathrm{\Omega }z$$

第六步，构造SQ

$$SQ=\frac{1}{\sqrt{2}}\left[\begin{array}{cc}I& -iI\\ -I& -iI\end{array}\right]$$

然后得到原BSE方程的归一化的本征右矢

$$v=SQLz{\mathrm{\Omega }}^{-1/2}$$

本征右矢的构成是

$$v=\left[\begin{array}{cc}X& Y^{\ast }\\ Y& X^{\ast }\end{array}\right],X^{†}X-Y^{†}Y=1$$

与之正交归一的本征左矢构成是

$$v_l=\left[\begin{array}{cc}X& -Y^{\ast }\\ -Y& X^{\ast }\end{array}\right],v_l^{†}v=I$$

光谱

在tda近似下，吸收谱是按下式计算的

$$\begin{array}{rl}{ϵ}_2(\omega )& \stackrel{\text{length}}{=}\sum _S\frac{4{\pi }^2}{N_{k}V}{|\sum _{ai\mathit{k}}⟨i\mathit{k}|\overrightarrow{r}|a\mathit{k}⟩X_{ai\mathit{k}}^S|}^2\delta (\omega -{\mathrm{\Omega }}_S)\\ & \stackrel{\text{vel\_div}}{=}\sum _S\frac{4{\pi }^2}{N_{k}V}{\left|\sum _{ai\mathit{k}}\frac{⟨i\mathit{k}|\overrightarrow{v}|a\mathit{k}⟩}{E_a-E_i}{X}_{ai\mathit{k}}^S\right|}^2\delta (\omega -{\mathrm{\Omega }}_S)\\ & \stackrel{\text{vel}}{=}\sum _S\frac{4{\pi }^2}{N_{k}V{\omega }^2}{|\sum _{ai\mathit{k}}⟨i\mathit{k}|\overrightarrow{v}|a\mathit{k}⟩X_{ai\mathit{k}}^S|}^2\delta (\omega -{\mathrm{\Omega }}_S)\end{array}$$

根据矩阵的谱分解，full应该是按下式推广为

$${ϵ}_2(\omega )=\sum _S\frac{4{\pi }^2}{N_{k}V{\omega }^2}(\sum _{vck}⟨v\mathbf{k}|\overrightarrow{v}|c\mathbf{k}⟩X_{vck}^S+⟨c\mathbf{k}|\overrightarrow{v}|v\mathbf{k}⟩Y_{vck}^S){(\sum _{vck}⟨v\mathbf{k}|\overrightarrow{v}|c\mathbf{k}⟩X_{vck}^S-⟨c\mathbf{k}|\overrightarrow{v}|v\mathbf{k}⟩Y_{vck}^S)}^{\ast }\delta (\omega -{\mathrm{\Omega }}_n)$$

但是经过测试，上式和tda的结果相差很大。目前程序中使用的还是

$${ϵ}_2(\omega )=\sum _S\frac{4{\pi }^2}{N_{k}V{\omega }^2}{|\sum _{vck}⟨v\mathbf{k}|\overrightarrow{v}|c\mathbf{k}⟩X_{vck}^S+⟨c\mathbf{k}|\overrightarrow{v}|v\mathbf{k}⟩Y_{vck}^S|}^2\delta (\omega -{\mathrm{\Omega }}_n)$$