ヘマタイト (α-Fe₂O₃) 中の電子・正孔ポーラロンの局在状態計算：自己束縛と2サイト構造

ヘマタイト (α-Fe₂O₃) は可視光を吸収し、安価で化学的にも安定な酸化物半導体で、光電気化学的な水分解の光アノード材料として盛んに研究されています。一方で、キャリアの移動度が非常に小さく、光で生成した電子や正孔の多くが表面に届く前に再結合してしまうことが、性能を制約する要因となっています。移動度が小さいのは、余剰の電子や正孔が周囲の格子を歪ませ、その歪みに自らを捕らえた「小ポーラロン (small polaron)」を Fe サイトに形成して、隣のサイトへの熱活性化ホッピングによって移動するためと考えられています。第一原理計算でポーラロンを扱うには、まず Fe³⁺ (3d⁵、高スピン) が反強磁性的に配列した母体の絶縁体状態を安定に求め、そのうえで余剰キャリアが局在した状態を探す必要があります。本解析では、第一原理計算ソフトウェア Advance/PHASE (ver.5.0) の DFT+U 機能を用いて、α-Fe₂O₃ 中の余剰電子と正孔の局在状態を計算しました。d 軌道の占有行列の初期値を工夫して反強磁性の母体を安定に求めたうえで、非局在の余剰電子が構造緩和の過程で自発的に自己束縛する様子を示し、得られる電子ポーラロンの形が Hubbard $U$ の値によって 1 サイト型と 2 サイト (Fe–Fe 対) 型の間で変わることを明らかにします。さらに、2 サイト・ポーラロンが移動する過程のエネルギー障壁を求め、ハイブリッド汎関数による文献値と比較します。正孔については、Fe 3d 軌道に $U$ (4.3 eV) を加える DFT+U では、対称性を破った初期状態から出発しても局在状態が得られないことを示し、手法の適用範囲を議論します。

Keywords: 第一原理計算, DFTシミュレーション, 小ポーラロン, DFT+U, ヘマタイトα-Fe₂O₃, 反強磁性, 自己束縛エネルギー, ホッピング伝導, 占有行列, 電子・正孔

計算方法：反強磁性母体中のポーラロン

反強磁性母体と DFT+U

α-Fe₂O₃ はコランダム型構造をとり、Fe³⁺ は高スピンの 3d⁵ 配置にあります。c 面内で稜を共有する Fe どうしはスピンの向きがそろった 2 重層 (bilayer) をつくり、この 2 重層が c 軸方向に交互に反平行に積み重なった反強磁性秩序が基底状態です [1] (Néel 温度は約 955 K [2])。GGA では Fe 3d 電子の局在が不十分で、バンドギャップが大きく過小評価されるため、Fe 3d 軌道に Hubbard $U$ 補正を加えた DFT+U を用いました [3] (以下、有効値 $U_{\text{eff}}$ を単に $U$ と書きます)。

反強磁性の Fe³⁺ 母体では、DFT+U の自己無撞着 (SCF) 計算の結果が初期状態に大きく左右されます。d 軌道の占有行列をゼロから始めると、最初の反復では Hubbard ポテンシャルがすべての d 軌道に一様な斥力として働くため、計算が低スピンの金属的な状態に落ち込むことがあります。そこで本解析では、各 Fe の d 軌道の占有行列を、サイト対称性に合わせた高スピン d⁵ の配置で初期化し、SCF 計算の初期段階ではこの占有行列を保ち、計算が落ち着いてから自由に変化させました。

余剰キャリアの導入と局在の判定

余剰電子 (正孔) は、系に電子を 1 個加える (取り除く) ことで導入し、系全体の電荷 −1 (+1) は一様な背景電荷で中和しました。以下、ポーラロンを局在させる Fe を Fe0 と呼びます (図1)。上向きスピンの Fe0 に入った余剰電子は、下向き (少数) スピンの d 軌道を占めて Fe²⁺ (3d⁶) をつくります。一方、正孔は上向き (多数) スピンの d 電子を 1 個取り去り、Fe⁴⁺ (3d⁴) をつくります。いずれの場合も全スピン $n_\uparrow - n_\downarrow$ は −1 となるため、この値に固定して計算しました。

母体の Fe はもともと大きな磁気モーメント (約 4 μ_B) を持つため、スピン密度の分布だけではポーラロンの位置を特定することはできません。そこで、DFT+U で計算される各 Fe の d 軌道の占有行列から、少数スピンの d 占有の増加量 $\delta n$ を求め、局在の指標としました。

$$\delta n_i = \mathrm{tr}\,\mathbf{n}_i^{\text{min}} - \underset{j}{\mathrm{median}}\; \mathrm{tr}\,\mathbf{n}_j^{\text{min}}$$

ここで $\mathbf{n}_i^{\text{min}}$ はサイト $i$ の少数スピンの d 占有行列 (5×5) です。$\mathrm{tr}$ は行列のトレース（対角成分の和）であり、対象サイトの少数スピン d 軌道の合計電子数を表します。また、$\mathrm{median}$ は中央値であり、同じスピン副格子の Fe についてとります。正孔では、多数スピンの d 占有の減少量を同様に評価します。Fe 3d と O 2p の混成や、占有行列を定義する射影子の取り方の影響により、1 つの Fe に局在した電子でも $\delta n$ は 1 にはなりませんが、余剰電子がどの Fe にどれだけ分布しているかは明確に読み取れます。この指標は SCF 計算の過程で直接得られるため、構造緩和の各ステップで局在の様子を追うこともできます。

自己束縛エネルギー

ポーラロンの束縛の強さは、自己束縛エネルギー $E_{\text{ST}}$ で評価します。

$$E_{\text{ST}} = E_{\text{deloc}} - E_{\text{loc}}$$

ここで $E_{\text{deloc}}$ は、緩和した中性母体の構造のまま余剰電子が全 Fe に広がった非局在状態の全エネルギー、$E_{\text{loc}}$ は構造緩和した局在状態の全エネルギーです。エネルギーの比較には、すべて同じ k 点の組で計算した値を用います。

計算モデルと計算条件

計算モデルには、コランダム型 α-Fe₂O₃ の六方晶の慣用単位胞 (空間群 $R\bar{3}c$、30 原子) を 2×2×1 に拡張した 120 原子 (Fe₄₈O₇₂) のスーパーセルを用いました (図1a)。格子定数は実験値 ($a$ = 5.035 Å、$c$ = 13.747 Å) に固定し、内部座標は $U$ = 4.3 eV で最適化しました。図1b に示すように、Fe0 は同じ 2 重層内の 3 つの Fe (原子 11、23、47) と 2.97 Å の距離で稜を共有しています。これらの Fe は Fe0 とスピンの向きが同じで、面内でのホッピングの相手になります。一方、c 軸方向に面を共有する Fe (2.90 Å) はスピンの向きが反対で、ホッピングの経路としては不利です。実験でも、面内の電気伝導度は c 軸方向より約 4 桁大きいことが知られています [4]。

図1. (a) 計算に用いた 2×2×1 スーパーセル (120 原子)。青と橙はスピンの向きが逆の Fe、緑はポーラロンを局在させる Fe0、灰色は O です。右側の矢印は、Fe 2 重層ごとのスピンの向きを示します。(b) Fe0 を含む Fe 2 重層を c 軸方向から見た図。塗りつぶしと白抜きは高さの異なる Fe を表します。Fe0 は 3 つの Fe (原子 11、23、47) と 2.97 Å で稜を共有しています。

主な計算条件を表1にまとめます。$U$ は主に 4.3 eV とし、ポーラロンの形の $U$ 依存性を調べるため、一部の計算は 5.0 eV でも行いました。

表1. 計算モデル・計算条件

原子構造	2×2×1スーパーセル / 120 原子 (Fe48O72)
格子定数	\(a\) = 5.035 Å、\(c\) = 13.747 Å (実験値に固定)
擬ポテンシャル	PAW
交換相関汎関数	GGA (PBE) ＋ Hubbard \(U\) (Fe 3d)、\(U_{\text{eff}}\) = 4.3 eV (一部 5.0 eV)
射影子 (DFT+U)	\(l\) = 2、半径 2.30 bohr
カットオフ (波動関数 / 電荷密度)	36 / 324 Rydberg
k点サンプリング	バルク (30 原子)：6×6×2、スーパーセル：2×2×1 (Monkhorst–Pack、c 軸方向に 1/2 ずらし。図4 の緩和の経過のみ面内方向にもずらした k 点)
磁性	反強磁性 (図1a)、全スピン \(n_\uparrow-n_\downarrow\) を固定 (中性 0、電子・正孔 −1)
系の電荷	−1 (余剰電子 1 個) / +1 (正孔 1 個)
構造最適化	内部座標の全原子緩和、最大力 \(1\times10^{-4}\) hartree/bohr

計算結果と考察

反強磁性母体と Hubbard $U$

まず 30 原子のバルクで $U$ を 3.0〜6.0 eV の範囲で変え、バンドギャップと Fe の磁気モーメントを評価しました (図2)。

図2. バルク α-Fe₂O₃ の (a) バンドギャップ と (b) Fe の磁気モーメント (半径 2.30 bohr の球内で積分) の $U$ 依存性。緑の帯は実験の光学ギャップ (2.0〜2.2 eV [5])、灰色の帯はポーラロンの計算に用いた $U$ の範囲です。$U$ = 3.0 eV の白抜きの点は、高スピンの占有行列から計算を始めたにもかかわらず低スピンの金属状態に収束した結果です。$U$ = 4.0 eV で収束した状態から計算を始めると、高スピンの反強磁性絶縁体が得られます。

バンドギャップは $U$ とともに単調に増加し、$U$ = 4.3 eV で 1.92 eV、5.0 eV で 2.12 eV と、実験の光学ギャップに近い値になります。Fe の磁気モーメントは 3.99〜4.19 μ_B で、いずれの $U$ でも全 Fe が同じ大きさの高スピン状態をとります (球内で積分した値のため、実験値の 4.6〜4.9 μ_B [3, 5] より小さくなります)。$U$ = 4.3 eV は Fe₂O₃ のクラスター計算から自己無撞着に求められた値で、DFT+U によるポーラロンの先行研究 [3] でも用いられています。本解析の主な計算はこの値で行いました。

$U$ = 3.0 eV では、高スピンの占有行列から始めても、Fe のモーメントが 1.1〜2.1 μ_B に減った低スピンの金属状態に収束しました。ところが、$U$ = 4.0 eV で収束した状態から計算を始めると、ギャップ 1.54 eV の高スピン反強磁性絶縁体が得られ、そのエネルギーは低スピン状態より 13.6 eV (Fe 1 個あたり 1.13 eV) も低くなりました。低スピン状態は基底状態ではなく、SCF 計算がエネルギーの高い状態に落ち込んだ結果です。$U$ が小さく d 電子を局在させる効果が弱い場合でも、$U$ の大きい計算で得た高スピン状態から出発すれば、正しい基底状態が得られます。

余剰電子の自発的な自己束縛

次に、余剰電子が全 Fe に広がった非局在状態 (緩和した母体の構造、すべての Fe で $\delta n$ < 0.02) から出発し、そのまま構造緩和を行いました (図3)。格子の歪みや占有行列の偏りなど、局在のきっかけは一切与えていません。

図3. 非局在の余剰電子から出発した構造緩和。(a) 全エネルギー (第 1 ステップ基準)、(b) 電子が局在した Fe の $\delta n$。$U$ = 4.3 eV (青) では約 130 ステップ目、5.0 eV (橙) では約 25 ステップ目に、電子が 1 つの Fe に自発的に局在しました。

$U$ = 4.3 eV では、電子は約 130 ステップの間、いくつかの Fe に部分的に偏る ($\delta n$ は最大 0.17) 状態を行き来したのち、1 つの Fe に急速に局在し、$\delta n$ = 0.27 に達しました。この間に全エネルギーは 73.1 meV 低下しています。$U$ = 5.0 eV では局在がより早く起こり、$\delta n$ = 0.33、エネルギーの低下は 158.7 meV でした。

$U$ = 4.3 eV の終状態は、Fe0 の八面体をあらかじめ膨張させた構造から緩和して得た局在状態と、等価なサイトどうしで比べて 0.06 meV の差で一致しました。局在のきっかけを与えなくても同じ局在状態に到達したことは、非局在状態が余剰電子にとって不安定であり、格子の歪みを伴う局在状態の方が安定であることを直接示しています。この計算は C3 対称性を保ったまま行ったため、電子は 3 回回転軸上の Fe に、C3 軸方向に伸びた a₁ 軌道 (d_z² 型) を占める形で局在しました。この 1 サイト・ポーラロンでは、局在した Fe の Fe–O の平均結合長が 0.044 Å (2.2%) 伸び、面内の隣接 Fe との距離は 2.970 Å から 2.94 Å へわずかに縮んでいます (表2)。

ポーラロンの形：1 サイトか 2 サイトか

C3 対称性を保った計算では、電子は 3 つの隣接 Fe に対して等価な a₁ 軌道に入るしかありません。そこで、対称性の制約を外し、C3 対称性を破った状態を調べました。同じ構造のまま、a₁ 軌道を隣接 Fe (原子 11) の方向に傾けた状態は、C3 対称な状態より 6.4 meV 低くなります。この状態から構造緩和を行うと、$U$ = 4.3 eV では電子が原子 11 にも移り始め、最終的に Fe0 と原子 11 が電子をほぼ等分に分け合う状態 ($\delta n$ = 0.16 / 0.15) に到達しました (図4)。この間に全エネルギーは約 54 meV 低下しています。原子 23 の方向に傾いた状態から出発した独立な計算も、Fe0 と原子 23 の対に同じ形で緩和しました。2 つの終状態の構造は C3 回転で互いに 0.0013 Å 以内で重なり、エネルギーの差は 0.1 meV でした。

図4. 2 サイト・ポーラロンの形成 ($U$ = 4.3 eV、対称性なし)。(a) Fe0 と相手の Fe の $\delta n$、(b) 全エネルギー (原子 11 の方向に傾いた 1 サイト状態の第 1 ステップ基準)。実線は原子 11 の方向に傾いた状態から、破線は原子 23 の方向に傾いた別の軌道状態から出発した計算で、破線の出発点のエネルギーは実線より 35 meV 高くなっています。どちらも、電子を 2 つの Fe で分け合う同じ状態に到達しました。(a) には緩和の前半を示し、最終状態の値 (0.16 / 0.15) を点線で示しました。

$U$ = 5.0 eV で同じ初期状態から緩和した場合も含め、得られた 3 種類の電子ポーラロンを図5と表2で比較します。

図5. Fe0 を含む Fe 2 重層における余剰電子の分布 ($\delta n$)。(a) 1 サイト・ポーラロン (a₁、C3 対称、$U$ = 4.3 eV)、(b) 2 サイト・ポーラロン (Fe–Fe 対、4.3 eV)、(c) 1 サイト・ポーラロン (原子 11 の側に傾く、5.0 eV)。中央が Fe0、真上が原子 11 で、$\delta n$ ≥ 0.03 の Fe に値を、Fe0 と原子 11 の間に距離を示しました。

表2. 電子ポーラロンの電子分布と構造 (母体：Fe–O 平均結合長 2.029 Å、面内 Fe–Fe 距離 2.970 Å)

状態	\(U\) (eV)	\(\delta n\) (Fe0 / 相手)	磁気モーメントの減少 (μB)	Fe–O の伸び (Å)	Fe0–相手の距離 (Å)
1 サイト (a1、C3 対称)	4.3	0.27 / 0.02	0.41 / 0.03	+0.044 / 0.000	2.938
2 サイト (Fe–Fe 対)	4.3	0.16 / 0.15	0.24 / 0.22	+0.028 / +0.026	2.788
1 サイト (相手側に傾く)	5.0	0.28 / 0.06	0.41 / 0.10	+0.046 / +0.011	2.819
文献：HSE (Fock 交換 12%) [1]	—	2 サイト	0.23 / 0.24	+0.026 / +0.028	—
文献：PBE+U [3]	4.3	1 サイト	—	—	2.92 (面内平均 2.99)

磁気モーメントの減少は、半径 2.30 bohr の球内で積分したモーメントが、同じ計算の中の他の Fe (同じスピン副格子の中央値) よりどれだけ小さいかを示します (文献 [1] の値は Hirshfeld 解析によるスピンモーメントの変化)。Fe–O の伸びは、Fe0 と相手の Fe のそれぞれについて、6 本の Fe–O 結合の平均長の母体からの変化を示します。

$U$ = 4.3 eV で得られた状態のうち最も安定な 2 サイト・ポーラロンでは、2 つの Fe のモーメントがそれぞれ 0.24 μ_B と 0.22 μ_B 減少し、Fe–O 結合は平均で +0.028 Å と +0.026 Å 伸び、対をなす 2 つの Fe の距離は 2.970 Å から 2.788 Å へ 6.1% 縮んでいます。これらの値は、Fock 交換を 12% 含むハイブリッド汎関数 (HSE) を用いて同じ大きさのセルで計算された Ahart らの 2 サイト・ポーラロン (モーメントの変化 0.23 / 0.24 μ_B、Fe–O の伸び +0.026 / +0.028 Å) [1] とよく一致します。Ahart らは、隣接サイト間の電子カップリングが再配置エネルギーの半分を上回るため、1 サイトに局在した電子は不安定になり、2 サイトに広がると説明しています [1, 2]。$U$ = 4.3 eV の DFT+U でも、同じ 2 サイト型のポーラロンが得られることがわかります。

一方、$U$ = 5.0 eV では、同じ初期状態から緩和しても 2 サイト型にはならず、電子の大部分が Fe0 に残る 1 サイト型 ($\delta n$ = 0.28 / 0.06) に収束しました。Fe0 の Fe–O 結合は平均で +0.046 Å 伸び、Fe0 と原子 11 の距離は 2.819 Å に縮んでいます。電子の大部分が 1 つの Fe に局在する点では、$U$ = 4.3 eV の DFT+U 計算で報告されている 1 サイト・ポーラロン [3] と同じ型です。ただし、相手の Fe との距離の縮み (5.1%) は文献の値 (面内平均 2.99 Å に対して 2.92 Å、2.3%) より大きく、後で示すように相手の Fe への移動の障壁も非常に小さいことから、2 サイト的な性格を強く残した状態といえます。なお、DFT+U における $U$ の効き方は占有行列を定義する射影子の取り方に依存するため、同じ $U$ の値でも d 電子を局在させる効果は計算コードによって異なります。このように、DFT+U で得られる電子ポーラロンの形は $U$ の値に敏感であり、得られた状態が 1 サイト型と 2 サイト型のどちらであるかを、各 Fe の電子分布で確かめることが重要です。

状態密度

図6に、非局在状態と 2 種類のポーラロンの状態密度を示します。

図6. 状態密度 (価電子帯上端 VBM を基準)。上向きがスピン上向き、下向きがスピン下向きの成分です。点線は VBM と伝導帯下端 (CBM) の位置を示します。(a) 非局在の余剰電子は伝導帯の底に入り、ギャップ内準位は現れません。(b) 1 サイト・ポーラロン、(c) 2 サイト・ポーラロン (いずれも $U$ = 4.3 eV) では、スピン下向きの占有準位がギャップ内に分離します。

ポーラロンの占有準位は、1 サイト・ポーラロンでは CBM の 0.16 eV 下、2 サイト・ポーラロンでは 0.33 eV 下にあり、より安定な 2 サイト・ポーラロンの方が深い準位をつくります ($U$ = 5.0 eV の 1 サイト・ポーラロンでは 0.43 eV 下)。この準位がスピン下向きの成分にだけ現れることは、上向きスピンの Fe³⁺ に入った余剰電子が少数スピンの d 軌道を占めて Fe²⁺ となることに対応しています。

自己束縛エネルギー

各状態の自己束縛エネルギーを表3にまとめます。

表3. 電子ポーラロンの自己束縛エネルギー

状態	\(U\) (eV)	\(E_{\text{ST}}\) (meV)
1 サイト (a1、C3 対称)	4.3	66
2 サイト (Fe–Fe 対)	4.3	126
1 サイト (a1、C3 対称)	5.0	150
1 サイト (相手側に傾く)	5.0	192

$U$ = 4.3 eV で得られた状態のうち最も安定な 2 サイト・ポーラロンの $E_{\text{ST}}$ は 126 meV で、C3 対称な 1 サイト・ポーラロン (66 meV) より 60 meV 大きくなります。この差のうち 6 meV は a₁ 軌道の傾き、54 meV は 2 サイト化による安定化です。$U$ = 5.0 eV では、相手側に傾いた 1 サイト・ポーラロンの $E_{\text{ST}}$ が 192 meV となり、$U$ を大きくするほど束縛は強くなります。なお、図3 の緩和エネルギー (73.1 / 158.7 meV) は、C3 対称性を課した非局在状態を基準にした値です。この基準は対称性を課さない非局在状態より 7 meV ($U$ = 4.3 eV) と 9 meV (5.0 eV) 高いため、図3 の値は表3 の C3 対称な 1 サイト・ポーラロンの値よりその分だけ大きくなっています。また、非局在状態の構造に実験の内部座標を用いた場合との差 (母体の構造緩和の寄与) は 4 meV 程度です。

ホッピング経路

$U$ = 4.3 eV で最も安定な状態は 2 サイト・ポーラロンなので、輸送の素過程は「1 サイトに局在した電子が隣の Fe に跳ぶ」過程ではなく、Fe–Fe 対から別の Fe–Fe 対への移動になります。ここでは、その中で移動距離が最も短い過程として、一方の Fe (Fe0) を保ったまま相手の Fe が原子 11 から原子 23 に入れ替わる移動 (ピボット) を調べました。これは、Ahart らがスピン拘束 DFT で扱った、隣接する 2 サイト状態どうしの間の移動 [2] と同じ過程です。計算では、Fe0–原子 11 の対と Fe0–原子 23 の対の緩和構造の間を線形補間し、各構造で電子状態を求めました (図7a)。比較のため、$U$ = 5.0 eV の 1 サイト・ポーラロンについても、Fe0 から原子 11 への線形補間経路を計算しました (図7b)。こちらの経路は、Fe0 と原子 11 の中点を中心とする反転で両端が互いに移り合うため、補間座標 $x$ と $1-x$ のエネルギーは等しくなります。

図7. 線形補間経路上のエネルギー (上段) と電子の分布 $\delta n$ (下段)。(a) $U$ = 4.3 eV の 2 サイト・ポーラロンのピボット移動 (Fe0–原子 11 の対 → Fe0–原子 23 の対)、(b) $U$ = 5.0 eV の 1 サイト・ポーラロンの Fe0 から原子 11 への移動。灰色の線は 2 状態モデルのフィット、緑の破線はハイブリッド汎関数を用いたスピン拘束 DFT による同じ移動の活性化エネルギー [2] です。(b) の白抜きの点は反転対称性から得られる点、灰色の×は原子 11 側の電子状態から出発して別の軌道状態に収束した計算 (参考) です。

ピボット移動のエネルギーは中点 ($x$ = 0.5) で最大となり、線形補間経路上の障壁は 68 meV でした。中点に近づくにつれて電子は Fe0 に寄り ($\delta n$ = 0.20)、中点で相手の Fe が原子 11 から原子 23 に入れ替わります (図7a 下段)。この障壁は、ハイブリッド汎関数を用いたスピン拘束 DFT による同じ移動の活性化エネルギー 67 meV [2] とよく一致します。始状態と終状態をそれぞれ放物線で表す 2 状態モデル (Marcus 理論) をエネルギー曲線にフィットすると、再配置エネルギーは $\lambda$ = 293 meV (文献 [2] では 363 meV) となりました。同じフィットから得られる電子カップリングは 5 meV と小さく、文献 [2] でもこの移動のカップリングは 26 meV と小さな値です。文献 [2] は、2 つの状態で電子の入る Fe の軌道の向きが 120° 異なることを、カップリングが小さい理由に挙げています。カップリングが小さいことは、この移動が非断熱的な性格を持つことを示しています。

この障壁から、ホッピングの速度と移動度を見積もりました。

$$k = \nu \exp\left(-\dfrac{E_a}{k_{\mathrm{B}} T}\right), \qquad D = \dfrac{z a^2 k}{4}, \qquad \mu = \dfrac{e D}{k_{\mathrm{B}} T}$$

ここで $E_a$ は障壁、$\nu$ = 1.72×10¹³ s⁻¹ は Fe–O 伸縮振動の振動数 (実験値、文献 [2] による)、$a$ = 2.52 Å は 1 回の移動での 2 サイト・ポーラロンの中心の変位、$z$ = 4 は等価な移動の数、$D$ は面内 (2 次元) の拡散係数です。室温 (300 K) では $k$ = 1.2×10¹² s⁻¹、$\mu$ ≈ 0.03 cm²/(V·s) となり、ドープした試料などで報告されている実験の移動度 (およそ 10⁻²〜10⁻¹ cm²/(V·s) [2, 3]) と同程度になります。ただし、この値は線形補間経路の障壁 (真の障壁の上限) と断熱極限 (透過係数 1) を仮定した見積もりです。フィットで得た電子カップリング (5 meV) を用いて Landau–Zener 理論で透過係数を見積もると 0.04 程度となり、非断熱性を考慮した移動度は 1 桁以上小さくなります。文献 [2] では、この移動の電子カップリングを 26 meV (透過係数 0.6) と求め、同じ移動だけによる移動度を 0.020 cm²/(V·s)、中心間の距離が 5.04 Å の移動も含めた電子の移動度を 0.098 cm²/(V·s) と報告しています。

一方、$U$ = 5.0 eV の 1 サイト・ポーラロンでは、Fe0 から原子 11 への移動の障壁は線形補間経路でも 4 meV しかなく、再配置エネルギーも 43 meV と小さな値でした (図7b)。障壁が室温の熱エネルギー (26 meV) より十分小さいため、電子は Fe0 と原子 11 の間をほぼ自由に行き来でき、$U$ = 5.0 eV の状態も実質的には Fe–Fe 対に広がったポーラロンとして振る舞います。$U$ = 4.3 eV の DFT+U 計算による 1 サイト・ポーラロンでは、線形補間による障壁が 0.17 eV と報告されており [3]、本計算の $U$ = 5.0 eV の状態は、この 1 サイト・ポーラロンよりも 2 サイト・ポーラロンに近いといえます。長距離の輸送には、$U$ = 4.3 eV の場合と同様に、Fe–Fe 対の向きが変わる過程が必要になると考えられます。

正孔ポーラロン

正孔についても、Fe0 を Fe⁴⁺ とする局在状態を探索しました (表4)。

表4. 正孔の局在を試みた計算 ($U$ = 4.3 eV)

計算	初期構造と電子状態	\(\Delta E\) (meV)	結果
(a) Fe4+ の初期状態	Fe0 の FeO6 八面体を 2.5% 収縮、Fe0 を d4 配置とする (C3 対称性を保つ計算)	+121	ほぼ非局在の状態 (Fe0 の多数スピン d 占有の減少 0.05)
(b) O− の初期状態	O0 の周囲の Fe–O 結合を 3% 伸長、O0 に正孔を配置	+154	O0 に局在せず (O0 のスピン −0.04)
(c) (a) からの構造緩和	—	−8.9	八面体の収縮が解消し、正孔は全体に非局在
(d) Fe4+ の初期状態、対称性なし	(a) と同じ構造、対称性を課さずに d4 配置とする (C3 対称性を破る配置を含む 5 通り)	+121	(a) と同じ状態 (エネルギー差 0.4 meV 以内)
非局在の参照状態	母体の構造 (実験の内部座標)	0	—

$\Delta E$ は非局在の参照状態に対するエネルギー差です。

Fe0 の八面体を収縮させ、Fe0 の d 電子を 4 個 (d⁴) とした初期状態から C3 対称性を保って計算すると、正孔はほぼ非局在の状態に戻りました。Fe0 の多数スピン d 占有の減少は 0.05 にとどまり、ギャップ内に空準位も現れません。そのエネルギーは非局在の参照状態より 121 meV 高く、正孔は八面体の収縮を支えるほど強く束縛されていません。この構造から緩和すると八面体の収縮は元に戻り、正孔は全体に一様に広がりました (表4(c) のエネルギーの低下は、正孔が非局在のまま結晶の構造が緩和したことによるものです)。対称性を課さずに、正孔を入れる軌道を変えた 5 通りの d⁴ 配置 (二重に縮退した e 軌道の一方だけを空にした、C3 対称性を破る配置を含む) から始めた場合も、すべて同じ状態に戻りました (表4(d))。対称性を課さない構造緩和でも、正孔は局在しませんでした。O の 2p 軌道に正孔を置いた O⁻ 型の初期状態から始めた場合も、O 2p にも $U$ (6 eV) を加えた場合も、正孔は O0 に残りませんでした。

非局在の正孔のスピン密度は、約半分が O 2p 上にあります。Fe 3d だけに $U$ を加える DFT+U では、このように O 2p の性格が強い正孔を Fe 上に束縛できないと考えられ、GGA+U の計算で非局在の正孔より安定な局在正孔が見つからなかったとする報告 [5] とも整合します。一方、ハイブリッド汎関数を用いた Ahart らは、正方晶的に歪んだ Fe の 1 サイトに局在した正孔ポーラロン (スピンモーメントの変化 0.66 μ_B) を報告しています [1, 2]。正孔ポーラロンを記述するには Fock 交換を含む汎関数などが必要と考えられ、ここに Fe 3d への $U$ だけを用いる DFT+U の適用限界があります。

まとめ

第一原理計算ソフトウェア Advance/PHASE の DFT+U 機能を用い、反強磁性絶縁体ヘマタイト α-Fe₂O₃ 中の余剰電子と正孔の局在状態を計算しました。d 軌道の占有行列の初期値を高スピン d⁵ の配置とし、$U$ が小さい場合には $U$ の大きい計算で収束した状態から出発することで、反強磁性の母体を安定に得られました。$U$ = 4.3 eV では、バンドギャップは 1.92 eV、Fe の磁気モーメントは 4.09 μ_B (球内積分) です。非局在の余剰電子から出発した構造緩和では、局在のきっかけを与えなくても、電子が 1 つの Fe に自発的に局在しました (自己束縛)。対称性の制約を外すと、$U$ = 4.3 eV では余剰電子が隣接する 2 つの Fe に分かれた 2 サイト (Fe–Fe 対) ポーラロンが、得られた状態の中で最も安定となり (自己束縛エネルギー 126 meV)、その構造はハイブリッド汎関数による文献の結果とよく一致しました。このポーラロンが一方の Fe を保ったまま相手の Fe を入れ替えて移動する過程の障壁は 68 meV で、ハイブリッド汎関数を用いたスピン拘束 DFT の値 (67 meV) と一致し、断熱極限で見積もった室温の移動度は約 0.03 cm²/(V·s) でした (非断熱性を考慮すると、これより小さくなります)。$U$ = 5.0 eV では電子の大部分が 1 つの Fe に残る 1 サイト型になりますが (自己束縛エネルギー 192 meV)、相手の Fe への移動の障壁は 4 meV しかなく、2 サイト的な性格を残しています。このように、DFT+U で得られる電子ポーラロンの形は $U$ に敏感ですが、占有行列から各 Fe の電子分布を追うことで、どの状態が得られたかを判定できます。一方、正孔は Fe 3d に $U$ = 4.3 eV を加えた計算では、対称性を破った初期状態から出発しても局在せず、ハイブリッド汎関数などが必要になるという DFT+U の適用限界も明らかになりました。

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

お問い合わせ

参考文献

1. C. S. Ahart, J. Blumberger, and K. M. Rosso, "Polaronic structure of excess electrons and holes for a series of bulk iron oxides", Phys. Chem. Chem. Phys. 22, 10699 (2020).
2. C. S. Ahart, K. M. Rosso, and J. Blumberger, "Electron and Hole Mobilities in Bulk Hematite from Spin-Constrained Density Functional Theory", J. Am. Chem. Soc. 144, 4623 (2022).
3. N. Adelstein, J. B. Neaton, M. Asta, and L. C. De Jonghe, "Density functional theory based calculation of small-polaron mobility in hematite", Phys. Rev. B 89, 245115 (2014).
4. N. Iordanova, M. Dupuis, and K. M. Rosso, "Charge transport in metal oxides: A theoretical study of hematite α-Fe₂O₃", J. Chem. Phys. 122, 144305 (2005).
5. J. Lee and S. Han, "Thermodynamics of native point defects in α-Fe₂O₃: an ab initio study", Phys. Chem. Chem. Phys. 15, 18906 (2013).

関連ページ

• 第一原理計算ソフトウェア Advance/PHASE
• 解析分野：ナノ・バイオ
• 産業分野：材料・化学