2次元材料の荷電欠陥の第一原理解析：Cu₂WS₄単層の空孔形成エネルギー#

PHASE logo

2次元材料の点欠陥は、半導体特性や触媒活性を大きく左右する重要な要素です。しかし、荷電欠陥を周期境界条件の第一原理計算で扱うと人工的な静電相互作用がエネルギーに混入し、2次元材料ではその誤差が真空層の厚さに依存して増大するため、2次元系特有の静電補正が不可欠となります。本解析では、第一原理計算ソフトウェアAdvance/PHASEと、2次元補正法（Freysoldt–Neugebauer法）を実装したsxdefectalign2dを連携させ、水素発生反応（HER）触媒として注目されるCu₂WS₄単層のCu空孔（V_Cu）とS空孔（V_S）について、荷電状態ごとの形成エネルギーと遷移準位を評価しました。補正に必要な誘電応答の算出から、真空層の厚さやスーパーセルサイズに対する収束性の検証までを一貫して行っています。

Keywords: 第一原理計算 (DFT), 2次元材料, 点欠陥, 荷電欠陥, 静電補正, Freysoldt–Neugebauer法, sxdefectalign2d, 欠陥形成エネルギー図, 遷移準位, Cu₂WS₄, 水素発生反応 (HER)

荷電欠陥の形成エネルギーと2次元補正#

欠陥の電子的性質を理解する上で基本となるのが、欠陥の荷電状態ごとの形成エネルギーをフェルミ準位の関数として表した「欠陥形成エネルギー図」と、その交点である「熱力学的遷移準位」です。周期境界条件を用いた計算では、欠陥電荷とその周期像、および中性化のための一様背景電荷との間に人工的な静電相互作用が生じます。3次元バルク結晶とは異なり、真空層を挟む2次元スラブモデルではこの誤差が真空層の厚さに比例して増大するため、スーパーセルを大きくして外挿する単純な方法でも正しい値が得られるとは限りません [1, 2]。そのため、スラブ内部と真空領域の誘電率の急峻な変化を考慮した2次元専用の補正手法が必要です。

荷電状態 q の欠陥 D^q の形成エネルギーは次式で与えられます [3, 4]。

$E_\text{f}[\text{D}^q] = E_\text{tot}[\text{D}^q] - E_\text{tot}[\text{host}] - \sum_i n_i \mu_i + q\,(E_\text{VBM} + E_\text{F}) + E_\text{corr}$

ここで $E_\text{tot}$ はスーパーセルの全エネルギー、 $n_i$ は欠陥形成に伴い外部とやり取りした原子の数（系に加えた場合は $+1$ 、取り除いた場合は $-1$ ）、 $\mu_i$ はその化学ポテンシャル、 $E_\text{VBM}$ は完全系の価電子帯上端、 $E_\text{F}$ はVBMから測ったフェルミ準位です。 $E_\text{corr}$ が人工的な静電相互作用を取り除く補正項です。

Freysoldt–Neugebauer法 [1] では、2次元材料を「厚さ $s$ 、誘電率 $\varepsilon$ の誘電体スラブ（外側は真空）」としてモデル化し、孤立欠陥と周期スーパーセルの静電エネルギー差から補正エネルギーを求めます。さらに、DFTで得た静電ポテンシャル差から長距離ポテンシャルを差し引いた短距離ポテンシャルが真空の上下で同じ一定値 $C$ になるようにモデル電荷の位置を調整し、ポテンシャル整合項 $q\,C$ を加えます。スラブの誘電率 $\varepsilon$ と厚さ $s$ は、スーパーセル全体の面内誘電率 $\varepsilon_\parallel$ と面直誘電率 $\varepsilon_\perp$ から次式で決定されます。

$\varepsilon = \dfrac{\varepsilon_\parallel - 1}{1 - \varepsilon_\perp^{-1}}, \qquad s = c\left[\dfrac{1}{1-\varepsilon_\parallel} + \dfrac{1}{1-\varepsilon_\perp^{-1}}\right]^{-1}$

本解析では、補正に必要な面内・面直の誘電応答（電子とイオンの寄与を含む静的な誘電応答）をAdvance/PHASEにより一貫して算出しました（詳細は解析事例「双極子補正と線形応答計算を用いた2次元材料の面内・面直誘電率の算出」を参照）。

計算モデルと計算条件#

Cu₂WS₄は、CuとWがそれぞれS原子の四面体に配位した層状化合物です。その単層（空間群 $P\bar{4}2m$ ）は、WとCuからなる正方格子の上下にS原子が配置された構造を持ちます（図1）。まず単位胞で格子定数を最適化し、PBEの結果から a = 5.49 Å として2×2スーパーセルを作成しました（参考：分散力補正を含む文献 [5] のPBE+D3値は5.43 Å）。次に、スーパーセル中央付近のCu原子、または上側のS原子を1つ欠損させて空孔モデルを構築しました。構造最適化にあたっては、空孔周辺の原子にわずかな初期変位を与えた上で、全原子の緩和を行っています。k点メッシュは3×3×1とし、価電子帯上端および伝導帯下端のk点が含まれるように設定しました。真空層を含むセルの高さ c は 15 Å を基本としましたが、補正の妥当性を検証するため、緩和後の構造を固定して c = 12.5, 17.5, 20 Å とした場合の計算も実施しました。さらに、スーパーセルサイズに対する計算結果の収束性を確認するため、同等のk点密度を用いた3×3スーパーセルでの計算も行っています（詳細は後述）。

Cu2WS4単層の2x2スーパーセルと空孔サイト
図1. Cu₂WS₄単層の2×2スーパーセル。(a) 上面図と、本解析で扱った空孔サイト（破線の円。V_S は上側のS原子）。(b) 側面図（1列分の原子）。

表1. 計算条件の概要

項目	設定
擬ポテンシャル	ウルトラソフト型（生成エンタルピー計算: PAW型）
交換相関汎関数	GGA-PBE
カットオフエネルギー	25 Ry (波動関数)/ 225 Ry (電荷密度)
欠陥モデル	2×2スーパーセル、セル高さ c = 15 Å（検証用に12.5, 17.5, 20 Å）。収束性の検証に3×3スーパーセル（c = 15 Å）
k点サンプリング	2×2：3×3×1、3×3：2×2×1（Γ点を含むメッシュ。どちらも単位胞の6×6×1メッシュに相当し、価電子帯上端（M点）と伝導帯下端（Γ点）を含む）
荷電状態	V_Cu: +1, 0, −1　V_S: +2, +1, 0, −1
スピン分極	電子数が奇数の V_Cu⁰, V_S⁺, V_S⁻ で考慮
構造最適化	全原子（格子定数は固定）、残留力 2×10⁻⁴ hartree/bohr 以下

荷電状態は、スーパーセルの電子数を増減させて計算しました。補正に必要な静電ポテンシャルは、PHASEから3次元分布（Gaussian cube形式）として出力し、補正計算に用いています。

補正計算は、PHASE-Viewerの荷電欠陥補正画面から直感的に実行できます（図2）。静電ポテンシャル、誘電体スラブのパラメータ、欠陥の位置と電荷を指定することでsxdefectalign2d [4] が実行され、補正エネルギーが算出されます。

PHASE-Viewerの荷電欠陥補正パネル
図2. PHASE-Viewerの荷電欠陥補正の画面（2次元材料用の slab (FNV) モード）。V_S²⁺（2×2スーパーセル、c = 15 Å）の設定例です。

なお、形成エネルギーの計算では、スーパーセルの全エネルギー、価電子帯上端のエネルギー $E_\text{VBM}$ 、ポテンシャル整合に用いる静電ポテンシャルのエネルギーの基準をそろえる必要があります。本解析では基準をそろえたうえで、完全系から電子をわずかに取り除いた計算の全エネルギー差（ヤナックの定理）から求めた $E_\text{VBM}$ が、固有値から求めた値と 0.02 eV の差で一致することを確かめました。

計算結果と考察#

基本物性と誘電応答#

PBEで求めた単層Cu₂WS₄のバンド構造は価電子帯上端がM点、伝導帯下端がΓ点にある間接遷移型で、バンドギャップは 1.60 eV です。文献のPBE値 1.60 eV とよく一致します（HSE06では2.49 eV [5]）。以下の欠陥形成エネルギー図は、このPBEのギャップの範囲で示します。

補正に用いる誘電応答を表2にまとめます。イオンの寄与は、面直方向の分極率 $\alpha_\perp$ に対して約1%、面内方向の $\alpha_\parallel$ に対して約7%です。これらを合わせた静的な応答を再現する誘電体スラブの厚さは $s$ = 4.88 Å、誘電率は $\varepsilon$ = 6.89 となりました。

表2. 単層Cu₂WS₄の誘電応答と、補正に用いたスラブモデルのパラメータ

量	電子のみ（イオン固定）	電子＋イオン（イオン緩和）
$\alpha_\parallel$ (Å)	2.129	2.286
$\alpha_\perp$ (Å)	0.328	0.332
$\varepsilon_\parallel$（c = 15 Å）	2.784	2.915
$\varepsilon_\perp$（c = 15 Å）	1.379	1.385
スラブの誘電率 $\varepsilon$	6.48	6.89
スラブの厚さ $s$ (Å)	4.88	4.88

※ 分極率は単位面積あたりの値（長さの次元）。 $\varepsilon_\parallel,\ \varepsilon_\perp$ はセル高さ c = 15 Åのスーパーセル全体の値です。分極率と誘電率の計算方法、収束性の確認は、解析事例「双極子補正と線形応答計算を用いた2次元材料の面内・面直誘電率の算出」にまとめています（本表の値は同事例の表5と同じです）。

ポテンシャル整合：DFTと誘電体モデルの比較#

図3に、V_S²⁺ と V_Cu⁻（2×2スーパーセル、c = 15 Å）について、sxdefectalign2dが出力する平面平均ポテンシャルを示します。欠陥系と完全系のDFTポテンシャルの差 $\Delta V_\text{DFT}$ （橙）は、背景電荷の影響で真空中に放物線状の大きな変化を示しますが、モデル電荷による長距離ポテンシャル $V_\text{model}$ （黒）がこれをよく再現し、差 $V_\text{sr}$ （緑）は真空領域で平坦になります。平坦部の値がポテンシャル整合項 $C$ です。

図3. sxdefectalign2dによるポテンシャル整合（2×2スーパーセル、c = 15 Å）。(a) V_S²⁺、(b) V_Cu⁻。灰色の帯は誘電体スラブ、縦の破線は孤立系として扱った領域の境界です。緑の点線が平坦部の値 $C$ 。

真空層の厚さに対する収束性：補正の有効性#

補正の有効性を検証するため、真空層の厚さ（セル高さ c）を変えて形成エネルギーを比較しました（図4）。補正を行わない場合、形成エネルギーは c の増加とともに直線的に増大し、V_S²⁺ では 1.42 eV もの大きな変動が生じます。これは2次元スラブ計算特有の周期境界による人工的な静電誤差であり、未補正の計算値からは真のエネルギーを評価できません。一方、sxdefectalign2dによる2次元補正を適用すると、どのセル高さでも形成エネルギーは 0.021 eV 以内の高い精度で一定値に収束しました。これにより、真空層の厚さによらない孤立欠陥の形成エネルギーが得られることが確認できました。

真空層の厚さに対する形成エネルギーの収束性
図4. セル高さ c（真空層の厚さ）に対する形成エネルギーの変化（2×2スーパーセル、k点3×3×1。構造は c = 15 Åで最適化したものに固定し、c = 15 Åでの補正後の値を基準とした相対値）。白丸：補正なし、四角：sxdefectalign2d補正あり。

表3. 荷電欠陥の補正エネルギー（2×2スーパーセル、k点3×3×1、c = 15 Å）

欠陥種	荷電状態 q	補正エネルギー $E_\text{corr}$ (eV)	ポテンシャル整合項 $C$ (eV)
V_Cu⁺	+1	0.46	−0.109
V_Cu⁻	−1	0.69	−0.114
V_S²⁺	+2	2.09	−0.101
V_S⁺	+1	0.46	−0.113
V_S⁻	−1	0.69	−0.116

表3に示すように、静電補正エネルギーは欠陥の電荷に応じて 0.46〜2.09 eV と極めて大きく、欠陥の安定性を正しく議論する上で本補正が不可欠であることが分かります。

化学ポテンシャルの安定領域#

空孔の形成エネルギーは、成長やアニール環境における各元素の化学ポテンシャル $\mu_i$ （原子の取り込みやすさ・放出しやすさ）に依存します。PBEで求めた生成エンタルピー（表4）をもとに、Cu₂WS₄単層がCuSやWS₂などの競合相に分解せず安定に存在できる熱力学的安定領域を求め、図5に示します。本解析では、安定領域の両端にあたる以下の2条件を評価しました。

点A（Cu-rich条件）：Cuが最も過剰な環境（金属CuおよびWS₂と共存）。
点B（Cu-poor / S-rich条件）：Cuが最も欠乏しSが豊富な環境（単体SおよびWS₂と共存）。

表4. 母相および競合相の生成エンタルピー（PBE、1式量あたり）

物質	生成エンタルピー $\Delta H$ (eV/f.u.)	備考
Cu₂WS₄（単層）	−3.39	解析対象（母相）
WS₂	−2.51	競合相（2H型バルク）
CuS	−0.41	競合相（コベリン）

化学ポテンシャルの安定領域
図5. Cu₂WS₄単層の化学ポテンシャル安定領域（緑）。各境界線は単体（Cu, S）や競合相（WS₂, CuS）の析出条件を表し、点A（Cu-rich）および点B（Cu-poor）で形成エネルギーを評価しました。

欠陥形成エネルギー図と遷移準位#

図6に、Cu-rich（点A）およびCu-poor（点B）の各条件における欠陥形成エネルギーをフェルミ準位の関数として示します。最も安定な荷電状態を結んだ線（実線）の傾きが変わる点が、熱力学的な電荷の遷移準位（表5）です。

欠陥形成エネルギー図
図6. Cu₂WS₄単層の空孔形成エネルギー vs. フェルミ準位（(a) Cu-rich、(b) Cu-poor）。白地がバンドギャップ（1.60 eV）で、太線が最も安定な荷電状態、細線が各荷電状態、丸印が遷移準位を表します。

表5. 主要な空孔の熱力学的遷移準位（VBM基準, eV）

欠陥種	遷移	遷移準位 (補正後)	補正なし	物理的特性
Cu空孔（V_Cu）	(0/−)	0.74 eV	0.05 eV	深いアクセプター準位（n型キャリアを補償）
S空孔（V_S）	—	バンドギャップ外	—	ギャップ内では中性（V_S⁰）が常に安定

Cu空孔（V_Cu）：
本解析で扱った空孔の中で最も生成しやすい欠陥です。中性空孔の形成エネルギーは 0.85〜1.30 eV と低く、アクセプターとして振る舞います。(0/−) 遷移準位はVBMから 0.74 eV とやや深い位置にあります。
特にCu-poor条件では、フェルミ準位が伝導帯（CBM）に近づくと V_Cu⁻ の形成エネルギーがほぼ 0 eV（−0.01 eV）まで著しく低下します。これは、ドナー添加により電子をドープしてn型化しようとしても、自発的に V_Cu⁻ が生成して電子を打ち消す（フェルミ準位をピニングする「補償中心」として働く）ことを示しています。一方、(0/−) 準位が深いため、熱励起によるp型伝導への寄与は限定的と考えられます。

S空孔（V_S）：
中性空孔の形成エネルギーは 2.43〜2.87 eV と V_Cu に比べて高く、生成しにくい欠陥です。荷電状態への遷移準位はいずれもバンドギャップ外（価電子帯内または伝導帯内）に位置するため、バンドギャップ内では中性状態（V_S⁰）が常に安定となります。このため V_S は、ギャップ内ではドナーとして働きません。

補正を行わない場合との比較#

図7に、補正なしと補正ありの形成エネルギー図を対比して示します（Cu-rich条件）。補正を行わない場合、V_Cu の (0/−) 遷移準位は 0.05 eV の位置に現れ、一見すると室温でも容易に正孔を供給する「浅いアクセプター」であるかのような誤った結論を導いてしまいます。さらに、その値は真空層の厚さに応じて大きく変動します。2次元補正を適用することで初めて、真空層の厚さに依存しない真の準位（0.74 eV：深いアクセプター）が得られます。このように、2次元材料の荷電欠陥解析では、定性的な物性判断を誤らないためにも適切な静電補正が不可欠です。

補正なしと補正ありの比較
図7. 補正なし（左）と補正あり（右）の形成エネルギー図の比較（Cu-rich、2×2スーパーセル、c = 15 Å）。静電補正により、見かけの浅い準位（0.05 eV）から真の深い準位（0.74 eV）へと正しく補正されます。

スーパーセルの大きさに対する収束性：3×3スーパーセルによる検証#

ここまでの結果は、扱いやすい計算規模である2×2スーパーセルによるものです。このセルサイズで十分な精度が得られているかを検証するため、より大きな3×3スーパーセル（完全系63原子）を用いて同一のk点密度・計算条件でサイズ収束性を確認しました（表6、図8）。

表6. スーパーセルサイズによる形成エネルギーと遷移準位の比較（Cu-rich、 $E_\text{F} = E_\text{VBM}$ ）

欠陥種	荷電状態	形成エネルギー (2×2セル / 28原子)	形成エネルギー (3×3セル / 63原子)	差 (3×3 − 2×2)	(0/−) 遷移準位
Cu空孔（V_Cu）	0（中性）	1.30 eV	1.29 eV	0.00 eV	2×2: 0.74 eV 3×3: 0.70 eV (差: −0.04 eV)
Cu空孔（V_Cu）	−1	2.04 eV	1.99 eV	−0.04 eV	2×2: 0.74 eV 3×3: 0.70 eV (差: −0.04 eV)
S空孔（V_S）	0（中性）	2.43 eV	2.37 eV	−0.06 eV	バンドギャップ外

2x2と3x3のスーパーセルの比較
図8. スーパーセルの大きさに対する収束性の検証。(a) 形成エネルギーの差（塗りつぶし：ギャップ内で安定な荷電状態、白抜き：ギャップ内で安定にならない荷電状態。灰色の帯は±0.1 eV）。(b) 3×3スーパーセル（実線）と2×2スーパーセル（破線）の形成エネルギー図の比較。

実用的な小規模セルでの高精度な収束：
バンドギャップ内で安定に存在する荷電状態（V_Cu⁰, V_Cu⁻, V_S⁰）について、形成エネルギーの差は 2×2 と 3×3 でわずか 0.06 eV 以内、熱力学的遷移準位 (0/−) の差も 0.04 eV とよく一致しました。
一般に荷電欠陥の第一原理計算では、セル間の人工的な相互作用を減らすために巨大なスーパーセルが必要とされがちです。しかし、sxdefectalign2dによる適切な2次元補正を適用することで、ギャップ内で安定な荷電状態については、計算負荷の小さい2×2スーパーセルでも十分に収束した結果が得られることを確認しました。これにより、計算コストを大幅に抑えつつ高精度な欠陥解析が可能となります。

※ なお、ギャップ内のどのフェルミ準位でも安定にならない荷電状態（V_Cu⁺ や V_S²⁺ など、正孔や電子が空孔に局在せずバンドの状態に広がる状態）では、形成エネルギーにセルサイズ依存性が見られますが、ギャップ内で安定にならないという結論は変わりません。

まとめ#

本解析では、第一原理計算ソフトウェアAdvance/PHASEと2次元静電補正ツールsxdefectalign2dを連携させ、水素発生反応触媒として注目されるCu₂WS₄単層の空孔形成エネルギーと遷移準位を評価しました。補正を行わない計算では真空層の厚さに応じて形成エネルギーが最大 1.42 eV も変動し、欠陥準位を見誤るリスクがありますが、2次元静電補正を施すことで形成エネルギーが真空層の厚さに依存せずに算出できることを実証しました。その結果、扱った空孔の中ではCu空孔（V_Cu）が最も生成しやすい欠陥であり、(0/−) 遷移準位が 0.74 eV の深いアクセプターとして働くこと、特にCu-poor環境ではn型化を自発的に妨げる補償中心となることが明らかになりました。一方、S空孔（V_S）は形成エネルギーが高く、ギャップ内では常に中性で安定です。さらに、本手法を用いることで、扱いやすい小規模なスーパーセルでも十分に収束した結果が得られるため、計算コストを大幅に抑えつつ信頼性の高い解析が可能です。

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

お問い合わせ

参考文献#

C. Freysoldt and J. Neugebauer, "First-principles calculations for charged defects at surfaces, interfaces, and two-dimensional materials in the presence of electric fields", Phys. Rev. B 97, 205425 (2018).
C. Freysoldt, J. Neugebauer, A. M. Z. Tan, and R. G. Hennig, "Limitations of empirical supercell extrapolation for calculations of point defects in bulk, at surfaces, and in two-dimensional materials", Phys. Rev. B 105, 014103 (2022).
C. Freysoldt, J. Neugebauer, and C. G. Van de Walle, "Fully Ab Initio Finite-Size Corrections for Charged-Defect Supercell Calculations", Phys. Rev. Lett. 102, 016402 (2009).
sxdefectalign2d（SPHInX add-ons）, https://sxrepo.mpie.de/.
T. Ouahrani, R. M. Boufatah, M. Benaissa, Á. Morales-García, M. Badawi, and D. Errandonea, "Effect of intrinsic point defects on the catalytic and electronic properties of Cu₂WS₄ single layer: Ab initio calculations", Phys. Rev. Mater. 7, 025403 (2023).

量	電子のみ（イオン固定）	電子＋イオン（イオン緩和）
\(\alpha_\parallel\) (Å)	2.129	2.286
\(\alpha_\perp\) (Å)	0.328	0.332
\(\varepsilon_\parallel\)（c = 15 Å）	2.784	2.915
\(\varepsilon_\perp\)（c = 15 Å）	1.379	1.385
スラブの誘電率 \(\varepsilon\)	6.48	6.89
スラブの厚さ \(s\) (Å)	4.88	4.88

2次元材料の荷電欠陥の第一原理解析：Cu2WS4単層の空孔形成エネルギー#