紧束缚模型与Wannier函数介绍
讨论固体物理,非常重要的条件就是晶格的周期性。在讨论后续的模型之前,我们先已经知道,周期性的哈密顿量得到的波函数应当为布洛赫波(Bloch waves),也是近自由电子模型、紧束缚模型共同的出发点:
ψkn(r)=eik⋅rukn(r)
其中k在第一布里渊区中,n是能带的标记,ukn(r)拥有晶格周期性:
ukn(r+任意晶格矢量R)=ukn(r).
可见布洛赫波弥散与整个晶体的空间,每个布洛赫波都有一个确定的参数k。如果对于一个高迁移率(high mobility)导电电子的体系,我们常用的近似叫做近自由电子模型,直接认为导电电子所感受到的背景电势(晶体赝势)可以忽略。电子波函数就是平面波,ukn=1。
但是对于晶格势能严重影响导电电子的体系,比如过度金属氧化物,我们用布洛赫波发展出了另一种从“局域电子”出发的模型:紧束缚模型(tight-binding systems)。
其实很容易理解:平面波弥散于整个空间,但是进行傅里叶变换,也就是不同的平面波积分,我们能够得到δ函数(局域化);类似地,布洛赫波近似于平面波弥散在整个空间,如果我们对其进行傅里叶变换,应该得到一个虽然不是完美δ函数,但也相当局域的函数。这个变换的结果就是Wannier瓦涅尔函数。可以容易想象,这种实质是傅里叶变换的纯数学操作自然不会丢失任何信息,也就是新的函数(Wannier函数)的正交性、完备性和变换前的布洛赫波一样。
数学的语言描述这种变换为:
(1)
∣ψRn⟩≡N1k∑B.Z.e−ik⋅R∣ψkn⟩
其中∣ψRn⟩就是Wannier态,下标从布洛赫态的k变成R,因为是傅里叶变换嘛,很好理解。
反过来的变换为:
(2)
∣ψkn⟩≡N1R∑eik⋅R∣ψRn⟩
需要注意,经过傅里叶变换,Wannier态在k空间上是弥散的(因而在坐标空间局域),而且Wannier态不再是哈密顿量的本征态(布洛赫波是,相加后就不是了)。
还有一个可以想象的极限状态:如果我们的晶体中原子距离非常遥远以至于可以看作孤立原子,那么局域的Wannier函数应当会退化为原子中的电子波函数。
上面说了个概念。下面的内容我们要玩真的,来看Wannier函数下,我们如何计算体系的哈密顿量。数学的手段是二次量子化的操作,不太清楚的可以直接看Condensed Matter Field Theory(2nd edition)的第二章第一节。
晶体中单电子的Hamilton量为:
H^0=2mp^2+V(r),
在二次量子化的表象下,波函数成为了粒子数升降算符(aσ†(r)),也要进入哈密顿量中:
(3)
H^0=∫draσ†(r)[2mp^2+V(r)]aσ(r),
其中σ代表电子的自旋,此文用不到电子自旋,但是为了理论的完备,还是写出来。
在实空间中求哈密顿量本征态,得到的就是上面所示的布洛赫态∣ψkn⟩,按照二次量子化,我们定义产生算符akσ†(是aσ†(r)的一个实例)使得
∣ψkn⟩=akσ†∣0⟩,
同时方便定义产生算符的共轭算符:湮灭算符akσ。在这种二次量子化的表象下,Hamilton量被写为(所有单个电子求和)
(4)
H^0=k∑ϵkakσ†akσ.
这里已经用到了性质:布洛赫态是该哈密顿量本征态,因此有确定的能量ϵk。
下面就是表象的变换,从布洛赫转换为Wannier,表象变换操作是大家在量子力学习题中应该见过多次的常用方法(运用完备性I=∣ψRn⟩⟨ψRn∣),下面直接给出结果,不清楚的看Condensed Matter Field Theory书公式(2.23)前后:
(5)
akσ†=N1R∑eik⋅Riaiσ†,aiσ†=N1k∑B.Z.e−ik⋅Riakσ†
其中aiσ†就是我们需要的Wannier态的产生算符∣ψRin⟩=aiσ†∣0⟩,i下标其实是简写了的Ri。将(5)代入(4)完成变换:
(6)
H^0=N1ii′∑k∑eik⋅(Ri−Ri′)ϵkaiσ†ai′σ≡ii′∑aiσ†tii′ai′σ,
其中tii′=N−1∑ii′∑keik⋅(Ri−Ri′)ϵk。这就是目前不考虑其它相互作用,紧束缚模型得到的哈密顿量。
在我们使用紧束缚模型时,常常用到的近似是:tii′仅存在于最近邻原子Ri,Ri′之间(因为Wannier函数局域,忽略更远近邻的影响)。
简单讨论一下紧束缚模型(6)公式的含义:
-
对于i=i′也就对同一位置的电子,哈密顿量的对角位,tii=N∑kϵk≡ϵ,理解上大致对应一个局域电子本身的能量;
-
i=i′最近邻电子之间,哈密顿量的非对角位,计算得到的tii′=N−1∑ii′∑keik⋅(Ri−Ri′)ϵk一般为一个负值(在紧束缚近似下,tii′=−t,t>0),理解为局域的电子向其最近邻位置有一定跳跃的可能,这种跳跃增加了局域电子的活动空间,从而降低了系统的能量。
这里大家可能有点迷糊:tii′=N−1∑ii′∑keik⋅(Ri−Ri′)ϵk这个公式看起来非常的不直观,不知道它代表了什么样的物理意义,怎么就莫名奇妙约等于一个负的常数−t了呢?如果我们换一种写法,不是用(4)的从布洛赫态展开哈密顿量,而是用(3)的更为一般的实坐标态下的哈密顿量写,就会更清楚。我们使用aσ†(r)→aiσ†的表象变换(推导过程与(5)原理一致):
aσ†(r)=i∑ψRi∗(r)aiσ†
将(3)变为:
H^0=ii′∑∫drψRi∗(r)aiσ†[2mp^2+V(r)]ψRi′(r)ai′σ,
也就是:tii′=⟨ψRi∣H^0∣∣∣ψRi′⟩=ti′i∗。这种量子力学中的矩阵元大家应该熟悉,体现的是最近邻两个Wannier函数的重合部分对能量的贡献,也应该理解为什么我们自然而然将它近似为常数。
紧束缚近似写为:
(7)
tii′=⎩⎪⎨⎪⎧ϵ−t0i=i′ii′ nearest neighboursotherwise
提醒一下,到这里我们讨论的全是哈密顿量H^0,本质是单电子哈密顿量,而没有考虑电子-电子的库伦相互作用。讨论相互作用的影响的紧束缚模型不在本文讨论范围内。
石墨烯结构

This figure is generated by TikZ/LaTeX(From郑奇靖博文)
石墨烯的蜂窝结构可谓非常经典了。只是需要说明下它的晶胞:原胞呈现菱形(灰色虚线),每个原胞内有两个碳原子A和B(红,蓝);每个原胞内有三条不同的最近邻连接线。a1,a2为两个原胞矢量。
晶格矢量R(x,y)=xa1+ya2,其中x,y为整数。
这种结构的成因:在石墨烯的平面内,sp2电子杂化,三条夹角120°的成键方向,产生每个碳原子连接三个碳原子的蜂窝结构。但是还剩下一个pz轨道垂直于平面,各个原子的pz轨道电子互相有一点的影响,产生π键。而我们的紧束缚模型,正好可以很好的用于解释这个π键电子的行为!
每个原胞中三条不同的最近邻线,以图中标注原子A起始为例,假设AB原子所在的晶胞位置为R(0,0):
- 向右连接本晶胞内R(0,0)的B,
- 向左上连接晶胞R(0,−1)中的B,
- 向左下连接晶胞R(−1,0)中的B。
这三个最近邻的作用应当是完全相同的,因为晶格的对称性。(旋转120°,240°完全一样)
其实你或者也可以从原子B起始,也写出3条最近邻线,这其实就是上面三条线晶格平移的结果:
- 向左连接本晶胞内R(0,0)的A,
- 向右下连接晶胞R(0,1)中的A,
- 向右上连接晶胞R(1,0)中的A。
从A开始或从B开始,二者等价。
紧束缚近似应用于石墨烯𝜋键
我们将(6)(7)应用于刚刚介绍的石墨烯,也就是计算石墨烯体系的能量。更为确切地说,是计算石墨烯(紧束缚近似下)的能带,即能量E关于动量空间k的函数。我们知道,其实Wannier态其实在动量空间是弥散的,所以最后的结果,我们不应该保留局域的Wannier态,而是仍然要傅里叶变换回布洛赫波的表象(拥有确定的k可以作为自变量)。
注意一个小问题,石墨烯的原胞内有两个碳原子,需要区分讨论,因此(7)中本来只有一个的aiσ†这里也要分化出两个ai1σ†和ai2σ†,下标1,2分别代表同一个原胞i中的原子A和B。两者对应的Wannier函数(各自的局域电子)自然也要求互相正交。(也就是说,对于每个晶格,并非只有一个Wannier态,而是2个原子分别有一个)Wannier函数变多了一倍,布洛赫波的个数自然也要扩大一倍:ak1σ†和ak2σ†。
由于这个变化,我们的哈密顿量的写法相对(6)要有小小的修改,看看ai1σ†和ai2σ†怎么加入进去。整个原胞内的哈密顿量有2个原子+3条最近邻组合。首先看2个原子上的电子本身都拥有能量,对应着哈密顿量的对角项,也就是:
(8)
H^diag=i∑ϵ(ai1σ†ai1σ+ai2σ†ai2σ)=ϵi∑kk′∑N1e−ik⋅Riak1σ†eik′⋅Riak′1σ+(1→2)=ϵkk′∑i∑N1ei(k′−k)⋅Riak1σ†ak′1σ+(1→2)=ϵkk′∑δ(k′−k)ak1σ†ak′1σ+(1→2)=ϵk∑ak1σ†ak1σ+(1→2)
其中第二行用到了(5)Wannier态和布洛赫态的转换;第三行用了∑iN1ei(k′−k)⋅Ri=δ(k′−k)。
剩下的还有3个最近邻作用,也就是最近邻hopping使得能量降低,对应的是哈密顿量中的非对角项。需要说明的是,虽说是3个最近邻作用,但是书写时我们经常把同一个作用的两写法A->B,与B->A都写出来,二者各乘0.5再加起来,这样写最为对称简洁。所以3个作用我们会写出对称的6项。
(9)
H^offdiag=====i∑−t(ai1σ†ai2σ+ai1σ†a(Ri+R(0,−1))2σ+ai1σ†a(Ri+R(−1,0))2σ)+i∑−t(ai2σ†ai1σ+ai2σ†a(Ri+R(0,1))1σ+ai2σ†a(Ri+R(1,0))1σ)−ti∑kk′∑N1[ei(k′−k)⋅Ri+ei[k′⋅(Ri−a2)−k⋅Ri]+ei[k′⋅(Ri−a1)−k⋅Ri]]ak1σ†ak′2σ+h.c.−tkk′∑i∑N1ei(k′−k)⋅Ri(1+e−ik′⋅a2+e−ik′⋅a1)ak1σ†ak′2σ+h.c.−tkk′∑δ(k′−k)(1+e−ik′⋅a2+e−ik′⋅a1)ak1σ†ak′2σ+h.c.−tk∑(1+e−ik⋅a2+e−ik⋅a1)ak1σ†ak2σ+h.c.
其中h.c.代表了Hermitian conjugate(厄米共轭)。
综合(8)(9),我们可以把总哈密顿量对于每一个k拆开写成H^(k),并写成矩阵形式:
(10)
H^(k)=[ak1σ†ak2σ†][ϵ−t(1+eik⋅a2+eik⋅a1)−t(1+e−ik⋅a2+e−ik⋅a1)ϵ][ak1σak2σ]
对其进行对角化,而本征值就是能量。求本征值的操作过于基础,此处忽略。直接得到结果:
E(k)=ϵ±t1+e−ik⋅a1+eik⋅a1+e−ik⋅a2+eik⋅a2+(e−ik⋅a1+e−ik⋅a2)(eik⋅a1+eik⋅a2)=ϵ±t3+2cos(k⋅a1)+2cos(k⋅a2)+2cos(k⋅(a1−a2)).
再为了方便后面绘图,代入xy坐标下的原胞矢量:
a1a2=2a(3,3)=2a(3,−3),
其中a就是晶格常数。我们就可以把k也写成分量,则有:
(11)
E(kx,ky)=±t3+2cos(3kya)+4cos(23kxa)cos(23kya).
我们忽略了能量基准ϵ。
下面便是使用Plotly对能带的绘图: