|
FireflyはGAMESS(US)の流れを汲む非経験的分子軌道計算プログラムで、PC
GAMESSを改称してFireflyと名づけられた。Fireflyの公式ホームページは こちらです。 質量は、周波数と振動の通常のモードだけに影響を及ぼします。 ・線形分子（１８０度の結合角）・ ・原子を固定する・ 1）デカルト座標で固定する $STATPTグループでIFREEZ(1)=1,2,3・・・$ENDで固定するデカルト座標を指定する。 対称性を使用する場合、COORD=CARTを指定すると対称操作が機能しないのでCOORD=UNIQUEを指定する。また IFREEZ(1)=・・・で指定する座標は対称操作によって生成される内部座標を指定する。したがって、インプットは２原子でも対称操作によって４原子の座標が生成され、そのうち3番目と4番目の原子を固定したい場合はIFREEZ(1)=7,8,9,10,11,12と指定する。どのような内部座標が生成するかは$CONTRLに、EXETYP=CHECKを指定して計算を流すことで確認できる。 構造はデカルト座標で入力するが、原子間距離、結合角度、二面体角など構造の一部を固定し他の構造を最適化したい場合は、＄CONTRLグループでNZVARで内部座標の数を指定。$ZMATグループでIFZMAT(1)、FVALUE(1)で固定する内部座標と値を指定する。 内部座標の指定は、例えば1,1,2は原子１と２の間の結合距離を、2,2,3,4は原子2,3,4で構成される結合角を、3,1,2,3,4は原子1,2,3,4で構成される二面体角を表す。先頭の１,２,３はそれぞれ結合距離、結合角、二面体角を表し、次の数字は原子の番号を表す。（この方法は「すぐできる量子科学計算」の128,129ページ参照） インプットの例 $CONTRL ICHARG=-1 MULT=1 SCFTYP=RHF RUNTYP=OPTIMIZE MAXIT=200 NZVAR=39 DFTTYP=B3LYP $END $SYSTEM TIMLIM=600000 MWORDS=10 $END $STATPT NSTEP=100 OPTTOL=0.0001 $END $SCF DIRSCF=.T. DAMP=.T. NPUNCH=0 $END $BASIS GBASIS=N31 NGAUSS=6 NDFUNC=1 DIFFSP=.T. $END $GUESS GUESS=HUCKEL $END $ZMAT DLC=.T. AUTO=.T. IFZMAT(1)=3,5,6,7,9 FVALUE(1)=10 $END $PCM SOLVNT=H2O $END $DATA Winmostar C1 O 8.0 -2.359744 -1.284661 -0.627286 O 8.0 2.922621 -0.893244 -0.425420 ・・・・ 上の例では原子5,6,7,9の二面対角が１０度に固定され他の構造は最適化されます。 2) Zマトリックスで固定する $CONTRLグループで $COORD= ZMTMPC NZVAR= (内部座標の数）を指定。 $STATPTグループでIFREEZ(1)=1,2,3・・・$ENDで固定するZマトリックスを指定する。 インプットの例 $CONTRL ICHARG=-1 MULT=1 SCFTYP=RHF RUNTYP=OPTIMIZE COORD=ZMTMPC MAXIT=200 NZVAR=39 DFTTYP=B3LYP $END $SYSTEM TIMLIM=600000 MWORDS=10 $END $STATPT NSTEP=100 OPTTOL=0.0001 IFREEZ(1)=21 $END $SCF DIRSCF=.T. DAMP=.T. NPUNCH=0 $END $BASIS GBASIS=N31 NGAUSS=6 NDFUNC=1 DIFFSP=.T. $END $GUESS GUESS=HUCKEL $END $PCM SOLVNT=H2O $END $DATA Winmostar C1 O 0 0 0 0 0 0 0 0 0 O 4.552125 1 0 0 0 0 1 0 0 O 2.253671 1 98.83861 1 0 0 1 2 0 O 1.506343 1 101.25025 1 3.82766 1 2 1 3 N 2.58399 1 97.03605 1 8.91471 1 3 1 2 C 1.49697 1 66.07976 1 -7.17025 1 5 3 1 C 1.531562 1 108.2754 1 131.11899 1 6 5 3 C 1.253044 1 72.38646 1 2.5658 1 1 2 3 C 1.48973 1 108.421 1 132 1 7 6 5 ・・・・ 上の例では２１番目の内部座標である９番目の炭素の二面対角１３２°が固定されます。 ・外部の基底系を使用する・ EMSL Basis Set
Exchangeのホームページで入手する。URLは https://bse.pnl.gov/bse/portal
Fireflyと完全な互換性のあるのは“GAMESS-US”フォーマットのもの 入手した基底系をメモ帳などに貼り付け.LIBの拡張子を持つ名前でインプットファイルと同じディレクトリに保存する。デフォルトではBASIS.LIB。 貼り付けたデータの一部を変更する。例えば元のデータが・・・ $DATA HYDROGEN S 2 1 5.4471780 0.1562850 2 0.8245470 0.9046910 S 1 1 0.1831920 1.0000000 OXYGEN L 3 1 8.5190000 -0.1455100 0.1100700 2 2.0730000 0.0828600 0.3496900 3 0.6471000 0.7432500 0.4809300 L 1 1 0.2000000 0.2847200 0.3072700 であるなら、最初を空行にして、元素名を短縮したものに変え、スペースの後に基底系を識別するためにname(8文字以下）を追加する。 （空行） H SBKJC(←これがname) S 2 1 5.4471780 0.1562850 2 0.8245470 0.9046910 S 1 1 0.1831920 1.0000000 (一行空ける） O SBKJC(←各元素ごとにnameを記入する） L 3 1 8.5190000 -0.1455100 0.1100700 2 2.0730000 0.0828600 0.3496900 3 0.6471000 0.7432500 0.4809300 L 1 1 0.2000000 0.2847200 0.3072700 インプット情報で $BASIS EXTFIL= .T. GBASIS=SBKJC $END のように .T. で外部ファイルの使用を指定し、GBASIS= の次に自分がつけたnameを記入する。それぞれの基底系のために個別のファイルを作ることもできるし、同じファイルの中に幾つかの基底系を含ませることもできる（nameで識別される）。 ・単位について・ Fireflyで出力されるエネルギーの値はハートリーの単位で、これを他の単位に換算する。 ・有効内殻ポテンシャル（ECP)について・
１ハートリー＝ 27.2114 eV ＝2625.4919 kj/mol (1eV ＝ 96.485kj/mol) ＝627.5112 kcal/mol (1eV = 23.0606kcal/mol) 電子数の多い原子の場合、すべての原子軌道を記述するには多くの基底関数が必要だがECPは内殻の電子については価電子に対する有効ポテンシャルに置換えることにより基底関数の数を激減させることができる。また重原子になると重要になる相対論効果を取り入れたものが提供されている。 FireflyにはSBKJC（Li～Rn）やHW（Na～Xe）がセットされていて、例えばSBKJCを使いたい時は $CONTRL ECP=SBKJC $END と指定するだけでいい。合わせて $BASIS GBASIS=SBKJC と指定する。ただしこの指定はすべての原子に適用される。 特定の原子にのみＥＣＰを適用したい時は$CONTRL ECP=READとして$ECPグループを追加する。$ECPグループにはすべての原子を含める必要がある。$DATAファイルに存在する原子の順番通りにそれぞれの原子のECPパラメーターを入力する。同じ種類の原子については最初にものにパラメーターを入力し、後のものについては（元素記号）-ECPと記入すれば最初のパラメーターが適用される。ECPを適用しない原子に対しては（元素記号）-ECP NONEと記入する。 例えばCH2O2の炭素原子と酸素原子にECPを使いたい時は下記の通りになる。 $ECP C-ECP GEN 2 1 1 ----- CARBON U(P) ----- -0.89371 1 8.56468 2 ----- CARBON U(S)-U(P) ----- 1.92926 0 2.81497 14.88199 2 8.11296 H-ECP NONE O-ECP GEN 2 1 1 ----- OXYGEN U(P) ----- -0.92550 1 16.11718 2 ----- OXYGEN U(S)-U(P) ----- 1.96069 0 5.05348 29.13442 2 15.95333 O-ECP H-ECP $END ・基底関数重ね合わせ誤差（BSSE)・ 「基底関数重ね合わせ誤差は,A‐B相互作用系の量子化学計算において不十分な基底関数を適用したときに,A,Bが互いに相手の基底関数を使って過度に安定化してしまうことにより生じる誤差です。・・・結合エネルギーの見積もりには基底関数車ね合わせ誤差は必須となります。」（「すぐできる量子科学計算」より） FireflyのマニュアルにはBSSEを補正する方法が三つ説明されています。その中でゴーストアトム（電荷は無く、基底関数だけを持つ）を使って結合エネルギーを求める方法を使ってみました。水二量体（A+B)の場合で説明します。 ① まず水分子一つを最適化します。そのエネルギーをE(H2O)とする ② 次に水分子をもう一つ適当な位置に追加して最適化します。（補正前の二量体エネルギー） ③ ②で最適化した構造のままAをゴーストアトムに指定してENERGY計算します。-E(B)ab- Bをゴーストアトムにして同じく計算します。-E(A)ab- ④ 次に②で最適化した構造のままA、B単独でそれぞれのエネルギー計算をします。-E(A)a- -E(B)b- ⑤ ③と④で得た四つのエネルギー値から下記の計算で誤差の大きさが分かります。 誤差＝E(A)ab - E(A)a + E(B)ab - E(B)b ⑥ 補正前の結合エネルギーは②でえた 補正前の二量体エネルギーから①で得たE(H2O) の二倍を引くことで得られます。 ⑦ 補正後の結合エネルギーは補正前の結合エネルギーから誤差分を差し引くことで求めることができます。 ゴーストアトムの指定はインプットファイルで＄DATAにある原子の電荷の前に負符号を付ければできます。下記は最初の水分子の三原子をゴーストアトムにしています。 O -8.0 0.000000 0.000000 0.000000 H -1.0 0.943392 0.000000 0.000000 H -1.0 -0.284431 0.900058 0.000000 O 8.0 -2.369099 -1.844959 0.000000 H 1.0 -1.561301 -1.349841 0.000000 H 1.0 -2.134704 -2.757777 0.000000 ・遷移金属を含む系の計算・ ① Fireflyでは初期軌道としてＨＵＣＫＥＬ（Ｒｎまで）とＨＣＯＲＥ（すべての元素）を選択できる。＄GUESSグループのGUESSというオプションで指定する。通常はＨＵＣＫＥＬがベターだが、遷移金属を含む系ではＨＣＯＲＥがより良い結果を与える。これはＨＵＣＫＥＬは常に遷移金属のECPを適切に処理するとは限らない為。 ② 遷移金属に限らず大きい基底関数を使用する場合、ある関数が別の関数の線形結合で書くことができる場合がある。これはdiffuse関数が使用される時によく起こる。この状況は数値不安定をもたらしＳＣＦの収束が低下する。このような線形依存性が検出されると「A PARTIAL LINEAR DEPENDENCE IN YOUR ATOMIC BASIS」というエラーメッセージが出力される。 この場合Fireflyのマニュアルによると良い出発点としてより正確な積分を行うことが推奨されている。$CONTRLグループにINTTYP=HONDO, ICUT=11, ITOL=30 と入力する。 また線形依存性が1.0Ｄ-7以下の固有値をもたらす場合収束に達しない場合があるので、そのような関数を手動で削除することによって線形依存性を減らすことができる。出力される関数のリストの中で最小の固有値を持つ関数を確認し、その関数を削除した基底系で計算を行うことで問題が解決することがある。時には問題となる関数をすべて削除するまで同じ手順を繰り返す必要がある。例えば関数のリストの最後が「0.8695E-07 - C 1 S, SHELL 4, AO 10」であれば、最小の固有値は一番目の炭素原子のSHELL 4に由来することが分かる。PUNCHファイルの$ DATAブロックからこの関数を削除して、残りの$ DATAブロック全体をインプットファイルに貼り付けて計算を行う。このためには・外部の基底系を使用する・ で説明した外部ファイルを使用する必要がある。 しかしこの方法は基底関数を削除することでエネルギー、勾配、その他の特性を変えることになるので注意が必要。 ・溶媒効果・ 溶媒を連続誘電体として近似する方法として下記の二つがある。 ①ＳＣＲＦ（自己無撞着反応場） 溶質を溶媒中の空孔の中に置き、溶質の電荷とそれにより発生する溶媒の静電ポテンシャルの相互作用により安定化される状態を求める。＄ＳＣＲＦグループでＤＩＥＬＥＣ＝＊＊（溶媒の誘電率） ＲＡＤＩＵＳ＝＊＊（空孔の球の半径 単位はÅ）空孔球の大きさは分子の端の原子のファンデルワールス半径の1.2倍の範囲が収まるように設定する。 ②ＰＣＭ 溶質分子の各原子を中心とする球（半径はその原子のファンデルワールス半径の１．２倍）の結合によって作られる範囲を空孔とする。＄ＰＣＭグループでＳＯＬＶＮＴ＝＊＊（溶媒の種類）と指定する。溶媒の種類としては下記が格納されている。 WATER (or H2O) CH3OH C2H5OH CLFORM (or CHCl3) METHYCL (or CH2Cl2) 12DCLET (or C2H4Cl2) CTCL (or CCl4) BENZENE (or C6H6) TOLUENE (or C6H5CH3) CLBENZ (or C6H5Cl) NITMET (or CH3NO2) NEPTANE (or C7H16) CYCHEX (or C6H12) ANILINE (or C6H5NH2) ACETONE (or CH3COCH3) THF DMSO (or DMETSOX) ・・・ＩＮＰＵＴ ＰＤＦより・・・ ＳＯＬＶＮＴ＝ＩＮＰＵＴと入力して溶媒の詳細をを指定する高度なキーワードも備えている。 ・ＭＥＰ計算・ 特定の二原子の位置関係（距離）を段階的に変化させながら、各々の段階でそれ以外の原子団を構造最適化してエネルギー最小値を求める。(下記の方法は角度の変化には無効です） $CONTRL RUNTYP=RSURFACE NZVAR=＊＊（Ｚマトリックスの座標数） $STATPT IFREEZ(1)=＊＊（Ｚマトリックスの何番目を固定するかを指定） $SURF NDISP1=＊＊（何回計算するかを指定）DISP1=＊＊（座標値の変化量を指定） VECT1(このカッコ内にIFREEZで指定した値を入力）＝1 ORIG1=＊＊（座標 の初期値） ・基底関数の種類・ 出力結果には下図のように基底関数の種類が示される。（太枠内）  Sはｓ軌道、X、Y、Zはｐ軌道、XX、YY、ZZ、XY、YZ、ZXはｄ軌道、XXX、YYY、ZZZ、XXY、XYY、YYZ、YZZ、ZZX、ZXX、XYZはｆ軌道 ・対称操作の利用・ 対称グループは$
DATAグループで指定できる。例えば水分子をＣ2Ｖと指定する時は以下の通り入力する。 標準では主回転軸はＺ軸、主回転軸に垂直な軸はＸ軸、σｖ平面はＸＺ面、σｈ平面はＸＹ面。 ・SCF計算が収束しない時・ |
