やってみよう量子化学計算! WinmostarでやるMOPAC入門
はじめに
いぇーいみんな量子化学計算してるー？してない？してるわけねーだろって？そうだよね。つよいサーバーとgaussianのライセンスいるもんね…
でも精度は低いものの、ノートPCレベルのスペックで量子化学計算ができるフリーソフトがあるのです。その名もWinmostar! このパッケージにはMOPAC6というソフトが含まれており、ボタン一つでぽちっと計算できちゃうのです。というわけでお前も量子化学計算やるんだよ！！！！
ソフトウェアの導入
Winmostar
https://winmostar.com/jp/dlFreeForm.php
このリンクに自分のメールアドレスを入力するとDLアドレスとライセンスが送られてきます。もし大学のアドレスを持っているなら、学生版にチェックして大学のアドレスを入力することで、ちょっとだけいい機能のものを使うことができますよ。MOPACは無償版でも使えます。
Chemsketch
https://www.acdlabs.com/resources/freeware/chemsketch/download.php
これは構造式を描画するソフトウェアです。こちらで構造を作っておいて、Winmostarで読み込むとはかどります。もしChemdrawや別のソフトウェアが持っている場合はそちらでもokです。
さっそくやってみよう！
1.構造式を作ろう！
まずは元となる構造式をChemsketchで書きます。
汚ったねぇなぁ！でも大丈夫！このボタンならね
右上の方にあるこのボタンをぽちっと押すと、
綺麗なポルフィリンになりました。これをFileから名前を付けて保存で、.mol形式にしてしまいましょう。
左上のFile→Saveas→ドロップダウンで.molを選択して保存します。
2.Winmostarで構造を整えよう！
早速Winmostarを起動して、先ほどのmolファイルをドラッグアンドドロップしてみましょう。
おや…？水素原子がありませんね。実は先ほどのmolファイルには水素原子が含まれていませんでした。水素原子を付けるには、chemsketch上でTool→Add Explicit Hydrogenを選ぶまたは、
Winmostar上で編集→水素を付加→すべての原子に付加をしましょう。
水素がついたら、右上にあるほうきマークをクリックして、構造をきれいにします。
これをやっておくと計算の収束が早くなります。
3.計算しよう！
では早速計算してみましょう。
MOP6W70実行を選択して、適当な名前で保存しちゃいましょう。何やら黒い画面に大量の文字が流れていきます。
CYCLEは計算の回数、TIMEは1回にかかった時間、LEFTは計算の残り時間(3600秒以上かかる場合は計算ができなかったものとして勝手に終わります。)、GRADはエネルギー勾配で、この値が小さいほど収束が近いです。最後にHEATは生成熱です。これが最小値になるように計算していきます。
上の方のCYCLE12ではGRADが4桁なのに対し、CYCLE40では0.144と相当小さくなっていますね。HEATも478から242まで下がっています。GRADがある値以下になったとき、構造最適化が終わったと判断して計算が終わります。
4.エネルギーと分子軌道を見てみよう！
計算が終わったら、この図のように分子軌道、電子密度(mgf)を選び、(名前).mgfを開きましょう。
今回はpolという名前を付けたので、pol.mgfを選びます。
左側がエネルギーダイアグラムです。わかりやすいようにeVにチェックしてあります。a.u.はeVの27.2114分の1の値です。量子化学計算では都合のいい値なので、a.u.で結果が出ることがあります。その場合は27.2114倍してくださいね。
さて、左上にHOMO:57と出ていますので、57番の軌道がHOMOであることが分かります。ということは58番がLUMOですね。WinmostarはやさしいのでHOMO-LUMOギャップまで計算してくれています。
早速HOMOを見てみましょう。Selected MOを57番にして、いざDraw!
こんな結果になりました。なんかその…よくわからないですね…。このDrawでは、電子密度が一定のラインを表示してくれるモードなので、電子密度がちょうど0.03の境界線を表示しています。これだと軌道が小さすぎて、よく見えないので、Isosurface Valueを0.01にしてみましょう。
分子全体に共役系が伸びている様子がよくわかりますね！また、表示をmeshからsmoothに変えることで、
きれいな図になりました！せっかくなのでLUMOも見てみましょう。LUMOは58番だったので、Selected MOを58にしていざDraw!
LUMOも無事描くことに成功しました！やったね！ちなみにHOMOとLUMOをよく見比べると、節の数がLUMOの方が多いことに気が付きます。
そのうち量子化学の記事でも触れますが、節(赤と青の境界線のことです)が、HOMOは8か所、LUMOは10か所存在します。一般に節の数が多ければ高いエネルギーを持つということが知られており、ポルフィリンでも確かに成り立っていることが分かりますね。(σh面も節だけどわかる人はわかるよね…)
おわり。
さて手ほどきは済んだのでみんな好きな分子を計算にかけてみよう！なお、このMOPACは半経験的な方法で、Gaussianがとっている非経験的(Ab initio)ではありませんので、変な分子には変な結果を返すことがあります。ご了承くださいまし。
また、精度もそーんなには良くなく、GaussianでB3LYP/3-21Gで計算した結果と比べるとこんな感じです。(マシンスペックの関係で3-21Gですが、論文に出すなら最低でもB3LYP/6-31G*で計算してください。)
MOPACの結果と軌道の形がちょこっとだけ違いますね。エネルギーも
| MOPAC | Gaussian | |
|---|---|---|
| HOMO | -7.793 | -5.335 |
| LUMO | -1.409 | -2.237 |
1eV程度のスケールで違いますね。1eVは化学の世界ではとても大きな違いです。なので、エネルギーを考えるときは参考程度に…。もちろん、Gaussianでもあっている保証はないですけどね。
ちょっと難しい話
現実の系では、溶媒や他の分子との相互作用で軌道のエネルギーが変化することがあります。Gaussianでは溶媒もモデリングして計算することはできますが、これもまた発展途上です。例えば、溶媒を「分極が起きる場」と考えて近似する分極連続場モデル(https://ja.wikipedia.org/wiki/%E5%88%86%E6%A5%B5%E9%80%A3%E7%B6%9A%E4%BD%93%E3%83%A2%E3%83%87%E3%83%AB)もありますが、あくまで溶媒は形を持った分子なので、それ相応の結果しか与えてくれません。また、QM/MM法（https://ja.wikipedia.org/wiki/QM/MM）と呼ばれる、量子化学計算と分子動力学法(剛体近似みたいなものです。)を行い、タンパク質などの巨大分子を頑張って計算したりしていました。
しかし、スーパーコンピューターの登場で、最近はCPUとメモリの暴力で大量の分子を第一原理的に計算をしようとする流れがあります。これにより、例えば固体表面での分子の挙動を計算することで固体触媒の反応機構を理解できたり、溶液中の分子の水素結合の挙動を見ることができるようになったりして、今まで手が付けられなかった分野にも計算化学が手を出せるようになるはずです。
最近、スーパーコンピューター京の後継機として、富岳が開発され、運用が始まっています。（https://japanese.engadget.com/2019/12/03/supercomputer/）富岳三十六景にあやかった名前も粋ですが、富岳がどんな成果を出してくれるかも楽しみですね。
2020.4.18
