遷移状態を探す時のトラブル

酢酸ビニルと塩化ビニリデンの遷移状態を探す時にちょっとトラブルが...。MOPACで遷移状態を求めるときに良く出くわすトラブルと解決方法をまとめておきます。

MOPAC7は知りませんが、MOPAC97では構造最適化の回数が500回になっています。遷移状態を求めようとしたとき、

MOPAC 7

CYCLE: 739 TIME: 0.00 TIME LEFT: 7096.0 GRAD.: 9.576 HEAT:-63.95119 
CYCLE: 740 TIME: 0.00 TIME LEFT: 7096.0 GRAD.: 8.015 HEAT:-63.97669 
CYCLE: 740 TIME: 0.00 TIME LEFT: 7096.0 GRAD.: 8.347 HEAT:-63.96488 
CYCLE: 741 TIME: 0.00 TIME LEFT: 7096.0 GRAD.: 8.730 HEAT:-63.96032 
TRUST RADIUS NOW LESS THAN 0.00100 OPTIMIZATION TERMINATING
GEOMETRY MAY NOT BE COMPLETELY OPTIMIZED
-------------------------------------------------------------------------------
PM3 PRECISE TS UHF XYZ
H9C6O2Cl2
19

MOPAC 97

CYCLE: 497 TIME: 0.000 TIME LEFT: 55.95M GRAD.: 0.103 HEAT:-64.24002 
CYCLE: 498 TIME: 1.000 TIME LEFT: 55.93M GRAD.: 0.100 HEAT:-64.24003 
CYCLE: 499 TIME: 0.000 TIME LEFT: 55.93M GRAD.: 0.141 HEAT:-64.24004 
CYCLE: 500 TIME: 0.000 TIME LEFT: 55.93M GRAD.: 0.104 HEAT:-64.24005 
EXCESS NUMBER OF OPTIMIZATION CYCLES CURRENT VALUE OF GRADIENT NORM = 0.103911
CURRENT VALUE OF GEOMETRY
PM3 PRECISE TS UHF XYZ

となって(MOPA7では計算が終了するのに）MOPAC97では計算が止まってしまう事があります。このようなときはキーワードに

とかを入れます。MOPACの計算が何日もかかっていた時の名残でしょう。

でも今回のVAc-VDCは遷移状態は求まります。その遷移状態の構造を使って振動計算を行うと、

A FAILURE HAS OCCURRED, TREAT RESULTS WITH CAUTION!! 
++++----**** FAILED TO ACHIEVE SCF. ****----++++ PM3 CALCULATION
MOPAC 97.00
2004/ 8/11


FOR SOME REASON THE SCF CALCULATION FAILED.

THE RESULTS WOULD BE MEANINGLESS, SO WILL NOT BE PRINTED.
TRY TO FIND THE REASON FOR THE FAILURE BY USING 'PL'.

CHECK YOUR GEOMETRY AND ALSO TRY USING SHIFT OR PULAY.

SCF CALCULATION FAILEDエラーが出て止まってしまいます。これは本当に構造が悪くてSCFが求まらない事もありますが、

とSCFの回数を増やしてあげるとこの場合は無事計算が終了してちゃんとした遷移状態が求まります。

0.出発の構造を変えてみる。（この場合は例えばVAc-ANから作ってみる。）
1.内部座標系とXYZ座標系を入れ替えてみる。
2.MOPAC7とMOPAC97を入れ替えて計算してみる。
3.分子の一部を動かして計算してみる。

などを試してみます。 MOPACのTSキーワードは結構強力なキーワードで今回のように遷移状態ライブラリーから出発した場合にはかなり楽に正しい遷移状態が求まります。

NORMAL COORDINATE ANALYSIS

Root No. 1 2 3 4 5 6 7 8

1 A 2 A 3 A 4 A 5 A 6 A 7 A 8 A 

-20.3 17.7 33.1 39.7 85.8 135.9 162.4 194.3

1 0.0017 -0.0013 0.0088 -0.0078 -0.0019 0.0155 0.0023 0.0010
2 -0.0003 0.0127 -0.0091 -0.0013 -0.0127 0.0192 0.0018 0.0167
3 -0.0010 0.0061 0.0041 -0.0126 -0.0111 -0.0509 -0.0120 0.0080
4 0.0017 -0.0015 0.0088 -0.0079 -0.0021 0.0155 0.0019 -0.0001
5 -0.0027 0.0181 -0.0294 -0.0243 -0.0273 0.0237 0.0009 0.0109

例えば虚の振動は一つだけでこのような解が求まる事があります。この虚の振動を見てみるとすぐに判りますが、これは反応方向へは向かっていません。（ねじれの振動）自分の経験上では正しい遷移状態だと虚の振動は-200より小さくなるので、それより大きいものは疑ってかかっています。

この他、虚の振動が複数存在する時などは明らかに計算失敗です。

0.出発の構造を変えてみる。
1.分子の一部を動かして計算してみる。

などは有効な事が多いです。MOPAC7のものと比べMOPAC97のTSキーワードは強力で、まず虚の振動が複数存在することは無いようです。

通常の計算であれば、最初の遷移状態を求めるときには
PM3 PRECISE TS UHF (XYZ)
などというキーワードを使い、求まった構造に対して
PM3 PRECISE FORCE UHF (XYZ)
とやって振動だけを計算します。ごくたまにこのFORCE計算をするときに構造も動いてしまう事があります。計算結果にfilename.arc(FOR012)などができていないかきちんとチェックをしましょう。特にMOPAC7ではそうした傾向が強いようなので、場合によると一度TSを求めた後、XYZと内部座標を入れ替えたりして数回TSを求めてからFORCE計算を行う。などということをやったりします。

自分はついつい振動方向だけをチェックしてIRC計算ははしょってしまうのですが、きちんとした計算は（特に初心者は）大事です。参考書を良く読んで正しい計算をしてください。

最初にPM3 PRECISE TS UHF (XYZ)のキーワードで計算してMOPAC7の場合はFOR012ができたら、その上の部分を消去してTSをFORCEに変えて計算すれば振動計算ができます。つまり入力データに電荷が入っていても問題なく計算できます。しかしMOPAC97ではこの電荷の部分は消去しないと計算できないようです。電荷の部分を一つ一つ消去するのはめんどくさいものです。構造の取り出しでラジカルもモノマーも選択しないでSearchボタンを押すとからっぽのテキストエリアが現れます。そこへ分子の構造の部分だけを（電荷の値ごと）コピーしてReadボタンを押します。その後ですぐにGetStructureボタンを押すと電荷が消去された構造が得られるのでそれをテキストエディターなどにコピーしてTSをFORCEに書き換えてセーブしてMOPACを計算するといいでしょう。TSの構造をチェックする意味でも有効です。

原子の種類を変えるときの注意。原子の種類を変えるときにはその新しい原子の標準結合長を使って残りの部分を動かします。この際遷移状態では多くがそうですがラジカルとモノマーの間に結合が無いためどちらか一方の分子は動かないという現象が出てしまいます。ご面倒でも分子をEditする際はラジカルとモノマーの間に結合を作ってください。（環の一部の原子の種類を変える場合は結合長は変更しません。）

アミドを計算する場合にはMMOKキーワードを入れてください。

マック版のMOPAC97(Chem3Dに付属のもの）は続けて何個も計算していると Not Enough Time to Complete Hessian というメッセージが出て止まってしまいます。その時にはアプリケーションを終了して立ち上げ直します。(単なるバグだと思います。）自分の私見ですがMOPAC97はChem3Dとの連携もおかしいし計算も遅いしSCFの回数などのオフセット値も適当ではないようです。余り使わない方が良さそうです。自分の使っているMOPAC７の最大原子数が４０という制限があるのでそれを超えるものは仕方なくMOPAC97を使っています。そのときのキーワードには、 T=100000 CYCLES=20000 ITRY=2000を加えています。

JAVAの発生させる座標がゼロに近い時、1.23xxxxE-4というように指数表記なってしまう事があります。このままMOPACを走らせるとエラーになるので（Eという原子は存在しないとかいうエラー）手動で0.000123xxxと書き直してください。

ピリカでその他の化学

遷移状態を探す時のトラブル