様々な結晶の陽電子寿命の第一原理計算#

PHASE logo

物質に打ち込まれた陽電子は、原子核の正電荷に反発されて格子間の空隙に溜まり、そこで電子と対消滅します。消滅までの時間（陽電子寿命）は陽電子が感じる電子密度で決まるため、原子空孔のような開いた欠陥があると寿命が延びます。これを利用した陽電子消滅法（PAT）は、非破壊で原子スケールの空孔型欠陥を検出できる数少ない手法ですが、測定された寿命から欠陥の種類とサイズを特定するには、「欠陥のない完全結晶なら何 ps になるはずか」を信頼できる精度で計算できることが前提になります。本解析では、第一原理計算ソフトウェアAdvance/PHASEに実装された二成分密度汎関数法を用い、金属 9 種・バンドギャップ系 7 種の計 16 種類の完全結晶についてバルク陽電子寿命を計算しました。電子には一般化勾配近似（GGA-PBE）、電子–陽電子相関には局所密度近似（Boroński–Nieminen）を用い、バンドギャップのある系には高周波誘電率 $\epsilon_\infty$ による補正を適用した結果、平均絶対誤差 4.5 %、RMSD 9.8 psに収まりました。あわせて、PHASEが出力する陽電子密度・電子密度・電子–陽電子対密度の cube ファイルを解析し、寿命の長短が結晶構造のどの特徴で決まっているのかを可視化しています。

Keywords: 第一原理計算 (DFT), 二成分密度汎関数法, 陽電子消滅法 (PAT), 陽電子寿命, 空孔型欠陥, Boroński–Nieminen, 増強因子, 高周波誘電率補正, cubeファイル解析

計算手法#

二成分密度汎関数法による陽電子寿命#

結晶中の陽電子は、通常は 1 個だけが存在する希薄な極限にあります。この場合、陽電子は電子系に影響を与えず、電子密度 $n_e(\mathbf{r})$ が作るポテンシャル中の 1 粒子問題として解くことができます（バルクモード、PHASEでは positron = BULK）。陽電子波動関数 $\psi_+$ が満たすのは次の式です。

$\left[-\tfrac{1}{2}\nabla^2 - V_\text{Coul}(\mathbf{r}) + V_\text{corr}[n_e(\mathbf{r})]\right]\psi_+(\mathbf{r}) = \varepsilon_+ \psi_+(\mathbf{r})$

ここで $V_\text{Coul}$ は電子とイオンが作る静電ポテンシャルを電子について書いたもので、陽電子は正電荷なので符号が反転して式中では $-V_\text{Coul}$ となります（イオン核から反発され、電子雲に引かれる）。 $V_\text{corr}$ は電子–陽電子相関ポテンシャルで、本計算では局所密度近似（Boroński–Nieminen [1]）を用います。この定式化の詳細は文献 [2] のレビューに、二成分DFTの自己無撞着な取り扱いは文献 [3] にまとめられています。得られた陽電子密度 $n_+ = |\psi_+|^2$ から、消滅率 $\lambda$ と寿命 $\tau$ が求まります。

$\lambda = \dfrac{1}{\tau} = \pi r_0^2 c \displaystyle\int n_+(\mathbf{r})\, n_e(\mathbf{r})\, \gamma[n_e(\mathbf{r})]\, d\mathbf{r}$

$r_0$ は古典電子半径、 $c$ は光速です。積分の中身が本解析の主役で、陽電子密度と電子密度が重なっている場所でだけ消滅が起こることを表しています。 $\gamma$ は増強因子と呼ばれ、陽電子がクーロン引力で周囲の電子を引き寄せる効果（電子相関）を表します。PHASEに実装されている LDA 増強因子は、電子密度パラメータ $r_s = (3/4\pi n_e)^{1/3}$ を使って次のように書けます。

$\gamma(r_s) = 1 + 1.23\,r_s + 0.9889\,r_s^{3/2} - 1.482\,r_s^{2} + 0.3956\,r_s^{5/2} + \dfrac{r_s^{3}}{6}$

これは Boroński–Nieminen [1] が与えた補間形式と同じ関数形で、低密度極限（ $r_s \to \infty$ ）で最終項 $r_s^3/6$ が支配し、消滅率が正しくポジトロニウムの値に漸近するように作られています [3]。

バンドギャップ系に対する誘電率補正#

上式の最終項 $r_s^3/6$ は、陽電子のまわりに電子が完全に集まる（完全遮蔽）という前提に対応し、低密度極限でポジトロニウム形成の描像を再現するための項です。ところがバンドギャップのある物質では、電子は束縛されていて完全には集まれません。Puskaら [4] はこの不完全遮蔽を高周波（光学）誘電率 $\epsilon_\infty$ で近似し、最終項を次のように置き換えることを提案しました。

$\dfrac{r_s^{3}}{6} \;\longrightarrow\; \dfrac{r_s^{3}}{6}\left(1 - \dfrac{1}{\epsilon_\infty}\right)$

PHASEでは positron_convergence{} ブロックに epsilon_ele を書くだけでこの補正が有効になります（既定はOFF＝金属扱い、 $\epsilon_\infty = \infty$ ）。本解析では 7 種のギャップ系すべてにこの補正を適用し、比較のため補正なしの計算も行いました。

計算対象と計算条件#

対象は文献 [5] と同じ 16 系です。格子定数はすべて実験値（室温）を用い、構造最適化は行っていません。擬ポテンシャルは、陽電子が滞在する格子間領域の電子密度を正しく再現することが重要です。消滅率への寄与が大きい半芯電子は、可能な限り価電子として扱いました（3d〜5d 遷移金属の d 電子、Ga 3d、In 4d、Cd 4d、Te 4d、さらに Nb の 4s4p と K の 3p）。一方 Ge と As については 3d をコアに残す擬ポテンシャルしか利用できず、コア電子との消滅の扱いが他系とは異なります。

表1. 計算条件

項目	設定
計算モード	`positron = BULK`（零陽電子密度極限）
交換相関汎関数（電子）	GGA-PBE
電子–陽電子相関	LDA（Boroński–Nieminen）
ギャップ系の補正	`epsilon_ele` = $\epsilon_\infty$（表2参照）
擬ポテンシャル	ノルム保存／ウルトラソフト
カットオフ（波動関数／電荷密度）	36–50 Ry ／ 144–405 Ry
k点サンプリング	金属 12×12×12（既約 231 点：fcc、140 点：bcc、Fe・Ni はスピン 2 成分で 2 倍）、バンドギャップ系 4×4×4–6×6×6（既約 19–60 点）、いずれも四面体法
スピン分極	Fe・Ni のみ考慮（強磁性）
陽電子固有値の収束判定	1×10⁻⁷ Ry

Fe と Ni は強磁性体であり、電子–陽電子相関汎関数は全電子密度 $n_\uparrow + n_\downarrow$ の汎関数として評価する必要があります。本計算ではスピン分極を明示的に考慮しています。

計算結果と考察#

陽電子寿命：実験値との比較#

16 系すべての計算結果を実験値に対してプロットしたものが図1です。文献 [1] と同じ形式（金属＝四角、バンドギャップ系＝丸、誤差 5 % と 10 % の補助線）で示しています。

図1. 16 種類の結晶のバルク陽電子寿命：Advance/PHASEによる計算値と実験値の比較。14 系が誤差 10 % 以内（うち 11 系は 5 % 以内）に収まり、最も外れる Ni でも 11.0 %、平均絶対誤差は 4.5 %、RMSD は 9.8 ps です。横棒は実験値の文献間ばらつき（表2の括弧内）を表します。左上は 90–135 ps の混み合った領域の拡大図です。

寿命は Ni の 97 ps から K の 376 ps まで約 3.9 倍の幅がありますが、計算値はその全域にわたって実験値に追随しています。実験値との差が最も大きいのは Ni（-11.0 %）、Cu（-10.8 %）、Fe（-9.1 %）といった 3d 遷移金属で、いずれも計算値が実験値より短くなる（＝消滅率を過大評価する）方向に系統的にずれています。これは LDA の増強因子が 3d 遷移金属で消滅率をやや過大評価するという、この手法に固有の既知の傾向で、同種のスキームによる文献の計算値にも同じ方向・同じ程度のずれが見られます。ただし Cu の実験値は文献間で 110.7–122 ps とばらついており、この程度の差は実験値の不確かさと同程度である点にも注意が必要です。

表2. 全16系の計算結果と実験値

系	構造	$\epsilon_\infty$	本計算 (ps)	実験値 (ps)
Al	fcc	—（金属）	167.39	165.0 (160.7–166.0)
Cu	fcc	—（金属）	107.02	120.0 (110.7–122.0)
Ag	fcc	—（金属）	123.94	130.0 (120.0–133.0)
Au	fcc	—（金属）	110.78	116.0 (113.0–123.0)
Ni	fcc（強磁性）	—（金属）	97.00	109.0 (105.0–110.0)
Fe	bcc（強磁性）	—（金属）	100.88	111.0 (106.0–114.0)
Nb	bcc	—（金属）	123.95	120.0 (119.0–125.0)
Na	bcc	—（金属）	321.28	338.0
K	bcc	—（金属）	375.78	397.0
Ge	ダイヤモンド構造	16.00	222.72	228.0 (220.5–230.0)
Si	ダイヤモンド構造	12.00	223.07	218.0 (216.0–222.0)
GaAs	閃亜鉛鉱型	10.89	226.85	231.6
InP	閃亜鉛鉱型	9.61	239.87	240.0
CdTe	閃亜鉛鉱型	7.10	279.85	285.0 (280.0–291.0)
C	ダイヤモンド構造	5.70	97.58	98.0 (97.5–115.0)
HfO₂	単斜晶 P2₁/c	4.00	155.52	170.0

格子定数は実験値（室温）。陽電子寿命の実験値は [6,7] から。括弧内は文献間のばらつきで、系によっては 10 % 以上あるため、計算値との差は幅をもって見る必要があります。

なお同種の相関スキームを用いた他コードの計算値 [6,7] とも、15 系の平均絶対差 3.1 ps で一致しています。

バンドギャップ系における誘電率補正の効果#

バンドギャップ系 7 種について、 $\epsilon_\infty$ 補正の有無を比較したものが図2(a)と表3です。

誘電率補正の効果
図2. (a) $\epsilon_\infty$ 補正による実験値からのずれの変化。補正なしでは 7 系すべてが実験値を 1.4–15.1 % 過小評価しますが、補正を入れると HfO₂ を除く 6 系が ±2.5 % 以内に収まります。(b) 補正による寿命の増分は、cube ファイルから導かれる解析式で 0.3 ps 以内に再現できます。

表3. $\epsilon_\infty$ 補正の有無による陽電子寿命の比較

系	$\epsilon_\infty$	補正なし (ps)	補正あり (ps)	実験値 (ps)
Ge	16.00	216.48 -5.1 %	222.72 -2.3 %	228.0
Si	12.00	215.05 -1.4 %	223.07 +2.3 %	218.0
GaAs	10.89	217.74 -6.0 %	226.85 -2.1 %	231.6
InP	9.61	228.46 -4.8 %	239.87 -0.1 %	240.0
CdTe	7.10	259.44 -9.0 %	279.85 -1.8 %	285.0
C	5.70	94.36 -3.7 %	97.58 -0.4 %	98.0
HfO₂	4.00	144.31 -15.1 %	155.52 -8.5 %	170.0

補正の効き幅は $\epsilon_\infty$ が小さいほど、また寿命が長いほど大きく、CdTe で +20.4 ps に達します。この振る舞いには簡単な説明がつきます。補正で削られるのは増強因子の $r_s^3/6$ 項だけですが、消滅率の被積分関数の中ではこの項は

$n_e(\mathbf{r})\cdot\dfrac{r_s(\mathbf{r})^3}{6} = \dfrac{1}{6}\cdot\dfrac{3}{4\pi} = \dfrac{1}{8\pi}$

と、 $r_s$ の定義（ $n_e r_s^3 = 3/4\pi$ ）から位置にも物質にもよらない定数になります。したがって陽電子密度の規格化 $\int n_+\,d\mathbf{r} = 1$ と合わせると、補正による消滅率の減少分は

$\Delta\lambda = -\,\pi r_0^2 c \cdot \dfrac{1}{8\pi\,\epsilon_\infty}$

と物質によらず $1/\epsilon_\infty$ だけで決まり、寿命の増分は $\Delta\tau \simeq \tau^2\,|\Delta\lambda|$ すなわち寿命の 2 乗に比例します。CdTe の補正量が大きいのは、 $\epsilon_\infty$ が小さいことに加えて寿命そのものが長いためです。

この見積りは、実際の自己無撞着計算の結果を最大 0.27 ps の差で再現しています（図2(b)）。残差の出どころもはっきりしています。 $\epsilon_\infty$ 補正は陽電子ポテンシャルには入らないため、陽電子密度・電子密度は補正の有無で完全に一致します（出力される cube ファイルがバイト単位で同一であることを確認しました）。上式は価電子消滅だけを扱っていますが、PHASEはコア電子との消滅にも同じ補正を適用するため、その分（寿命に対して 0.02〜0.12 %）が差として残ります。

cubeファイルからの解析：寿命の長短はどこで決まるか#

PHASEは陽電子計算の実行時に 3 つの cube ファイルを自動出力します。positron.cube（陽電子密度 $n_+$ ）、electron.cube（価電子密度 $n_e$ ）、そして ep_pair.cube（電子–陽電子対密度 $n_+ n_e \gamma$ 、すなわち消滅率の被積分関数そのもの）です。これらを使うと、寿命という 1 つの数値の背後にある空間分布を直接見ることができます。

図3は、寿命の短い順に並べた 6 種類の結晶について、陽電子密度の断面を示したものです。6 面すべて同じ 7.2 Å × 5.2 Å の窓・同じカラースケールで描いており、直接比較できます。

陽電子密度の断面図
図3. 結晶種による陽電子密度分布の違い（寿命の短い順）。青が濃いほど陽電子密度が高い領域です。細い等高線は価電子密度で、陽電子が電子雲の「谷」に入り込んでいることが分かります。

どの結晶でも、陽電子はイオン核の位置から完全に排除され（白く抜けている領域）、格子間の空隙に分布しています。違いはその空隙の広さです。ダイヤモンドは原子間距離が 1.54 Å と短く価電子密度が非常に高いため、陽電子が入り込める低密度領域がほとんどなく、 $\langle n_e \rangle_+ = 0.07194$ bohr⁻³ と 16 系で最も高い密度を感じます。その結果、寿命は 97.6 ps とギャップ系で最短（16 系全体でも Ni の 97.0 ps に次ぐ短さ）になります。同じ共有結合・同じダイヤモンド構造でも、格子定数が 1.5 倍大きい Si では格子間チャネルがはるかに広く、陽電子は電子密度の低い場所に留まれるため（ $0.01961$ bohr⁻³）、寿命は 223.1 ps と 2.3 倍に伸びます。Na は価電子が 1 原子あたり 1 個しかない低密度金属で、セル全体がほぼ空隙として働き、321.3 ps という長い寿命になります。HfO₂ はイオン結晶で、Hf⁴⁺ と O²⁻ の非対称な配置が作る空隙に陽電子が集まります。陽電子が実質的に占める体積はセル体積の 60 % 程度で、図3の 6 系のなかでは最も局在的です。

図4(a)は、Si の共有結合を [111] 方向に横切ったときの 3 つの分布です。

cubeファイルからの定量解析
図4. (a) Si の [111] 方向に沿った価電子密度・陽電子密度・消滅密度。3 つの量はそれぞれ別の場所でピークを持ちます。(b) 全16系について、価電子消滅だけの寿命 $\tau_\text{val}$ を、陽電子が実際に感じる価電子密度 $\langle n_e \rangle_+ = \int n_+ n_e\, d\mathbf{r}$ に対してプロットしたもの。

価電子密度は結合中心で最大になり、陽電子密度はその正反対に格子間サイトで最大になります。そして両者に増強因子を掛けた消滅密度は、そのどちらとも違う分布を示します。Si の場合、結合中心の価電子密度は格子間サイトの 25.7 倍ある一方、陽電子密度は逆に格子間サイトが結合中心の 5.5 倍です。両者に増強因子 $\gamma$ を掛け合わせた消滅密度は、この 2 か所でほぼ同じ大きさ（比 0.95）になります。一方、イオン核の位置では陽電子が完全に排除されるため消滅密度は 3 桁以上小さくなります。つまり「陽電子がいる場所」でも「電子がいる場所」でもなく、その重なりが寿命を決めています。

この描像を 16 系すべてで定量化したのが図4(b)です。陽電子密度で重み付けした平均電子密度 $\langle n_e \rangle_+$ を横軸にとると、価電子消滅だけの寿命はきれいな冪則に乗ります。fcc・bcc・強磁性体を含む金属 9 系が $\tau_\text{val} \propto \langle n_e \rangle_+^{-0.46}$ （ $R^2 = 0.993$ ）、共有結合・イオン結合を含むギャップ系 7 系が $\langle n_e \rangle_+^{-0.62}$ （ $R^2 = 0.989$ ）です。フィットからの残差は金属で最大 9.8 %（Al）、ギャップ系で最大 8.0 %（HfO₂）です。結晶構造や化学結合の違いは、その大半が「陽電子がどこに落ち着き、そこの電子密度がいくらか」という 1 つの量に集約されて寿命に効いていることが分かります。

ただし 2 つの系列は重なりません。金属のフィットを同じ密度まで外挿して比べると、ダイヤモンドを除くギャップ系 6 系は寿命が 16–28 % 長く、ダイヤモンドだけは逆に 6 % 短くなっています。この隔たりのうち $\epsilon_\infty$ 補正で説明できるのは一部で、補正を外して同じ比較をすると系列間のずれは 17 % から 11 % に縮まります。残りは、平均密度が同じでも電子密度の空間的な偏り方（したがって実効的な増強因子）が金属と共有結合系で異なることによります。「陽電子が感じる平均密度」は寿命の第一近似としては非常に強力ですが、それだけで全系を 1 本の線に載せられるわけではありません。

まとめ#

第一原理計算ソフトウェアAdvance/PHASEの二成分密度汎関数法を用いて、金属 9 種・バンドギャップ系 7 種の計 16 種類の完全結晶のバルク陽電子寿命を計算しました。11 系が実験値の ±5 % 以内、14 系が ±10 % 以内（平均絶対誤差 4.5 %、RMSD 9.8 ps）に収まり、金属から半導体・酸化物まで一貫した精度が得られること、ギャップ系では高周波誘電率による補正が 1.4–15.1 % 効くことを確認しました。cube ファイルの解析からは、陽電子が例外なくイオン核から排除されて格子間の空隙に分布し、消滅がその陽電子密度と電子密度の重なりで決まること、16 系の寿命が「陽電子が実際に感じる価電子密度」に対して系列ごとに冪則に乗ること（金属 $R^2 = 0.993$ 、バンドギャップ系 $R^2 = 0.989$ ）が分かりました。完全結晶の寿命をこの精度で予測できることは、陽電子消滅測定の解釈に直接役立ちます。測定で得られた寿命が本計算のバルク値より長ければ空孔型欠陥の存在が示され、その差の大きさから欠陥のサイズを推定できるためです。

本解析の詳細や、研究への適用可能性に関するご相談はこちら

お問い合わせ

参考文献#

E. Boroński and R. M. Nieminen, "Electron-positron density-functional theory", Phys. Rev. B 34, 3820 (1986).
M. J. Puska and R. M. Nieminen, "Theory of positrons in solids and on solid surfaces", Rev. Mod. Phys. 66, 841 (1994).
M. J. Puska, A. P. Seitsonen, and R. M. Nieminen, "Electron-positron Car-Parrinello methods: Self-consistent treatment of charge densities and ionic relaxations", Phys. Rev. B 52, 10947 (1995).
M. J. Puska, S. Mäkinen, M. Manninen, and R. M. Nieminen, "Screening of positrons in semiconductors and insulators", Phys. Rev. B 39, 7666 (1989).
M. Saito, A. Nakamoto, T. Yamasaki, and M. Okamoto, "Two-component density functional calculations on lifetimes of positrons in a variety of crystals", The ISSP International Workshop and Symposium on "Foundations and Applications of the Density Functional Theory", The Institute for Solid State Physics, the University of Tokyo, Aug. 1 (2007).
J. M. Campillo Robles, E. Ogando, and F. Plazaola, "Positron lifetime calculation for the elements of the periodic table", J. Phys.: Condens. Matter 19, 176222 (2007).
W. Zhang, B. Gu, J. Liu, and B. Ye, "Accurate theoretical prediction on positron lifetime of bulk materials", Comput. Mater. Sci. 105, 32 (2015).