两带超导体边界效应及其对二硼化镁薄膜临界温度的影响研究

Boundary Effect of Two-Band Superconductors and Its Impact on the Critical Temperature of MgB₂ Thin Film

摘要: 基于两带Bogoliubov-de Gennes理论，我们针对两带超导体–绝缘体(或真空)之间的边界效应展开研究。在两带Ginzburg-Landau自由能表达式基础上引入新的表面项，得到其对超导序参量影响的特征尺度。在考虑边界效应的影响下，我们对二硼化镁超导薄膜临界温度随厚度的变化关系进行了理论计算，所得结果与实验数据相符。

Abstract: Based on two-band Bogoliubov-de Gennes theory, we study the boundary effect of an interface between two-band superconductor and insulator (or vacuum). New boundary terms are introduced into two-band Ginzburg-Landau free energy, and the characteristic length scale of boundary effect can be estimated. Taking into account this boundary effect, we study the critical temperature dependence of film thickness for magnesium diboride. Our numerical results are in good agreement with the experimental data measured in this material.

1. 引言

2001年日本科学家秋光纯等报道了MgB₂具有超导电性，其临界转变温度 $T_{c0}$ 高达40 K [1]，该发现迅速激起了全世界范围内关于该超导材料的研究热潮。MgB₂是由六方紧密排列的镁原子层和类石墨结构的硼原子层交替组成，其中硼原子对二硼化镁超导电性起着关键作用。研究发现MgB₂属于常规超导体，符合传统的电声耦合机制 [2]，其费米面是由两个二维的σ带和两个三维的π带组成 [3] [4]。随后，有关该超导体两个能隙的图像被比热、核磁共振、角分辨光电子谱和电子隧道等一系列实验广泛证实 [5] [6] [7] [8] [9]。

MgB₂超导薄膜对于电子、高场和频射腔等应用具有重要意义 [10] [11] [12]，而制备出高质量的超导薄膜是研究其物理特性及实际应用的重要前提。目前人们已对MgB₂超导薄膜的沉积技术进行了大量研究，并积累了非常丰富的制备经验，其中混合物理化学气相沉积法(HPCVD)是制备MgB₂超导薄膜较为有效的技术 [13] [14] [15]。例如，冯庆荣团队利用HPCVD技术在SiC衬底上制备出一系列厚度为8 μm到8 nm的MgB₂超导薄膜 [16] [17]，发现当薄膜厚度从30 nm减小到8 nm时，临界温度下降了13%。潘杰云等人在MgO衬底上生长出厚度在40 nm到10 nm的MgB₂薄膜，并研究了超导临界温度和厚度间的关系 [18]，发现在10 nm的薄膜中， $T_c$ 急剧下降到34 K。于是人们提出了各种不同的机制来解释MgB₂超导薄膜中临界温度被抑制的问题，如超导序参量的量子相位波动 [19] 和无序引发的库仑相互作用增强等 [20]。但是到目前为止，对于上述实验数据的解释仍然没有达成共识。

在本文中，我们将针对两带超导体–绝缘体(或真空)之间的边界效应展开研究，并试图探究两带超导体边界效应对MgB₂薄膜临界温度的影响。我们首先在超导体界面为超导序参量引入合适的边界条件，然后基于两带Bogoliubov-de Gennes理论得到了MgB₂薄膜表面项的特征尺度。然后我们根据两带Ginzburg-Landau (GL)理论计算了MgB₂薄膜临界温度与厚度之间的依赖关系，并将理论计算结果与实验数据进行比较，发现理论分析与实验数据一致。我们的工作表明边界效应是理解MgB₂薄膜超导电性的一个重要因素。

本文主要由四部分组成，在接下来的一部分中，我们回顾了两带Bogoliubov-de Gennes理论和GL方程，对自由能中引入的新表面项进行讨论。第三部分是MgB₂临界温度的理论计算。第四部分是文章的结论。

2. 两带Ginzburg-Landau理论和超导体–绝缘体界面的边界条件

根据Zhitomirsky和Dao的理论工作 [21]，我们可以写出两带超导体的哈密顿量为：

$H={\displaystyle {\sum }_{i\sigma }{c}_{i\sigma }^{†}\left(r\right)\hat{h}\left(r\right)c_{i\sigma }\left(r\right)}-{\displaystyle {\sum }_{i{i}^{\prime }}{g}_{i{i}^{\prime }}{c}_{i\uparrow }^{†}\left(r\right)c_{i\downarrow }^{†}\left(r\right)c_{i^{\prime }\downarrow }\left(r\right)c_{i^{\prime }\uparrow }\left(r\right)}$ (1)

方程中 $i,i^{\prime }=1,2$ 指能隙序数， $\sigma =\uparrow ,\downarrow$ 是自旋态。 $\hat{h}\left(r\right)$ 表示正常金属单粒子哈密顿量， $g_{i{i}^{\prime }}$ 为散射振幅且 $g_{12}=g_{21}$。

为简化哈密顿量的形式，引入能隙函数

${\Delta }_i\left(r\right)=-{\displaystyle {\sum }_{i^{\prime }}{g}_{i{i}^{\prime }}\langle c_{i^{\prime }\downarrow }\left(r\right)c_{i^{\prime }\uparrow }\left(r\right)\rangle }$ (2)

在平均场近似下，我们由此可以得到有效哈密顿量的形式为

$H_{eff}={\displaystyle {\sum }_{i\sigma }{c}_{i\sigma }^{†}\left(r\right)\hat{h}\left(r\right)c_{i\sigma }\left(r\right)}+{\displaystyle {\sum }_i\left[{\Delta }_i\left(r\right)c_{i\uparrow }^{†}\left(r\right)c_{i\downarrow }^{†}\left(r\right)+\text{H}\text{.c}\text{.}\right]}$ (3)

通过引入描述准粒子激发的产生湮灭算符 $b^{†},b$，我们可以利用Bogoliubov变换对有效的哈密顿量形式进行对角化

$c_{i\uparrow }\left(r\right)={\displaystyle {\sum }_k\left[u_{ik}\left(r\right)b_{ik\uparrow }-v_{ik}^{*}\left(r\right)b_{ik\downarrow }^{†}\right]}$ (4)

$c_{i\downarrow }\left(r\right)={\displaystyle {\sum }_k\left[u_{ik}\left(r\right)b_{ik\downarrow }+v_{ik}^{*}\left(r\right)b_{ik\uparrow }^{†}\right]}$ (5)

方程中 $k$ 代表波矢量。对角化后的哈密顿量可写为

$H_{eff}=E_g+{\displaystyle {\sum }_{ik\sigma }{E}_{ik}{b}_{ik\sigma }^{†}{b}_{ik\sigma }}$ (6)

$E_g$ 是基态的能量， $E_{ik}$ 是激发态的能量。

利用对易关系 $\left[c_{\sigma }\left(r\right),H_{eff}\right]$ 和公式(4)~(6)，我们可以推导出两带Bogoliubov-de Gennes方程

$\left(\begin{array}{cc}\hat{h}\left(r\right)& {\Delta }_i\left(r\right)\\ {\Delta }_i^{*}\left(r\right)& -{\hat{h}}^{*}\left(r\right)\end{array}\right)\left(\begin{array}{c}{u}_{ik}\left(r\right)\\ v_{ik}\left(r\right)\end{array}\right)=E_{ik}\left(\begin{array}{c}{u}_{ik}\left(r\right)\\ v_{ik}\left(r\right)\end{array}\right)$ (7)

由公式(6)，我们同时也能得到 $\langle b_{ik\uparrow }^{†}{b}_{ik\uparrow }\rangle =f\left(E_{ik}\right)$，其中 $f\left(E_{ik}\right)={\left[1+\exp \left(E_{ik}/k_{B}T\right)\right]}^{-1}$。结合公式(2)，可得能隙函数的自洽方程为

${\Delta }_i\left(r\right)={\displaystyle {\sum }_{i^{\prime }k}{g}_{i{i}^{\prime }}{v}_{i^{\prime }k}^{*}\left(r\right)u_{i^{\prime }k}\left(r\right)\times \left[1-2f\left(E_{i^{\prime }k}\right)\right]}$ (8)

与单带类似 [22]，当 ${\Delta }_i$ 为小量时，我们可以通过公式(7)和(8)得到能隙自洽方程的线性表达式

${\Delta }_i\left(r\right)={\displaystyle {\sum }_i{\displaystyle \int K_{i{i}^{\prime }}\left(r,r^{\prime }\right){\Delta }_{i^{\prime }}\left(r^{\prime }\right)\text{d}{r}^{\prime }}}$ (9)

其中

$K_{i{i}^{\prime }}\left(r,r^{\prime }\right)=g_{i{i}^{\prime }}{k}_{B}T{\displaystyle {\sum }_{k{k}^{\prime }\omega }\frac{{\Phi }_{i^{\prime }k}^{*}\left(r^{\prime }\right){\Phi }_{i^{\prime }{k}^{\prime }}^{*}\left(r^{\prime }\right){\Phi }_{i^{\prime }k}\left(r\right){\Phi }_{i^{\prime }{k}^{\prime }}\left(r\right)}{\left({\epsilon }_{i^{\prime }k}-i\hslash \omega \right)\left({\epsilon }_{i^{\prime }{k}^{\prime }}+i\hslash \omega \right)}}$ (10)

${\Phi }_{i^{\prime }k}\left(r\right)$ 是正常金属电子本征函数， $\hat{h}{\Phi }_{i^{\prime }k}={\epsilon }_{i^{\prime }k}{\Phi }_{i^{\prime }k}$。 $\omega =\left(2\upsilon +1\right)\pi k_{B}T/\hslash$ 代表频率， $\upsilon$ 是整数。

我们结合块体材料 $K_{i{i}^{\prime }}\left(r,r^{\prime }\right)$ 的精确表达形式并在能隙方程中加入非线性项，由公式(9)可以得到两带GL方程 [18]

${\alpha }_1\left(T\right){\Delta }_1+{\beta }_1{\left|{\Delta }_1\right|}^2{\Delta }_1-{\gamma }_1{\nabla }^2{\Delta }_1-R_{12}{\Delta }_2=0$ (11)

${\alpha }_2\left(T\right){\Delta }_2+{\beta }_2{\left|{\Delta }_2\right|}^2{\Delta }_2-{\gamma }_2{\nabla }^2{\Delta }_2-R_{12}{\Delta }_1=0$ (12)

方程中参数取值如下

${\alpha }_{1,2}=N_{1,2}\left[\frac{{\lambda }_{22,11}}{\lambda }-\frac{1}{{\lambda }_{\max }}-\ln \left(\frac{T_{c0}}{T}\right)\right]$， ${\beta }_i=\frac{7\zeta \left(3\right)N_i}{16{\pi }^2{\left(k_B{T}_{c0}\right)}^2}$ (13)

${\gamma }_i=\frac{7\zeta \left(3\right){\hslash }^2{N}_i{v}_{Fi}^2}{16{\pi }^2{\left(k_B{T}_{c0}\right)}^2}$， $R_{12}=\frac{N_1{\lambda }_{12}}{\lambda }=\frac{N_2{\lambda }_{21}}{\lambda }$ (14)

其中 ${\lambda }_{i{i}^{\prime }}=g_{i{i}^{\prime }}{N}_{i^{\prime }}$， $N_{i^{\prime }}$ 为两个能带各自费米能级态密度。 $\lambda ={\lambda }_{11}{\lambda }_{22}-{\lambda }_{12}{\lambda }_{21}$ 以及 ${\lambda }_{\max }=\frac{1}{2}\left[\left({\lambda }_{11}+{\lambda }_{22}\right)+\sqrt{{\left({\lambda }_{11}-{\lambda }_{22}\right)}^2+4{\lambda }_{12}{\lambda }_{21}}\right]$ 为 $\lambda$ 矩阵的最大特征值。 $T_{c0}$ 是块体材料的临界温度， $v_{Fi}$ 表示每个能带的平均费米速度。

在各向同性的均匀体系中，我们可以略去两带GL方程的梯度项。当 $T=T_{c0}$ 时，由(11)、(12)式得能隙自洽方程

$\left(\begin{array}{cc}{\lambda }_{11}& {\lambda }_{12}\\ {\lambda }_{21}& {\lambda }_{22}\end{array}\right)\left(\begin{array}{c}{\Delta }_1\\ {\Delta }_2\end{array}\right)={\lambda }_{\max }\left(\begin{array}{c}{\Delta }_1\\ {\Delta }_2\end{array}\right)$ (15)

与此同时，根据公式(11)和(12)我们可以写出弱耦合的两带GL自由能表达式

$F_V=F_1+F_2+F_{12}$ (16)

其中

$F_i={\displaystyle {\int }_V\left[{\alpha }_i\left(T\right){\left|{\Delta }_i\right|}^2+\left({\beta }_i/2\right){\left|{\Delta }_i\right|}^4+{\gamma }_i{\left|\nabla {\Delta }_i\right|}^2\right]\text{d}r}$ (17)

以及

$F_{12}={\displaystyle {\int }_V\left(-R_{12}{\Delta }_1^{*}{\Delta }_2+\text{c}\text{.c}\text{.}\right)\text{d}r}$ (18)

考虑超导体–绝缘体(或真空)间存在界面，我们添加新表面项

$F_S={\displaystyle {\int }_S\left({\displaystyle {\sum }_{i{i}^{\prime }}{\gamma }_i{D}_{i{i}^{\prime }}{\Delta }_i^{*}{\Delta }_{i^{\prime }}}\right)\text{d}s}$ (19)

该积分在超导体表面S进行，其中 $D_{i{i}^{\prime }}$ 为常数。

于是总的自由能可表示为

$F=F_V+F_S$ (20)

对 ${\Delta }_i^{*}$ 做变分，得

$\begin{array}{c}\delta F={\displaystyle {\int }_V\left({\alpha }_1{\Delta }_1+{\beta }_1{\left|{\Delta }_1\right|}^2{\Delta }_1-{\gamma }_1{\nabla }^2{\Delta }_1-R_{12}{\Delta }_2\right)\delta {\Delta }_1^{*}\text{d}r}\\ +{\displaystyle {\int }_V\left({\alpha }_2{\Delta }_2+{\beta }_2{\left|{\Delta }_2\right|}^2{\Delta }_2-{\gamma }_2{\nabla }^2{\Delta }_2-R_{12}{\Delta }_1\right)\delta {\Delta }_2^{*}\text{d}r}\\ +{\displaystyle {\int }_S{\gamma }_1\left(\nabla {\Delta }_1\cdot \hat{s}+D_{11}{\Delta }_1+D_{12}{\Delta }_2\right)\delta {\Delta }_1^{*}\text{d}s}\\ +{\displaystyle {\int }_S{\gamma }_2\left(\nabla {\Delta }_2\cdot \hat{s}+D_{21}{\Delta }_1+D_{22}{\Delta }_2\right)\delta {\Delta }_2^{*}\text{d}s}\end{array}$ (21)

当满足条件 $\delta F/\delta {\Delta }_i^{\text{*}}=0$ 时，我们不仅可以得到两带GL方程(11)和(12)，而且可得边界条件

${\nabla {\Delta }_i\cdot \hat{s}|}_S=-{\displaystyle {\sum }_{i^{\prime }}{D}_{i{i}^{\prime }}{\Delta }_{i^{\prime }}}$ (22)

对于通常应用于超导体–绝缘体(或真空)界面的Neumann边界条件，我们有 $D_{i{i}^{\prime }}=0$，但在此条件下超导临界温度与薄膜厚度的依赖关系无法得到合理解释。我们基于de Gennes [22] 关于单带超导体边界条件的理论计算，推广得到两带超导体边界条件中 $D_{i{i}^{\prime }}$ 的表达式为

$D_{i{i}^{\prime }}=\frac{a{N}_{i^{\prime }}}{{\gamma }_i{\lambda }_{\max }}$ (23)

我们认为由此表面项引起的边界效应将为解释 $T_c$ 被抑制的行为提供一种可能。

3. MgB₂超导薄膜临界温度的计算

在本节中，我们将根据边界效应解释临界温度 $T_c$ 与MgB₂薄膜厚度之间的依赖关系。设薄膜厚度为d，其边界位置由 $z=-d/2$ 和 $z=d/2$ 决定。

根据公式(11)和(12)，两带GL方程可写为

$\left(\begin{array}{cc}{\hat{H}}_{11}& {\hat{H}}_{12}\\ {\hat{H}}_{21}& {\hat{H}}_{22}\end{array}\right)\left(\begin{array}{c}{\Delta }_1\left(r\right)\\ {\Delta }_2\left(r\right)\end{array}\right)=0$ (24)

其中

${\hat{H}}_{i{i}^{\prime }}=-{\gamma }_i{\nabla }^2+{\alpha }_i\left(T\right)$ (25)

和

${\hat{H}}_{12}={\hat{H}}_{21}=-R_{12}$ (26)

由于临界温度附近序参量较小，我们忽略了包含 $\beta$ 的高阶项贡献。

与单带超导体类似 [23]，我们设超导薄膜能隙方程的解为

$\left(\begin{array}{c}{\Delta }_1\left(z\right)\\ {\Delta }_2\left(z\right)\end{array}\right)=\left(\begin{array}{c}{\mu }_1\\ {\mu }_2\end{array}\right)\cos \left(kz\right)$ (27)

为常数。由边界条件

${\frac{\text{d}{\Delta }_i}{\text{d}z}|}_{z=\pm d/2}=\mp {\displaystyle {\sum }_{i^{\prime }}{D}_{i{i}^{\prime }}{\Delta }_{i^{\prime }}}$ (28)

可得

$\left(\begin{array}{cc}k\tan \left(\frac{kd}{2}\right)-D_{11}& -D_{12}\\ -D_{21}& k\tan \left(\frac{kd}{2}\right)-D_{22}\end{array}\right)\left(\begin{array}{c}{\mu }_1\\ {\mu }_2\end{array}\right)=0$ (29)

由

$\left|\begin{array}{cc}k\tan \left(\frac{kd}{2}\right)-D_{11}& -D_{12}\\ -D_{21}& k\tan \left(\frac{kd}{2}\right)-D_{22}\end{array}\right|=0$ (30)

得波矢k满足方程

$\left[k\tan \left(\frac{kd}{2}\right)-D_{11}\right]\left[k\tan \left(\frac{kd}{2}\right)-D_{22}\right]=D_{12}{D}_{21}$ (31)

同时，从公式(24)和(27)出发，两带GL方程可改写为

$\left(\begin{array}{cc}{H}_{11}& H_{12}\\ H_{21}& H_{22}\end{array}\right)\left(\begin{array}{c}{\mu }_1\\ {\mu }_2\end{array}\right)=0$ (32)

其中

$H_{ii}={\gamma }_i{k}^2+{\alpha }_i\left(T\right)$ (33)

和

$H_{12}=H_{21}=-R_{12}$ (34)

在 $T=T_c$ 时，决定两带超导薄膜临界温度的条件为

$H_{11}{H}_{22}-H_{12}{H}_{21}=0$ (35)

对于典型的两带超导体MgB₂，已知 $T_{c0}\approx 40\text{K}$ [1]，平均晶格常数 $a\approx 0.4\text{nm}$ [24]。 $\sigma$ 带费米面的态密度是 $N_1=0.16{\text{eV}}^{-1}$， $\pi$ 带费米面的态密度是 $N_2=0.25{\text{eV}}^{-1}$ [2]。由比热数据的理论拟合得 ${\lambda }_{11}=g_{11}{N}_1=0.4$，比值 $g_{11}:g_{22}:g_{12}=1:0.3:0.2$ [21]，据此我们可以计算出 $\lambda =0.066$， ${\lambda }_{\max }=0.44$， $R_{12}=0.30{\text{eV}}^{-1}$。结合费米平均速度 $v_{F1}=3.7\times {10}^{14}\text{nm}\cdot {\text{s}}^{-1}$， $v_{F2}=4.5\times {10}^{14}\text{nm}\cdot {\text{s}}^{-1}$ [2]，公式(14)给出 ${\gamma }_1=39{\text{nm}}^2\cdot {\text{eV}}^{-1}$， ${\gamma }_2=87{\text{nm}}^2\cdot {\text{eV}}^{-1}$。将所得参数带入公式(23)，最终得特征尺度 $D_{11}={\left(268\text{nm}\right)}^{-1}$， $D_{12}={\left(171\text{nm}\right)}^{-1}$， $D_{21}={\left(598\text{nm}\right)}^{-1}$， $D_{22}={\left(383\text{nm}\right)}^{-1}$。基于公式(31)和(35)，我们给出了超导临界温度随薄膜厚度的变化趋势，如图1所示。从图中可以看到临界温度随薄膜厚度的减小而逐渐降低，当膜厚度降低到10 nm时，临界温度 $T_c\approx 33\text{K}$。值得注意的是，该理论结果与不同衬底生长的MgB₂实验数据均相符。

Figure 1.The critical temperature as a function of MgB₂ film thickness

图1. 超导临界温度随MgB₂薄膜厚度的变化曲线

4. 结论

在本文中，我们对两带超导体–绝缘体(或真空)之间的边界效应展开细致研究，发现对于典型的两带超导体二硼化镁Neumann边界条件无法解释薄膜厚度与临界温度的依赖关系。于是我们在两带GL理论中引入了新表面项，并给出表面项参数的微观表达式，该理论能够很好地解释当薄膜厚度减少时临界温度被显著抑制的实验现象。研究结果表明，在MgB₂等多带超导薄膜的研究中由新表面项引起的边界效应可能产生不可忽略的物理效应。