84
第４章 最急エントロピー勾配量子熱力学に基づく化学吸着
§4-1.
最急エントロピー勾配量子熱力学§4-1-1. 基礎と先行研究[1]
最急エントロピー勾配量子熱力学（steepest-entropy-ascent quantum thermodynamics，
SEAQT）は熱力学アンサンブルに基づく第一原理的な枠組みである．平衡から遠い非平 衡状態であっても取り扱うことができ，力学的な手法（例えば，分子動力学）よりも，
はるかに低い計算コストで非平衡過程の挙動を予測することができる．今日に至るまで，
枠組みの開発，マルチスケールな時間・空間における反応系・非反応系への適用および 実験との比較による検証が続けられてきた[2-25]．GrmelaとBerettaの議論[6,7,26-28]に 基づくと，非平衡の枠組みは，一般に，不可逆な緩和過程と可逆なダイナミクスの結合 である．SEAQT では，不可逆な緩和過程は SEA 原理（the principle of steepest entropy ascent）[2,3,29]に，すなわち，エネルギーや粒子数で制約された，状態空間におけるエ ントロピーの勾配ダイナミクスに従うものと考える．これはMEP原理（the principle of
maximum entropy production）[30]と等価なものである．最近，Liとvon Spakovskyによ
り開発された状態密度法（density of state method）[12]によって，SEAQTの計算可能な 適用範囲は無限次元の状態空間にまで拡大した．その結果，実際的な工学問題のモデリ ングに対して有用なツールとなった．SEAQT の工学応用の一例に，固体酸化物形燃料 電池のカソードにおける酸素およびクロム酸化物還元のモデリングがある[18,19]．この モデリングでは，物質拡散・熱拡散および電気化学反応・化学反応の連成が，マルチス ケールな時間・空間にわたって，SEAQTという一つの枠組みの中で取り扱われている．
この他にも，いくつかの実際的な系へのSEAQTの適用が報告されている[11-17,20]．
§4-1-2. 運動方程式[12]
状態の発展の仕方は運動方程式によって記述される．SEAQT において，状態は密度 行列𝜌̂によって表現され，運動方程式の一般形は，
𝑑𝜌̂
𝑑𝑡 = 1
𝑖ℏ[𝜌̂, 𝐻̂] + 1
𝜏(𝜌̂)𝐷̂(𝜌̂) (4.1)
ここで，右辺第1項は量子リウヴィル方程式（時間依存シュレディンガー方程式）に従 う可逆な状態発展を表し，𝐻̂はハミルトン演算子，ℏはディラック定数である．右辺第2 項はSEA原理から導出される不可逆な状態発展を表し，𝐷̂は散逸演算子，𝜏は緩和時間 である．散逸演算子はヒルベルト空間における制約された勾配から導出されるが，この とき，ヒルベルト空間の幾何学的特徴を表す計量テンソルが特定されている必要がある．
第４章 最急エントロピー勾配量子熱力学に基づく化学吸着
85
特別な場合として，粒子が独立に分布しているdilute Boltzmann gas状態を考えるとき，
状態は1粒子密度行列（1粒子固有状態の基底において対角行列）で表される．さらに，
ヒルベルト空間の異なる次元が一様となるFisher-Rao計量を採用すると，式(4.1)の可逆 項はなくなり，不可逆な緩和過程だけが残る．このような可逆項と不可逆項が連成しな い形式の運動方程式によって，孤立な化学反応系を取り扱う．いま，密度行列は対角で あるから，状態は対角要素で{𝑝𝑘}と表される．これは，系のエネルギー固有値{𝜖𝑘}に対 する確率分布である．すなわち，エネルギー固有値𝜖𝑘をもつ固有状態 |𝜙𝑘⟩が確率𝑝𝑘で 現れるような熱力学的状態について考えている．固有状態のラベル𝑘に代わって，エネ ルギー準位のラベル𝑗を用い，𝜖𝑗での縮重度を𝑛𝑗，𝜖𝑗をもつ固有状態の確率の和を𝑝𝑗とす る．系の熱力学的状態{𝑝𝑗}によって，系の熱力学的状態量であるエネルギーおよびフォ ン・ノイマンエントロピー[31,32]を表すと，それぞれ，
𝐸 = 〈𝑒〉 = ∑ 𝑝𝑗
𝑗
𝜖𝑗 (4.2)
𝑆 = 〈𝑠〉 = ∑ −𝑝𝑗ln (𝑝𝑗 𝑛𝑗)
𝑗
(4.3)
ここで，〈⋯ 〉はアンサンブル平均を意味する．いま孤立系を考えているので，次式の質 量保存則，エネルギー保存則を満足する必要がある．
𝐼 = ∑ 𝑝𝑗
𝑗
= 1 (4.4)
𝐸 = ∑ 𝑝𝑗
𝑗
𝜖𝑗= constant (4.5)
散逸演算子𝐷̂はSEA原理から導出されるのであった．すなわち，状態（{𝑝𝑗}）空間にお いて，保存則による制約を満たす最大のエントロピー勾配の方向に沿って，時々刻々と 状態は発展する．状態空間における確率𝐼，エネルギー𝐸およびエントロピー𝑆の勾配を 𝒈𝐼，𝒈𝐸，𝒈𝑆とすると，𝒈𝐼と𝒈𝐸が張る多様体𝐿(𝒈𝐼, 𝒈𝐸)に垂直な𝒈𝑆の成分𝒈𝑆⊥𝐿(𝒈𝐼,𝒈𝐸)が，そ の状態発展の方向である．図4.1に3次元ではあるが，イメージ図を示す．
第４章 最急エントロピー勾配量子熱力学に基づく化学吸着
86
図4.1 状態発展の方向のイメージ図．
グラム行列式の比の形で，
𝒈𝑆⊥𝐿(𝒈𝐼,𝒈𝐸)=
|
𝒈𝑆 𝒈𝐼 𝒈𝐸
(𝒈𝑆, 𝒈𝐼) (𝒈𝐼, 𝒈𝐼) (𝒈𝐸, 𝒈𝐼) (𝒈𝑆, 𝒈𝐸) (𝒈𝐼, 𝒈𝐸) (𝒈𝐸, 𝒈𝐸)
|
|(𝒈𝐼, 𝒈𝐼) (𝒈𝐸, 𝒈𝐼) (𝒈𝐼, 𝒈𝐸) (𝒈𝐸, 𝒈𝐸)|
(4.6)
ここで，(⋯ , ⋯ )は内積を表す．したがって，Fisher-Rao計量で，dilute Boltzmann gas状 態の孤立系に対する運動方程式は次のように書き下せる[12]．
𝑑𝑝𝑗 𝑑𝑡 =1
𝜏
|
−𝑝𝑗ln (𝑝𝑗
𝑛𝑗) 𝑝𝑗 𝑝𝑗𝜖𝑗
〈𝑠〉 1 〈𝑒〉
〈𝑒𝑠〉 〈𝑒〉 〈𝑒2〉
|
|1 〈𝑒〉
〈𝑒〉 〈𝑒2〉|
(4.7)
ここで，〈𝑒𝑠〉，〈𝑒2〉はアンサンブル平均として，
〈𝑒𝑠〉 = ∑ −𝑝𝑗𝜖𝑗ln (𝑝𝑗 𝑛𝑗)
𝑗
(4.8)
〈𝑒2〉 = ∑ 𝑝𝑗𝜖𝑗2
𝑗
(4.9)
次のような行列式を定義すると，式(4.7)は簡潔な形になる．
第４章 最急エントロピー勾配量子熱力学に基づく化学吸着
87 𝐴1 = | 1 〈𝑒〉
〈𝑒〉 〈𝑒2〉| (4.10)
𝐴2= |〈𝑠〉 〈𝑒〉
〈𝑒𝑠〉 〈𝑒2〉| (4.11)
𝐴3= |〈𝑠〉 1
〈𝑒𝑠〉 〈𝑒〉| (4.12)
𝑑𝑝𝑗 𝑑𝑡 =1
𝜏(−𝑝𝑗ln (𝑝𝑗 𝑛𝑗) − 𝑝𝑗
𝐴2
𝐴1+ 𝑝𝑗𝜖𝑗𝐴3
𝐴1) (4.13)
一般に，緩和時間𝜏は状態の関数であり，運動方程式の解は実時間における時々刻々の 状態発展の経路を与える．一方で，この緩和時間を定数においた運動方程式の解は状態 発展の経路のみを与えて，途中の各状態を通過する実時間についての情報は与えない．
どちらの場合の経路も同一なものであることが，実時間から無次元時間への変数変換で 理解される[12]．すなわち，式(4.13)の𝜏を定数1として解くとき，𝑡は無次元時間である．
以上の運動方程式をODEソルバーで解くと，ある非平衡な初期状態から平衡な終状態 に至るまでの状態発展，すなわち，非平衡な途中状態すべてがわかる．平衡状態は次の カノニカル分布で与えられ，このときエントロピーは最大値をとる．
𝑝𝑗=
𝑛𝑗exp (− 𝜖𝑗 𝑘B𝑇) 𝑍
(4.14)
𝑍 = ∑ 𝑛𝑗exp (− 𝜖𝑗 𝑘B𝑇)
𝑗
(4.15)
ここで，𝑘Bはボルツマン定数，𝑇は平衡温度，𝑍は分配関数である．カノニカル分布以外 のすべての分布は，非平衡状態を与える．非平衡状態においては，部分系の仮平衡状態
（subsystem hypoequilibrium state）という概念が導入される[12]．系の固有状態を M個 のグループに分割する．このグループ一つ一つを部分系といい，部分系のエネルギー構 造を{𝜖𝑗,𝑚}，{𝑛𝑗,𝑚}，（𝑚 = 1,2, ⋯ , M，𝑗 = 1,2, ⋯ , W𝑚）と表す．𝑚は部分系のラベル，𝑗は 部分系内のエネルギー準位のラベルである．各部分系内でカノニカル分布が達成されて いるとき，すなわち，すべての𝑚に対して，
第４章 最急エントロピー勾配量子熱力学に基づく化学吸着
88 𝑝𝑗,𝑚= 𝑃𝑚
𝑛𝑗,𝑚exp (− 𝜖𝑗,𝑚 𝑘B𝑇𝑚) 𝑍𝑚
(4.16)
𝑍𝑚= ∑ 𝑛𝑗,𝑚exp (− 𝜖𝑗,𝑚 𝑘B𝑇𝑚)
𝑗
(4.17)
であるとき，系はM次の仮平衡にあるという．ここで，𝑇𝑚は部分系温度という．𝑃𝑚は 部分系内の確率の総和である．あらゆる非平衡状態は，ある次数の仮平衡状態にある．
最も平衡からかけ離れている状況では，Mは固有状態の数に等しい．式(4.13)の運動方 程式による状態発展において，一度，仮平衡に達した部分系は，その後も仮平衡状態に 保たれることがわかっている[12]．同様に，それぞれに仮平衡にあった2つの部分系同 士が平衡したならば，すなわち，その2つの部分系にわたってカノニカル分布が達成さ れたならば，その後もそれらの間の仮平衡状態は保たれる．そうして最終的に，系全体 での平衡に達する（図4.2）．ここで，時間発展𝑃𝑚(𝑡)および𝑇𝑚(𝑡)は部分系間の確率輸送，
熱輸送を説明している．
図4.2 部分系の仮平衡が進み，全系の平衡に至るイメージ．
SEAQT ではエネルギー固有構造によって系が定義され，エネルギー準位の数と等しい
次元をもつ状態空間で運動方程式を解く．いま取り扱いたいのは化学反応系であるが，
例えば，気体分子の並進運動のエネルギー準位は連続と見れるほどに密であり，状態空 間の次元の数は膨大になる．したがって，実際に計算を実行することは不可能である．
この問題を解消するため，Liとvon Spakovskyは状態密度法（density of state method）を 開発した[12]．この手法は，同様のエネルギー準位にあって同様の確率をもつ状態は，
同様の状態発展を示すという考えに従う．そこで，十分に近傍にあるエネルギー準位を
第４章 最急エントロピー勾配量子熱力学に基づく化学吸着
89
束ねて，束ねられた固有状態の数は縮重度によって考慮する擬似系（pseudosystem）が 導入される．先ず，近傍にあるエネルギー準位{𝜖𝑗}を束ねるための一定のエネルギー間 隔を定める．極めて高いエネルギー準位がもつ確率は極めて小さいので，適切にカット オフエネルギー𝜖cutを設定して，
𝑒𝑖 = 𝜖ground+ 𝑖
R(𝜖cut− 𝜖ground), 𝑖 = 1, ⋯ , R (4.18)
ここで，𝜖groundは基底エネルギー，𝑖はエネルギー区間のラベル，Rは区間の数である．
すなわち，𝑖番目の区間は[𝑒𝑖−1, 𝑒𝑖]であり，区間の長さは∆𝐸 = 𝑒𝑖− 𝑒𝑖−1．この区間ごとに エネルギー準位のラベル𝑗を付して，{𝜖𝑗𝑖}，{𝑛𝑗𝑖}，（𝑖 = 1,2, ⋯ , R，𝑗 = 1,2, ⋯ , W𝑖）と表す と，区間𝑖内で束ねられた固有値の数および平均値は，
𝑁𝑖 = ∑ 𝑛𝑗𝑖
W𝑖
𝑗=1
(4.19)
𝐸𝑖 = 1
𝑁𝑖∑ 𝑛𝑗𝑖𝜖𝑗𝑖
W𝑖
𝑗=1
(4.20)
と計算される．すなわち擬似系のエネルギー準位と縮重度は，{𝐸𝑖}，{𝑁𝑖}，（𝑖 = 1,2, ⋯ , R） となる．擬似系に対する運動方程式は式(4.2)，式(4.3)，式(4.8)～式(4.13)において，𝜖𝑗 → 𝐸𝑖，𝑛𝑗 → 𝑁𝑖，𝑝𝑗 → 𝑝𝑖と読み替えたものである．この擬似系に対する運動方程式と比較す るために，次にオリジナル系の運動方程式について考える．区間𝑖内のエネルギー準位 についての運動方程式を足し合わせると，
𝑑𝑃𝑖 𝑑𝑡 =1
𝜏(〈𝑠〉𝑖− 𝑃𝑖𝐴2
𝐴1+ 〈𝑒〉𝑖𝐴3
𝐴1) (4.21)
ここで，
𝑃𝑖 = ∑ 𝑝𝑗𝑖
W𝑖
𝑗=1
(4.22)
〈𝑠〉𝑖 = ∑ −𝑝𝑗𝑖ln (𝑝𝑗𝑖 𝑛𝑗𝑖)
W𝑖
𝑗=1
(4.23)
〈𝑒〉𝑖 = ∑ 𝑝𝑗𝑖𝜖𝑗𝑖
W𝑖
𝑗=1
(4.24)
第４章 最急エントロピー勾配量子熱力学に基づく化学吸着
90
い ま ， 仮 平 衡 に あ る 部 分 系 に つ い て 考 え て い る と す る ． 次 の 擬 似 連 続 の 条 件
（quasicontinuous condition）
∆𝐸 = 𝑒𝑖− 𝑒𝑖−1 ≪ 𝑘B𝑇𝑚 (4.25)
を満たすとき，オリジナル系の𝑃𝑖と擬似系の𝐸𝑖がもつ確率は，ほとんどの区間𝑖において 等価になることが証明されている[12]．したがって，擬似系をオリジナル系の近似とし て用いることができる．エネルギー構造が状態密度
𝐷(𝜖) = 𝐶0𝜖𝛼−1 (4.26)
で与えられるとき，式(4.19)，式(4.20)の計算は，次の積分によって実行される．
𝑁𝑖 = ∫ 𝐷(𝜖)𝑑𝜖
𝑒𝑖+∆𝐸 𝑒𝑖
= 𝐶0∫ 𝜖𝛼−1𝑑𝜖
𝑒𝑖+∆𝐸 𝑒𝑖
=𝐶0
𝛼 [(𝑒𝑖+ ∆𝐸)𝛼− 𝑒𝑖𝛼] (4.27)
𝐸𝑖 = 1
𝑁𝑖∫ 𝐷(𝜖)𝜖𝑑𝜖
𝑒𝑖+∆𝐸 𝑒𝑖
=𝐶0
𝑁𝑖∫ 𝜖𝛼𝑑𝜖
𝑒𝑖+∆𝐸 𝑒𝑖
= 𝛼
𝛼 + 1
(𝑒𝑖+ ∆𝐸)𝛼+1− 𝑒𝑖𝛼+1
(𝑒𝑖+ ∆𝐸)𝛼− 𝑒𝑖𝛼 (4.28) 以上の手続きをまとめる．まず，系のエネルギー固有構造から擬似系のエネルギー構造 を構築する．その擬似系に初期状態を設定して，運動方程式を解くことにより状態発展 を求める．特に，初期状態を適切な仮平衡状態に設定すれば，状態発展は部分系温度と 部分系確率により特徴づけられる．
第４章 最急エントロピー勾配量子熱力学に基づく化学吸着
91
§4-2.
吸着過程のモデリング§4-2-1. 系の定義
前節では，SEAQTの一般論，特にLiとvon Spakovskyによる定式化[12]について解 説した．本節では，半導体表面における化学吸着過程のSEAQTモデリングを提案する．
本研究では，次の反応メカニズムを対象とする．
S[A] + B(g) ↔ S[C] + D(g) (4.29)
ここで，S[⋯ ]は表面を，⋯ (g)は気体分子を表している．すなわち，分子Bが表面Aに 吸着して，表面Cとなり，このとき分子Dが脱離する．系は，分子Bが吸着する前の
部分系1（反応物）と，分子Bが吸着した後の部分系2（生成物）から構成される．部
分系1は表面Aにおける吸着子A1，A2，A3，…と分子Bを含み，部分系2 は表面C における吸着子C1，C2，C3，…と分子Dを含むように定義する．このとき，表面Aと 表面Cは，吸着子の違い以外は，同様であるものとする．部分系1に含まれる原子の種 類と数は，部分系 2 に含まれる原子の種類と数に当然等しい．図 4.3 に式(4.29)の反応 メカニズムのイメージ図を示す．
図4.3 反応メカニズムのイメージ図．
部分系1および部分系2のエネルギー固有構造を構築し，部分系1のエネルギー準位に 分布する確率と部分系2のエネルギー準位に分布する確率とを比較することで，反応の 進行度を理解することができる．そのために，§4-2-2では，分子および吸着子のエネルギ ー固有構造の構築方法を示し，§4-2-3では，部分系および全系のエネルギー固有構造の構築 方法を示す．系のエネルギー固有構造が決定すれば，初期状態から運動方程式を解くことで，
吸着反応の非平衡発展を求めることができる．