magattacaのブログ

日付以外誤報

雰囲気量子化学入門(前編) ~シュレーディンガー方程式からハートリー・フォック法まで〜

量子化学ってなんだか格好良くて憧れてしまいますよね!で、学生の頃疑問だったのが講義と実践の圧倒的解離。。。

講義ではいつも「シュレーディンガー方程式入門!」「水素原子解いちゃうよ!」で終わってしまうのに、学会や論文では、「ここはDFTでー、B3LYPでー」みたいな謎用語が繰り出される。。。、 「え!何それ??何この飛躍???」となっていました。

f:id:magattaca:20210117225311p:plain

で、数式わからないけど知ったかぶりたい!格好つけたい!というわけでそれっぽい用語(?)をひろってみました。

参考文献はこちら!本棚の奥から出てきた本です。

量子化学 -分子軌道法の理解のために (化学の指針シリーズ)

量子化学 -分子軌道法の理解のために (化学の指針シリーズ)

  • 作者:中嶋 隆人
  • 発売日: 2009/10/30
  • メディア: 単行本(ソフトカバー)

では早速、雰囲気量子化学入門!まずは前編!ハートリー・フォック法についてお勉強!

基本の復習

まず、基本の復習です。とりあえずシュレーディンガー方程式が解ければ、その分子がどんな感じのやつかわかるんだ、と!

f:id:magattaca:20210117225456p:plain

で、「ハミルトニアンが決まるのが大事」ということですが、 どうも「ハミルトニアンエルミート演算子」ということに関連しているらしい。

f:id:magattaca:20210117225517p:plain

固有値実数だから観測量として意味をもつ」、ということでしょうか?

これを踏まえてもう一度定常状態のシュレーディンガー方程式を見返します。こんな感じ?

f:id:magattaca:20210117225536p:plain

・・・エルミートってそんな物理化学的な意味合いにつながってたんですね。 線形代数の格好いい名前だけど、なんだかよくわからないやつくらいにしか思ってませんでした。。。

では、この大事なハミルトニアンをどう導くか?

古典的なハミルトン関数をつくっておいて演算子を使って書き直す」ことで導出できるそうです。 以下のような「量子化の手続き」と呼ばれる対応規則を用いればOK!!簡単!!

f:id:magattaca:20210117225556p:plain

分子のハミルトニアンの式は長いので省略します。(・・・LaTexにもう飽きた)

雰囲気量子化学の流れ

さて、本題。水素原子からDFTへの穴埋めです。

あやふやな雰囲気ですが、キーワードを拾っていくとこんな感じみたいです。

f:id:magattaca:20210118231543p:plain

多粒子問題のシュレーディンガー方程式を解けないので、近似を頑張って1粒子問題のハートリーフォック方程式までもっていった。 でも、どうしても誤差(電子相関)の問題が残った。解決のためにポスト・ハートリーフォック法が考えられたが、計算コストがとても大きくなった。

で、より計算コストの低い解決策が密度汎関数(DFT)で、「波動関数ではなく電子密度から出発する」という根本的な違いがある。 DFTが解くのはシュレーディンガー方程式そのものではなく等価な別のもの。原理的には厳密に電子相関を見積もることができるらしい。

ただDFTにも「汎関数の正確な形がわからない」という問題があり、近似が導入される。現在のDFT計算の多くはコーン・シャム近似に基づいており、 コーン・シャム法では汎関数の運動エネルギー項のためにコーン・シャム軌道を、また交換相関汎関数と呼ばれる項を導入した。 *1

で、この交換相関汎関数として最も有名なものにB3LYPがある。

やった!B3LYPでてきた!

さっぱり意味がわかりませんが、とりあえずこんな感じに追っていけば論文でよく見るアレにたどり着ける!

ハートリー・フォック法までの流れ

では、前半シュレーディンガー方程式〜ハートリー・フォック方程式までの流れをもう少し詳しく追って見ましょう。

こんな感じ。

f:id:magattaca:20210117225716p:plain

ボルン・オッペンハイマー近似と分子軌道

多原子分子のシュレーディンガー方程式は厳密には解けないので近似が必要です。 近似法の一つとして分子軌道法があり、その基礎としてボルン・オッペンハイマー近似(≒断熱近似)があります。 これは「電子の運動に対して原子核の運動を固定させて考えよう」というもので、原子核と電子を分離することで、 「原子核と電子の多粒子問題」を「電子のみに着目した問題」へと簡略化することができます。

「原子マジで重いしもう止めて良くない??」ってやつですね!

「電子のみ」となりましたが、依然として多電子系は3体以上の多体問題なのでさらに近似が必要です。 ここで導入されるのが分子軌道(Molecular orbital, MO)で、「一つの電子の座標だけを含む1電子軌道関数」です。 分子軌道の概念をもちいることで「1電子の問題」にまで近似することができます。

ちなみに、電子の座標には位置の座標だけでなく電子スピンの座標も含まれます。

MOが出てくると実験化学屋でも親しみを感じられますね!光れ!HOMO-LUMO!

パウリの排他原理とスレーター行列式

さて1電子の問題まで近似できましたが、ここでさらに「波動関数パウリの排他原理を満たさなければならない」という問題がのこってます。 パウリの排他原理は「2つの電子の座標を交換すれば、波動関数符号が変わる」というもので、「電子の交換に関して波動関数反対称である」といいます。

分子軌道法の最も簡単な形、ヒュッケル分子軌道法では「n電子系の波動関数をn個の分子軌道の積」として近似します。 ですが、このハートリー積と呼ばれる近似では、電子の座標を交換しても波動関数の符号が変わらずパウリの原理を満たせません。

ここで、スレーターは行列式の性質「任意の2つの行(or 列)を入れ替えると符号が変わる」に着目し、波動関数行列式で表しました(スレーター行列式)。 行列式波動関数 (Slater determinant)は自動的に反対称化されているためパウリの原理を満たすことができます。 *2

f:id:magattaca:20210117225851p:plain

変分法とハートリー・フォック方程式

さて、近似と波動関数の表し方が決まりました。あとは変分法によりハートリー・フォック方程式を得ることができます。

変分原理によると「近似のエネルギーは真のエネルギーと等しいか、それよりも高く」なります。 従って「最良の近似波動関数は近似のエネルギーを最小にするように選べば良い」ということになります。

以上を踏まえて、スレーター行列式について、波動関数を構成する軌道関数に関してエネルギー期待値を最小化すれば最も良い近似式となります。 ラグランジュの未定乗数法を使って規格直交条件のもとで変分法を適用することで導くことができ、これがハートリー・フォック方程式です。 *3

f:id:magattaca:20210117230002p:plain

ハートリー・フォック法関連用語

ハートリー・フォック(HF)方程式の導出までの流れを辿りました。すでに頭の中がぐちゃぐちゃなので、整理を兼ねて関連用語を拾ってみましょう。

  • 独立粒子モデル

一つの電子に着目して、その電子は原子核と他の電子の作る平均場の中を運動している」という描像で、HF法はこの考え方に従っています。 電子一つとする近似で、問題の簡略化に繋がっていますが、同時に「2つ以上の電子がお互いに衝突したり、散乱したりする効果を含まない」という点で誤差の要因となります。*4

HF法では経験的なパラメーターを用いないので、HF法に基づく分子軌道法ab initio分子軌道法非経験的分子軌道法と呼びます。

  • 軌道エネルギーハートリー・フォックエネルギー

HF法のエネルギー(ハートリー・フォックエネルギー)はエネルギーの期待値です。これに対してHF方程式のフォック演算子固有値軌道エネルギーとよばれます。 ヒュッケル法の場合と異なり、「HFエネルギー ≠ 全軌道エネルギーの和」であることに注意です。

  • 空間軌道スピン軌道

分子軌道の電子の座標には位置の座標電子スピンの座標があります。 前者による空間部分を空間軌道とよび、空間軌道にスピン関数をかけたものをスピン軌道といいます。 スピン関数は規格直交条件をつかって積分できるので、スピン軌道によるHF方程式から、空間軌道だけのHF方程式を導くこともできます。

f:id:magattaca:20210117230223p:plain

  • RHF法とUHF

制限付き(Restricted)HF法と非制限(Unrestricted)HF法の違いはスピンの取り扱いによります。
RHF法では「 αスピンの電子とβスピンの電子をペアにして、同じ空間軌道につめる」ことで得たスピン軌道を使います。 対して、UHF法では「ペアにせず、異なる空間軌道**につめる」ことで得られたものを使います。

ハートリー方程式を解くには?

さて、シュレーディンガー方程式の近似式であるハートリー・フォック方程式を得ることができました。これを解くにはどうすれば良いでしょうか?

ここまでで得たHF方程式は、フォック演算子の運動エネルギー項に微分を、クーロン演算子と交換演算子積分を含む微積分方程式であり、 多原子分子のような多中心の場合、数値的に解くことは困難です。 そこで行列形式ハートリー・フォック・ローターン方程式に誘導し、self-consistent field (SCF)の手続きと呼ばれる手法により解く、ということが行われます。

この行列問題への置き換えはDFT計算においても主流の解き方となっているそうです。

ハートリー・フォック・ローターン方程式

ハートリー・フォック・ローターン方程式を導くには、まず基底関数を用いて分子軌道をLCAO (Linear Combination of Atomic Orbital)展開します。

LCAO展開は、分子軌道を「N個の原子軌道の線型結合」で表すもので、「分子の軌道は原子軌道の形を色濃く残している」といったタイトバインディングの考えに基づく近似です。 この時の原子軌道を基底関数、展開係数を分子軌道係数と呼びます。

HF方程式導出の時と同様、変分原理に従って、近似エネルギーが最小となるような近似式が最良の波動関数となります。 LCAO展開では分子軌道が「線型結合で近似した試行関数」となっているので、リッツ (Ritz)の変分法を用い、 変分パラメーターである展開係数(分子軌道係数)を求めれば良いことになります。

Ritzの変分法では永年方程式と呼ばれる方程式を解くことで展開係数を得ます。永年方程式は行列式なので、ここで行列問題への置き換えができていることになります。

f:id:magattaca:20210117230445p:plain

ちなみにローターンのつづりはRoothaanです。ジョジョの擬音っぽい。ウリィィィィィ!!

self-consistent field (SCF)の手続き

さてハートリー・フォック・ローターン方程式に導くことで行列問題となりました。

ハートリー・フォック・ローターン方程式は「重なり行列をもった一般化固有値問題」の形を取っています。

単純な場合であれば、重なり行列(S)の単位行列への変換(直交化*5により、「重なり行列をもたない普通の固有値問題」へと変換することで、簡単に計算機で解くことができるようになります。

ですが、ハートリー・フォック・ローターン方程式の場合、フォック行列に分子軌道係数からなる密度行列が含まれているため単純ではなく、 「繰り返しの手続き」を用いて計算することになります。この手続きをself-consistent field (SCF)の手続き(自己無撞着場の手続き、つじつまのあった場の手続き)と呼びます。

f:id:magattaca:20210117230327p:plain

上図の手続きをn回繰り返し、n回目と(n-1)回目の係数行列 or 密度行列がほとんど同じものになれば、SCFが収束(converge)したことになり、 ハートリー・フォック・ローターン方程式が解けたことになります。

基底関数

無事ハートリー・フォック・ローターン方程式の解き方までわかりました。 話を戻して、先のLCAO展開で導入された基底関数の中身について触れておきます。これは後編のDFTにも関係してきます。

基底関数の形

まず、基底関数の形としては大きくスレーター型(Slater-type)とガウス(Gaussian-type)があります。 スレーター型は水素原子の波動関数と同じ形で、より正確ですが計算コストが高く、現在では半経験的分子軌道法計算以外ではあまり用いられていないそうです。 対するガウス型は正確さに劣るものの、2電子積分を解析的に求めることができ、計算コストでメリットがあるため、基底として多用されているようです。

f:id:magattaca:20210117230522p:plain

上図で式の記号はそれぞれ、主量子数(n)、方位量子数(l)、磁気量子数(m)、球面調和関数(Y)を表します。 また、R_Aは原子Aの座標ベクトルで、軌道指数(ζ)は関数の広がりを表します。

カスプ条件Wikipedia 加藤の定理)とは、「電子密度が原子核の位置において尖点(カスプ)を有することを述べる」ものです。 スレーター型、ガウス型、両者のグラフを見比べれば、前者はx=0で尖っていてカスプ条件を満たし、後者は滑らか(微分が0)で条件を満たさないことがわかります。

短縮ガウス型基底関数

さて、多用されるガウス型軌道(GTO)ですが、より少ない個数の関数系で精度の高い計算を実現するため、短縮ガウス型基底という形で用いられます。

以下のような形を取っており、軌道指数の異なる原子ガウス型関数の線型結合で与えられます。 「GTO複数組み合わせることで近似的により正確なSTOに近づける」ことが目指されています。

f:id:magattaca:20210117230552p:plain

代表的な基底関数系として以下のものがあります。

f:id:magattaca:20210117230613p:plain

いよいよ6-31Gがでてきました!よく見るけどよくわからないやつ!

格好いいのでもう少し用語を拾います。

pople型基底関数系

こういうX-YZGといった形で表されるsplit-valence基底系は、ジョン・ポープル(Pople)らのグループにより提案されたのでPople型基底関数系とも呼ばれるようです。
ハイフンの前半は内殻軌道、後半は原子価軌道(valence ornital)の表し方を示し、後半の数によって呼び方が変わります。*6
例えば、6-31GのようなX-YZGはダブルゼータ(DZ)型、6-311GのようなX-YZWGであればトリプルゼータ(TZ)型と呼ばれます。 *7

具体例として、炭素原子では6-31Gは以下に対応します。*8

  • 内殻軌道 = 1s → 個の原始(primitive)GTOからなる短縮(contracted)GTO(CGTO)で表す
  • 原子価軌道 = 2s2p → 3個の原始GTOを短縮した短縮GTOと、1個の原始GTO線型結合で表す

原子価軌道は短縮GTOと原始GTOの組み合わせとなっていますが、グラフ形状をみると後者の原始GTOのほうが動径方向に広がっており、 前者が軌道のより内側、後者がより外側も表現する役割を担っているようです。*9

f:id:magattaca:20210117230710p:plain

Pople型以外でよく使われる基底系としては、Huzinaga-Dunning系、Roos ANO系、Dunning cc-pV nZ系といったものがあるそうです。 最後のダニング(Dunning)らによる「Correlation-Consistent(電子相関と整合な)基底系」(ex. cc-pVDZ)も論文でよく見かける気がします。

分極関数、分散関数

さて、6-31Gがわかりました。では、変化形の6-31G(d)6-31+G(d)とは???

それぞれ、「(d)分極関数(polarization function)」、「+分散関数(diffuse function)」を表すそうです。

分極関数は「一つ多くをもつ補助関数角運動量が高い)」で、これを追加することで基底系に柔軟性が加わり、 分子中の原子の電子密度の分極を表現できるようになります。

例えば水素原子であれば、球対称な1sに分極p関数を追加することで、非対称な分子軌道も表せる、といったことになります。

分極関数の表記は追加した軌道をそのまま書くか、アスタリスクをつけます。*10

一方、分散関数(広がった関数)は「非常に平たいガウス関数」で、原子軌道の原子から遠く離れた尾部を正確に表現することができるようになります。 例えば、アニオンや励起状態の電子といった、より広がった電子分布を表現する際に有効な補助関数となるそうです。*11

分散関数の表記は「Gの前にプラス(+)」をつけます。「++」だとより軽い元素(HとHe)にも分散関数が追加されることを示します。

・・・表記1文字にそんな違いがあったとは。。。知りませんでした。

以上、基底関数系の用語紹介でした。

前編まとめ

以上、雰囲気量子化学入門(前編)では、DFTの知ったかぶりに向けてシュレーディンガー方程式の復習からハートリー・フォック方程式とその解き方までを辿って見ました。

ざっくりとまとめると、多体問題解くために、近似して、行列表現にしたよ!って感じでしょうか?それにしてもどうしてこんな発想がでてくるのでしょう?本当すごいですね。

最後に、分子軌道を表す要素として基底関数系についても触れました。よくみる格好いいアレの正体がちょっとわかって楽しかったです。

っていうか脳みそシンプルなので分極関数の形を見て、「え、これ混成軌道と同じじゃないの?そんなちょいちょい足さんでも、もう混ぜてしまえばええんやないの?」って思ってしまいました。たぶん正確さ(定量性?)とか違うんですよね???そもそも分子軌道は1電子の近似だから、化学結合原子価結合法とは別物なのでしょうか?さっぱりわからない。

あとPople型でゼータと呼ぶのがなぜかもわかりませんでした。唯一分かったのはエルミートには格好いいだけじゃない意味があったということ!

格好つけるために数式をLaTeXでコピペしてみましたが、意味はわからなかった!*12

ひどい誤解がたくさんあると思うのでご指摘いただければ幸いです。

参考資料

この記事は主に「量子化学 ー分子軌道法の理解のためにー 中島隆人 著」(裳華房)の内容を参考にしています。 かつて優秀な研究室の先輩にDFTについて知りたいといったら、「とても簡単にわかりやすく紹介されているから」とオススメされました。 私の能力では字面を追うだけで精一杯でしたが、数式の変換も記載してあり、コンパクトに幅広い情報がまとまっているので、 入門書としてとてもよいのではないかな?と思います。ご興味を持たれた方はぜひご一読を!そして、私に解説してください。お願いします。

また、インターネット上に公開してくださっている以下の資料がわかりやすく、参考にさせていただきました。

おまけ

念のため、観測量に関連して「演算子Aの期待値」の定義を復習します。ついでに記号が似てるのでブラケット表現も。 だいたいこんな感じ。

f:id:magattaca:20210117233249p:plain

*1: コメントでご指摘いただいた内容をもとに修正してみました 01/18 update

*2:スレーター行列式には「反対スピンを持つ電子が同じ軌道を占めることが可能になる」という欠点があり、1電子軌道関数に反します。 同じ向きのスピンを持つ電子については、行列式の「同じ成分を持つ行(or 列)が2つあると値が0になる」という性質から、同じ軌道を占めることができません。

*3:以下の図では正準(cononical)ハートリー・フォック方程式としています。これはHF方程式をラグランジュの未定乗数法で解く際に、未定乗数が対角行列となるようにユニタリ変換したものです。正準HF方程式はフォック演算子固有値方程式であり、固有値軌道エネルギー(εi)、対する固有関数がスピン軌道(φi)となります。「HF方程式の解をユニタリ変換してもHF方程式の解となる」こと、「フォック演算子がエルミート演算子であるため、対角化するようなユニタリ変換を選ぶことができる」ことから、このような選び方が可能になります。ここでもエルミートが活躍!(参考 Wikipedia-ハートリー=フォック方程式

*4:後編で出てくる電子相関の問題です

*5:重なり行列はエルミート行列なので、ユニタリー変換により単位行列に変換することができます

*6:この数は「異なる軌道指数を持つ関数のセットの数」を表すそうです。 参考 Chem-Station 計算科学:基底関数って何?

*7:Wikipediaに「ゼータ ζ はSTO基底の指数部を表わすために用いられることが多い」と書かれていますが、なぜそれでPople型の呼称に使われるのかはよくわかりません。どなたか教えてください

*8:使われている原子GTOの数をまとめると、1sに4個、2Sに4個、2pに4個なので、結局s軌道に10個、p軌道に4個使われていることになります。
また、短縮後の関数は1sに1個、2s2pに2個なので、結局s軌道に3個、p軌道に2個使われています。 これらをまとめて「(10s4p)→[3s2p]」といった表記もなされるそうです。

*9:参考 1. 計算科学.com Pople型基底関数
参考 2. Wikipedia 基底関数系(化学))
参考 3. PC CHEM BASIS.COM B3LYP/6-31G(d)とは何なのか?

*10:軌道をそのまま書くのとアスタリスクで微妙に意味が違うようなのですがよくわかりませんでした。

*11:参考 PC CHEM BASIS.COM 基底関数はどれを選べばいいの?

*12:LaTeXiTとかいうのすごい便利ですよ!