計算方法

2. 1. 全エネルギー計算

2. 1. 1. 概要

本研究では SiC 表面上におけるグラフェンの成長過程をエネルギー計算によって検討す る．具体的には（I）SiC表面からのSi原子の脱離エネルギーによってSiC表面熱分解にお ける表面モフォロジーの変化を，（II）C 原子の拡散ポテンシャルから SiC表面熱分解過程 で解放されたC原子の拡散過程を（第3 章），（III）SiC表面上に形成する可能性のある複 数のCクラスター構造の形成エネルギーを比較することで，実際に行われるCクラスタリ ング過程を検討する（第 4，5 章）．脱離エネルギー，拡散ポテンシャルおよび形成エネル ギーを計算するためには系の全エネルギーを計算する必要がある．系の全エネルギーEtotは 電子系の全エネルギーEelect（電子運動エネルギー，電子－電子，電子－原子核相互作用の 和）と原子核系の全エネルギーEnucl（原子核－原子核相互作用）との和である（式(2. 1)）¹．

 (2. 1)

式(2. 1)が示すように，系の全エネルギーEtot を計算するためには電子系の全エネルギー

Eelectと原子核系の全エネルギーEnuclを知る必要がある．本研究では，Tight-binding 全エネ

ルギー計算および第一原理計算により系の電子状態を計算し，得られた電子状態から電子 系の全エネルギーEelectを求めた．さらに，適切な原子間斥力ポテンシャルによって原子核 系の寄与Enuclを計算し，式(2. 1)から系の全エネルギーEtotを決定する^{2, 3}．系の全エネルギ ーEtotを基にして，上述の Si 原子脱離エネルギー，C 原子拡散ポテンシャル，C クラスタ ー形成エネルギーを計算し，グラフェン成長過程の議論を進める．

本研究では系の全エネルギーEtotを計算する際に系の構造最適化を行う．構造最適化とは，

一般的に未知となる応力緩和構造を決定する手続きであり，想定する系の最安定状態が求 まる．構造最適化を行うためには系内部の応力，すなわち原子間力を求める必要があるが，

これは系の全エネルギーの勾配から決定する．エネルギー計算を基礎とする一連の計算フ ローを図2. 1にまとめる．

33

図2. 1. 量子力学計算の計算フロー．

2. 1. 2. Tight-binding 全エネルギー計算

2. 1. 2. 1. Tight-binding 近似

半導体結晶中の各原子は閉殻していない価電子軌道を重ね合うことで共有結合を形成す る．このとき価電子軌道よりも内部に位置する閉殻軌道は共有結合の形成に寄与せず，ま た，価電子軌道の波動関数は共有結合形成の前後で変化しないものと考えられる（図2. 2）．

Tight-binding 近似では以上の推察に基づき，結晶中で共有結合をなす電子軌道を結合形成

前の孤立原子状態の価電子軌道（のブロッホ和）を線形結合することによって表現できる と仮定する．



図2. 2. 孤立原子の電子軌道と結晶中の電子軌道．

このような仮定のもとで共有結合性軌道を定式化する ⁴．今，結晶中の位置 rRn にあ る（孤立原子状態の）価電子軌道の波動関数を  = 1，2，3，…，X)とする．このと

34 き価電子軌道のブロッホ和は

Φ 1

√

⋅  (2. 2)

と定義される．式(2. 2)の総和は結晶内に並進配列されたユニットセル；n = 1，2，3，…，

N に関する和である．また，r，Rn，k はそれぞれユニットセル内の価電子軌道の位置，

結晶の並進ベクトルおよび電子状態を指定する波数ベクトルである．ブロッホ和k は結 晶中の電子状態を表す波動関数の基底となる．すなわち，Tight-binding 近似では結晶中で 共有結合をなす電子軌道の波動関数は各価電子軌道のブロッホ和k の線形結合として下 式で与えられる．

Ψ Φ  (2. 3)

ここで，総和は価電子軌道 = 1，2，3，…，X に関する和である．また添え字はバンド 指標であり，波数ベクトルkとともに結晶中の電子状態(, k)を指定する．式(2. 2)，(2. 3) が示すように，Tight-binding 近似では結晶中の共有結合に寄与する電子の波動関数k を 孤立原子の価電子軌道の波動関数  = 1，2，3，…，X)から構成する．ここで価電子軌 道の波動関数は既知であり，式の未知係数 aを求めることで価電子軌道を決定す る．

2. 1. 2. 2. 電子状態の支配方程式

結晶中の電子状態を知るためには，式(2. 3)における線形結合の係数 aを求めればよい．

係数aは1電子のシュレーディンガー方程式；式(2. 4)から求められる⁵．

Ψ Ψ  (2. 4)

ここでハミルトニアン演算子は

2  (2. 5)

である．また，は電子状態(, k)のエネルギー固有値である．式(2. 5)の は電子が結晶中

35

の他の電子や原子核から受ける相互作用を平均化したポテンシャルである．

式(2. 2)，(2. 3)を式(2. 4)に代入すると

1

√

⋅  

1

√

⋅ (2. 6)

式(2. 6)の左から価電子軌道 のブロッホ和 kの共役複素対をかけて全空間積分する．こ のとき左辺は

左辺 

⋅ _{ }

, 

∗

   (2. 7a)

また，右辺は

右辺

⋅ _{ }

∗

 

, 

(2. 7b)

となる．結晶の並進対称性より，式(2. 7a)について

1 _{⋅  }

∗

 

, 

 

^{⋅ } _^∗ _ (2. 8a)

一方，式(2. 7b)は

1 _{⋅  }

∗

 

, 

 

^{⋅ } _^∗ _ (2. 8b)

とすることができる．ここで，飛び移り積分hと重なり積分sを下式で定義する．

36

_{  }^∗ _  (2. 9)

_{  }^∗ _  (2. 10)

このとき式(2. 8a)は

^{⋅ } _^∗ _  

^{⋅ } _{ } (2. 11a)

さらに式(2. 8b)は

^{⋅ } _^∗ _ 

^{⋅ } _{ } (2. 11b)

となる．式(2. 11)を式(2. 6)に代入すれば，1電子のシュレーディンガー方程式；式(2. 4)は

^{⋅ } _{ } 

^{⋅ } _{ } (2. 12)

を得る．式(2. 12)は任意の価電子軌道；' = 1，2，3，…，Xについて成り立つので，未知 数ak ，kの固有方程式である．式(2. 12)を行列表示で表すと

 (2. 13)

ここで，ckはX次のベクトルであり，HkとSkはX次の正方行列である．各々の要素は

(2. 14)

37

, ^⋅  (2. 15)

, ^⋅ (2. 16)

である．このように，1 電子のシュレーディンガー方程式(2. 4)は固有方程式(2. 16)に帰着 する．この固有方程式を永年方程式と呼ぶ．

式(2. 15)，(2. 16)の飛び移り積分行列 H と重なり積分行列 S は原子軌道同士の相互作 用；式(2. 9)，(2. 10)から計算される．Tight-binding 全エネルギー計算では原子軌道同士の 相互作用を直接計算せず，先行研究によって定式化された簡単な関数から求める．したが

って Tight-binding 全エネルギー計算は一般的に計算コストが同じ量子力学計算である第一

原理計算に比べて少ない．求めた飛び移り積分行列Hと重なり積分行列Sから永年方程式 (2. 13)を導き，これを解くことで系の電子状態ak ,，kを求める．

2. 1. 2. 3. 電子系全エネルギーの計算

系の電子状態 ak ,，kが求まと，系の全エネルギーEtotに対する電子系の寄与 Eelectを 計算することができる．今，第1ブリユアン領域（1st B. Z.）内の全ての波数kに対して系 の電子状態 ak ，kが求まったとする．結晶中で共有結合に寄与する価電子は価電子バ ンド（V. B.）を占有するので，電子系全エネルギーEelectは式(2. 17a)によって求まる．

2 1

∈ . .

∈ . .

 (2. 17a)

式(2. 17)右辺の係数2は電子状態ak，kが異なるスピンを持つ2個の電子によって占有 されることに起因する．また第1ブリユアン領域（1st B. Z.）の波数kに対する総和は，結 晶中のすべてのユニットセルに対する足し合わせに相当する．結晶中の電子が第 1 ブリユ アン領域内で取ることのできる波数kの数はユニットセルの個数Nであるので，第1ブリ ユアン領域における波数 k の総和をユニットセル 1 個当たりへ換算するには式(2. 17)に示 すように係数 1/Nをかければよい．実際の計算では第1ブリユアン領域の総和は k空間に おける体積積分によって行われる．また，ユニットセル 1 個当たりへの補正は体積積分値 を第1ブリユアン領域の体積で割ることでなされる．このときの表式を式(2. 17b)に示す．

38 2

2 _{. .}

∈ . .

 (2. 17b)

ここで Vunitはユニットセルの体積，dkは k空間の体積要素である．このとき 2/Vunitが第 1ブリユアン領域の体積となる．

2. 1. 2. 4. NRL Tight-binding 全エネルギー計算スキーム

本研究で検討したSiC表面熱分解における表面モフォロジーの変化と C原子の拡散過程

（第 3 章）では，米国海軍研究試験所（Naval Research Laboratory; NRL）で開発された

NRL Tight-binding 全エネルギー計算を用いた ⁶．本手法では，与えられた計算スキームか

ら求められた系の電子状態ak NRL，k NRLから，原子核系全エネルギーEnuclを計算するこ となく，系の全エネルギーEtot が求まる．すなわち，系の全エネルギーEtot は電子系全エネ

ルギーEelectの表式(2. 17a)と同様，式(2. 18)で与えられる．

2 1

∈ . .

∈ . .

 (2. 18)

このような簡略化は，NRL Tight-binding 全エネルギー計算の計算スキームによって求まる エネルギー固有値k NRL に原子核系全エネルギーを含めることでなされる．すなわち，式 (2. 1)，(2. 17a)による一般的な全エネルギーの表式を以下のように表す．



2 1

∈ . .

∈ . .



2 1

∈ . .

∈ . .



2

∈ . .

∈ . .

(2. 19)

ここで，Neはユニットセル内の全価電子数である．式(2. 19)は固有エネルギーkを固有状

態(k)に関係しない量Enucl /Neで一様にシフトさせているため，NRL Tight-binding全エネ

ルギー計算により求めたエネルギー固有値k NRLは系のバンド構造を再現する．

次に，固有エネルギーk NRLを求めるため用いる飛び移り積分 h；式(2. 9)と重なり積 分s；式(2. 10)のモデリングについて述べる．NRL Tight-binding全エネルギー計算は価電

39

子軌道との相対位置r = (r + Rn)  rに対して，飛び移り積分 h(r)と重なり積分

s(r)を評価するためのモデル式およびパラメータセットを提供する環境依存型

Tight-binding スキームである．本研究では Bernstein らによりバルク SiC 結晶に対して最適化さ

れたパラメータセットを利用した ⁷．このパラメータセットは結晶中の共有結合に寄与す る原子軌道として，Si 原子の 3s，3px，3py，3pz軌道およびC原子の 2s，2px，2py，2pz軌 道を選択しており，これら8種類の原子軌道同士の飛び移り積分と重なり積分を図2. 3に 示す5種類の基本的な原子軌道の重ね方；ss，sp，ps，pp，ppに対して与える．



図 2. 3. 基本的な原子軌道の重ね方．(a) ss，(b) sp，(c) ps，(d) pp，(e) pp

図 2. 3に示した基本的な原子軌道の重ね方に対し，飛び移り積分と重なり積分を下式に よってモデリングする．まず，重ね方 = ss，sp，…の飛び移り積分hについて

40

 

       0 ∩ 

0 0 ∩ 

   exp _ 0

 (2. 20)



ここで r は相対位置r の絶対値である．また， J は原子軌道の周囲に存在する原子種 Jによる局所原子密度であり式(2. 21)で定義される．

_ exp _{ }  (2. 21)

式(2. 22)において，添え字jは原子軌道の周辺に存在する原子種Jの原子を表す．また，

r jは価電子軌道と原子j との間の距離である．式(2. 20)，(2. 21)における関数Fは原子 軌道同士の重なりを打ち切るカットオフ関数であり，式(2. 22)で与えられる．

1

1 exp 5

0

 (2. 22)

ここにRcutはカットオフ半径，l0はカットオフの急峻さを与えるスケーリング長である．

重なり積分のモデリングは飛び移り積分と類似の式(2. 23)で与えられる．

  ^{   exp }

   exp _  (2. 23)

ここで Aと Bはそれぞれ原子軌道，を持つ原子の原子種を表す．また，^は原子軌道 の正規直交性に対応したパラメータであり，式(2. 24)で与えられる．

1 , ,

0 ,  (2. 24)

41

以上の式(2. 20)－(2. 24)がNRL Tight-binding全エネルギー計算における原子軌道間の相 互作用の表式である．これらを用いて実際に原子軌道間相互作用を計算するためには式中 のパラメータを定める必要がある．上述のように，本研究ではBernsteinらが提案するパラ メータセットを用いた．NRL Tight-binding 全エネルギー計算の開発を行った NRL では，

SiC以外の様々な物質に対するパラメータセットを公開している⁸．

2. 1. 2. 5. SCED Tight-binding 全エネルギー計算スキーム

本研究で検討したSiC表面上におけるCクラスタリング過程（第4，5章）では，Yuら によって開発された自己無撞着環境依存型（Self-Consistent Environment Dependent; SCED）

Tight-binding 全エネルギー計算を用いた ⁹．本計算スキームは系の全エネルギーEtotを電子

系全エネルギーと原子核系全エネルギーの和として下式で与える．

2

∈ . .

∈ . .

1 2

1 2

1 2

(2. 25)

ここで第 1-3 項が電子系全エネルギー，第 4 項が原子核系全エネルギーである．またk SCED，Zi，Ni，Rijはそれぞれ固有状態(，k)の固有エネルギー，原子 i の持つ陽子数および 価電子数，原子i，j間の距離である．E0は素電荷eおよび真空の誘電率0から決まるクー ロン定数

4 (2. 26)

である．さらにVNは以下のように定まる．

_, 1 _, (2. 27a)

_{, ,} ∆ _, (2. 27b)

_, 1 1 _,

, (2. 27c)

ドキュメント内 九州大学大学院 (ページ 37-121)