うさぎだんご

うさぎだんご、ってなに?

電子状態計算イントロダクション【Wathematica2021アドベントカレンダー19日目】

はじめに

この記事はWathematicaアドベントカレンダー2021の19日目として執筆されたものです。この記事では量子化学計算というものについて軽くお話しすることにしました。これを読んだあなたにはぜひ「量子化学計算ではHF方程式を解いたり、DFTを用いた計算をしたりして電子状態を求め、化学反応や物性の予測をしてるらしい」とシッタカしてほしい、そんな気持ちで書いています。怖いので先に言っておきますが、不正確な記述が多いのでこの文章をあまりアテにしないでください。不正確な文言を電子の海に漂わせる趣味はまるでなく、書く以上は正しいことを書きたいので、誤りを見つけた場合はぜひコメントで優しく教えてください。よろしくお願いいたします。では、最初にモチベーションについてお話ししましょう。

電子状態計算のその前に

電子状態計算の目標

化学者の興味はだいたい以下の3点に絞られると思います。

  1. その化合物はどのような性質(物性)を持つのか
  2. ある化合物と別な化合物は混ぜるとどのような反応をするのか
  3. その化合物はいかなる反応の結果として得られるのか

とくにこのうちの1番や2番が大きなテーマです。1番に現れた物性というワードは物質の持つ物理的性質のことを指します。電気伝導率や熱伝導率、磁性、光学特性、などなど、物理の方程式に登場するけどその方程式からではわからず測定しなければならない物質固有の値のことですね。物性や化学反応は、その物質がどんな物質かによって異なる値を持つので、それを知ることができれば応用上大変よろしいというわけです。上で上げた3については合成屋さんの腕の見せ所のような気がしますが、この領域に関してはケム・インフォマティクスと呼ばれる分野からのアプローチがあるそうです。

さて、では1や2を決めている要素はなんなのか。それは電子です。実は化学はおおむねこいつがラスボスといっても過言ではありません。もちろん原子核が空間的にどのような分布をしているかは重要ですし、化学反応が始まるためには原子核がある程度のエネルギーを有していないといけないのですが、化学反応がどの原子上から始まるのかを決めるのは原子核の都合ではなくその周りを漂う電子の都合なので電子が決めているといっていいのです。……多分。このことについて少し説明しましょう。

第一原理的な視点から

用語でもやもやすると悪いので先に説明しますと、第一原理というのはシュレディンガー方程式のことです。第一原理という語は多分"一番偉い"程度の意味で、量子論におけるもっとも重要な基礎方程式がシュレディンガー方程式なのは間違いないことなので、それで第一原理と呼ばれるわけですね。ラテン語でかっこよくab initio量子化学計算と言ったりします。実験的に得られたパラメータを使わず、シュレディンガー方程式と初期条件、物理定数のみを頼りに計算するとてもシンプルな計算です。 さて、電子の都合で反応が決まると申し上げた件について説明したいわけです。ここでまず考えたいのは電子と原子核の質量比です。電子と原子核の質量はだいたい1000倍ほど違い、原子核のほうがずっと重いです。古典力学に明るくないあなたにはぜひドデカい岩と自分の体をひもで括り付けて動いてみてください。岩はほとんど動かないと思います。何が言いたいこと言うと、電子の状態の変化が起こるとき原子核はほぼ止まって見えるということです。このことを具体的に計算に取り込む場合、原子核は止まっていると仮定する近似を入れます。これをボルン-オッペンハイマー近似(BO近似)といいます。ちなみに、質量は全く違うけど持っている電荷は等しいから及ぼしあうクーロン相互作用は同程度というところがみそだと僕は勝手に思っています。また、原子核同士の反応(いわゆる核反応)は通常身近では起こらないので原子核の化学的な(?)変化は無視してよいでしょう。結局、分子と分子が近づいたとき、近づいた時の振る舞いを決めるのはだいたい電子だ、といえるわけです。ちなみに、超BO近似の化学も存在するらしいです。怖いね。というわけで、原子やそれの集合である分子についてとりあえず第一原理であるシュレディンガー方程式に従って計算すれば、振る舞い、分かっちゃうんじゃない~~? という気持ちになりますね。ところが、多体問題は解析的には解けないので数値計算するしかないわけです。数値計算をより効率的に行うために発展したのが電子状態計算という分野ということになるわけです。

第一原理計算で分かることとは?

数値的にであっても解くべきシュレディンガー方程式から何らかの数値解が得られたとします。この時電子のエネルギーが解として得られます。シュレディンガー方程式はそういう性質の方程式なので。では、得られた電子のエネルギーって何なんでしょう。エネルギーは実数値だから順序があって……? だから何が分かるんでしょう? 実は求めた電子のエネルギーと熱力学の結果を使って反応が自発的に進行するかどうかが分かります。

熱力学は高校でもやったかと思いますが、熱機関の考察から生まれた分野でありながら、マクロな系の平衡状態(マクロにみて変化がなくなった状態のこと)を正確に言い当てる恐ろしい分野です。この熱力学において定義される自由エネルギーという量を、電子のエネルギーや与えられた温度を用いて計算すると、反応前後での自由エネルギーの変化を見ることで反応が自発的に進行するかどうかが分かります。混ぜるだけで勝手に進むかどうかが分かるんですね。すごい。自由エネルギーが減少する方向に反応が進むということが熱力学によって知られているので、計算結果で自由エネルギー変化が負になっていれば反応は自発的に進行するし、逆なら進行しません。このことに関して分子などの化学反応の主体を擬人化して、エネルギーを減らす"ように"変化する、と表現することがよくあります。さらに、化学ではエネルギーが低い状態のことを"安定"と言ったりします。したがって、(自由)エネルギーがより低い場合により有利な生成物になっている可能性が高い、というふうに電子エネルギーは見ることができます。ちなみに、この電子状態計算において原子核は固定されていますが、この原子核の位置を変化させてエネルギーがどう変化することを調べ、極小に向かうように原子核の空間分布を変えていくことで構造最適化などができます。その他物性値を引き出す方法については勉強不足でよく知りません……ごめんなさい。

電子状態計算の超基礎的な話

シュレディンガー方程式の話

書いてるうちに僕の嫌いないわゆる大学の量子化学の講義みたいな感じになってきましたね。この点に関してはあきらめることにして、電子状態計算の超基礎的な話をしましょう。具体的には上でぼかした部分に関して説明します。シュレディンガー方程式の具体的な表式を下に示すと、

 \displaystyle
\hat{H}|\psi> = E|\psi>

となります。ただし、\hat{H}は系のハミルトニアンEは固有エネルギー、|\psi>状態ベクトルです。\hat{H}は線形演算子|\psi>はベクトル、Eは実数値なので、\hat{H}についての固有方程式ですね。知らない人のために簡単な解説をします。状態ベクトルというやつには考えたい系の情報すべてが詰まっています。ここからほしい情報を取り出したいわけです。取り出す操作に対応するのが「演算子をかける」ことです。状態ベクトル演算子をかけると演算子が対応する何らかの実数に化けます。上に示した式の左辺から右辺への変化を見てください。状態ベクトル演算子ハミルトニアンをかけた結果右辺ではハミルトニアンのあった場所にエネルギーがありますね。演算子は物理量(の種類)を表していて、出てきた実数はその物理量の具体的な値を示します。つまり、ある物理量に対応する演算子を作用させることは状態ベクトルから一つの物理量に関する情報を取り出す操作に対応します。

ハミルトニアンは系のエネルギーに対応する演算子なので状態ベクトルに作用させるとその状態のエネルギー固有値が求まります。これは状態ベクトルがある一つに決まっている場合であることに注意してください。実際には我々は系の状態を知らない、というかそれを知りたい。でもエネルギーも分からない、そういう状態です。ちなみに、多電子系では系の状態を知らないのならある意味では系のハミルトニアンもよく分からないという意味不明な状況に陥るはずです。ある電子が感じるポテンシャルはほかのすべての電子から生じるポテンシャルも含むためです。(ここ説明としてめっちゃ怪しい気がします)

リッツの変分原理

でも大丈夫。系の状態をとりあえず適当に決めてやって、それを出発点として系の基底状態(エネルギーが最も低い状態)を得る方法があります。そのことを保証するのがリッツの変分原理です。リッツの変分原理の主張は、知りたい系の基底状態のエネルギー期待値は、適当に決めたどんな状態(試行状態)のエネルギー期待値よりも小さい、というものです。リッツの変分原理があるおかげでめちゃくちゃやったとしても一番低いエネルギー期待値が出てきた試行状態が現実(の基底状態)に一番近いことになります。ちなみに基底状態であるというのが大事なポイントです。電子は光と相互作用してエネルギーを受け取ることがあります。その結果として電子は励起し、全体としては励起状態へと変化します。この励起状態が大事だったりすることもあるんですが、リッツの変分原理からわかるのは基底状態だけ、ということになります。

試行状態はいかにして作るか、ということになりますが、これは状態がベクトルであったことを思い出せば属するベクトル空間の基底をもってきて、それの線形結合で表せばよさそうです。残念ながらコンピューター上では有限個の基底しか扱えないので無限次元のヒルベルト空間の元である状態ベクトルを完全に再現できるわけではありません。が、ここは甘んじて受け入れましょう。実際には状態ベクトルシュレディンガー表現してやって得られた波動関数を適当な関数の組の線形結合で表すということをします。

ハミルトニアンの中身

具体的にハミルトニアンの中身を考えてみましょう。ハミルトニアンは一般に運動エネルギーとポテンシャルエネルギーからなります。運動エネルギーは、ボルン-オッペンハイマー近似から原子核は静止しているとして電子の分だけ考えます。ポテンシャルに関しては電荷の相互作用であるクーロンポテンシャルが効いていて、核-核間ポテンシャル、核-電子間ポテンシャル、電子-電子間ポテンシャルに分かれます。最初の核-核間ポテンシャルはボルン-オッペンハイマー近似の下では定数になるので無視しましょう。核-電子間ポテンシャルも電子の状態が決まれば決まる量なので、最も複雑なのは電子-電子間ポテンシャルです。

ここでお詫びと訂正なのですが、多電子系でのポテンシャルはクーロンポテンシャルだけではありません。交換相関ポテンシャルと呼ばれる純量子力学的なポテンシャルが働きます。これは電子-電子間ポテンシャルで、パウリの排他原理に基づく効果です。電子はフェルミオンと呼ばれる粒子の分類に含まれます。三次元空間ではフェルミオンのほかにはボソンというものがあります。同じ種類の粒子がたくさんある系では粒子同士は区別がつきません。そこで、それぞれの状態の電子のラベルを入れ替えることをしても状態は変化しないように思われます。ボソンはそのように、ラベルの入れ替え=交換に対して変わらないのですが、フェルミオンでは交換すると符号が反転するということが起きます。交換に対して反対称な粒子がフェルミオンというわけです。これはどういうことを示すかというと、同一の状態を持つフェルミオンは2つ以上存在できないということを示します。まさしくパウリの排他原理そのものです。

どうしてこれがポテンシャルになるかというと、電子が持つスピン自由度と関係があります。電子の状態は空間座標とスピンで指定されるわけですが、スピンが同一だった場合空間座標が一致することはできないのでそこに反発があるように見えるわけです。したがって、パウリの排他原理に起因するポテンシャルを交換相関ポテンシャルといい、クーロンポテンシャルと同じく重要なポテンシャルとなっているわけです。

おわりに

突然終わりになります。本当は第一原理のハミルトニアンにハートリー・フォック近似を導入してハートリー・フォック方程式を作ってそれの話をしたり、そのあとに密度汎関数理論の話を書いたりしたかったんですけど、19日目の記事なのにここを書いているのは19日の1時です……。まぁ18日しか執筆してないんですけどね!! 冒頭に述べたシッタカしてほしいことに到達できないまま記事が終了になってしまい大変悲しいので、後日別な記事で補足することにします。というわけで主張をまとめます。

  • 電子状態計算で知りたいことは物性やどんな反応をするかだよ
  • 第一原理というのはシュレディンガー方程式のことで、電子状態計算は実験によって得られたパラメータを使わないab initioな計算なんだ
  • 化学では、系のエネルギーが下がることを安定になると言うよ。電子状態計算の結果として得られたエネルギーを評価する際にはエネルギーの大小関係を見て、より小さいと反応の結果そっちになりそうだなって思うことが多いよ
  • 電子状態計算をサポートする強力な原理がリッツの変分原理だよ。いろんな状態を試した結果一番エネルギーが低いやつが基底状態にちかくなってるはずだよ
  • 多電子系では確かに電子間のクーロンポテンシャルも大変なんだけど、交換相関ポテンシャルっていう古典力学ではよく分からないやつが出てきてそっちも大変なんだよ

以上です。6,000文字に迫る長文となってしまいました。お疲れさまでした。