この記事では化合物の3次元構造データの一つであるxyzファイルから、Gaussianのインプットファイルを作成する方法を紹介します。
・xyzファイルを量子化学計算のインプットファイルに使用するための方法
・Gaussianのインプットファイルのルール
・実際のGaussianインプットファイル例
今回使用するxyzファイル
今回はベンゼンの3次元構造データに対する量子化学を行ってみましょう。
xyzファイルは以下のものです。
12 benzene C -0.00590 0.00560 -0.06240 C 0.03290 -0.66660 -1.27820 C 0.06640 -2.05580 -1.30340 C 0.06120 -2.77280 -0.11280 C 0.02240 -2.10050 1.10300 C -0.01110 -0.71140 1.12820 H -0.03200 1.08720 -0.04280 H 0.03700 -0.10840 -2.20520 H 0.09660 -2.57920 -2.25000 H 0.08730 -3.85430 -0.13240 H 0.01830 -2.65870 2.03000 H -0.04130 -0.18800 2.07480
xyzファイルの作り方については以前の記事をどうぞ。
Gaussianのインプットファイルの形式(benzene_single_point_e.com)
ベンゼンのxyzファイルを元に、計算レベルωB97X-D/Def2-SVP(P)でのDFT計算による一点計算(構造最適化はせず、xyzファイルの3次元構造における電子エネルギーの計算)を行うためのインプットファイル(benzene_single_point_e.com)を以下に示します。
各項目の説明は後ほど記載します。
%nprocshared=4 %chk=benzene_single_point_e.chk # wb97xd/Def2SVPP benzene 0 1 C -0.00590 0.00560 -0.06240 C 0.03290 -0.66660 -1.27820 C 0.06640 -2.05580 -1.30340 C 0.06120 -2.77280 -0.11280 C 0.02240 -2.10050 1.10300 C -0.01110 -0.71140 1.12820 H -0.03200 1.08720 -0.04280 H 0.03700 -0.10840 -2.20520 H 0.09660 -2.57920 -2.25000 H 0.08730 -3.85430 -0.13240 H 0.01830 -2.65870 2.03000 H -0.04130 -0.18800 2.07480 !分子指定セクションのあとは空行が1行必要です。付加情報セクションがない場合は2行空行にします。
gaussianのインプットファイルは
・Link0コマンド
・ルートセクション
・タイトルセクション
・分子指定セクション(ここでxyzファイルを使います)
・付加情報セクション
に分かれています。
またコメントアウト(プログラムに読まれないメモ書きのようなもの)は、!記号を先頭に書いて記載します。
上のインプットファイルでは22行目がコメントアウトなので、この行の記載がなくても計算は実行できます。
このブログでは、時々コメントアウトを入れて読者の方に理解しやすいようにしようと思いますが、コメントアウトのルールは使用するプログラムごとに異なるので注意してください。
それでは、xyzファイルから上記インプットファイルをどのようにして作成するのか、それぞれのセクションがどんな設定を行っているのかを見ていきましょう。
各セクションの記載方法
Link0コマンド
Link0コマンドはgaussianのインプットファイルの先頭に記載するコマンドで、計算ジョブの実行に関するシステムリソースやファイル設定を指定する役割があります。
「%」記号で始まる部分で、複数のコマンドを記載する際はその都度改行して「%」から始めます。
先ほど示したインプットファイル(benzene_single_point_e.com)では以下の部分がLink0コマンドです。
%nprocshared=4 %chk=benzene_single_point_e.chk
HPCシステムズのWebサイトにLink0コマンドの一覧が記載されているので、この中から基本的なもの2つを以下に示します。
使用するプロセッサ数の指定(%nprocshared)
%nprocshared=4
は4プロセッサを用いた並列計算を実施すること意味しています。
単純に数字を大きくしたほうが高速で計算が実行されますが、最大数はCPUのコア数に依存します。
スーパーコンピューターでは20-40程度にすることが多い部分ですが、プロセッサ数が多くなるほど比例的に計算が速くなるわけではなく、プロセッサ数上昇の効果はだんだん少なくなっていくことに注意が必要です。
計算結果が格納されるチェックポイントファイルの指定(%chk)
チェックポイントファイル(.chk)の名前を指定します。計算の途中経過や結果を保存するために使用されます。
今回の例では計算が実行されると、そのチェックポイントファイルの名前は「benzene_single_point_e.chk」とすることを指示しています。
GaussViewでインプットファイルを作成した場合はフルパスで記載されることがあり、環境によっては計算実行のエラーの原因となるため、慣れないうちは「/」記号が入っていないことを確認しておきましょう。
このチェックポイントファイルは、この計算結果を読み込んで次の計算ジョブを実行する際に読み込む形で利用することができます。
次の段階の計算で前のチェックポイントファイルを読み込むためには、%chkでファイル名を指定して読み込むことができます。
チェックポイントファイルの上書きを避けたいのであれば、以下のように%oldchkで以下のように記載することもできます。(管理人はこちらの方法を使っています。)
%chk=next_job.chk !次の計算のチェックポイントファイル名を指定 %oldchk=benzene_single_point_e.chk !前の計算結果のチェックポイントファイル名を指定
ルートセクション
ルートセクションは、計算手法、基底関数、計算タイプ、オプションなど、どのような計算を実行するのかを指示する部分です。
行の先頭に「#」記号をつけ、終わりには空行を1行入れるのがルールです。
benzene_single_point_e.comでは以下の部分がルートセクションです。
# wb97xd/Def2SVPP
これは、汎関数としてωB97X-Dを、基底関数としてDef2-SVP(P)を用いる計算レベルで計算を実行することを意味しています。
その他に記載がないため、上の計算レベルで入力構造の電子エネルギーを計算することとなります。
同じ計算レベルで構造最適化を行う場合は
# wb97xd/Def2SVPP Opt
のように、新しく「Opt」というキーワードを半角スペースの後に記載します。
計算レベルをB3LYP/6-31G(d)としたいのであれば
# B3LYP/6-31G(d) Opt
と記載すればよいです。
Link0コマンドのように、その都度「#」をつける必要はないことに注意しましょう。
キーワードは大文字と小文字の区別がないため、例えばB3LYPとb3lypは同じものと認識されます。
重ねてですが、すべて記載したら空行を入れましょう。
キーワード内な追加キーワード
ルートセクションの各キーワード(Opt等)では、それぞれのキーワード内で計算の詳細な条件をさらに変更するためのキーワード(一つ下の階層のキーワード)を入れることができます。
例えば構造最適化のOptキーワードでは
・遷移状態の構造最適化を行うための「TS」
・構造最適化における最大サイズ(1ステップで動く大きさのようなもの)を指定する「MaxStep=N(Nは任意の数字)」
・座標の固定化等を行うための「ModRedundant」
などがあります。
これらを適用するには元のキーワードに「=」をつけてその後に追加オプションを追加します。
遷移状態最適化では
Opt=TS
となります。
2つ以上の追加キーワードを用いるときはコンマ「,」で分けて小括弧「()」で括ります。
例えば以下のような感じです。
Opt=(TS,MaxStep=30)
タイトルセクション
計算ジョブの簡単な説明やメモを記述する部分です。分子の名前や計算レベルを書いておくとよいですが、その内容は計算の実行に直接影響しません。
benzene_single_point_e.comでは以下の部分が該当します。
benzene
このセクションも、終わりの目印として空行が1行必要です。
分子指定セクション
計算する分子の3次元構造や電荷、多重度を定義する部分です。
ここで記載した構造に対して、ルートセクションで指定した計算を実行するということになります。
benzene_single_point_e.comでは以下の部分です。
0 1 C -0.00590 0.00560 -0.06240 C 0.03290 -0.66660 -1.27820 C 0.06640 -2.05580 -1.30340 C 0.06120 -2.77280 -0.11280 C 0.02240 -2.10050 1.10300 C -0.01110 -0.71140 1.12820 H -0.03200 1.08720 -0.04280 H 0.03700 -0.10840 -2.20520 H 0.09660 -2.57920 -2.25000 H 0.08730 -3.85430 -0.13240 H 0.01830 -2.65870 2.03000 H -0.04130 -0.18800 2.07480
一行目の「0 1」はxyzファイルにない部分ですね。
この数字は左から順に
・系の総電荷(今回は±0)
・系の多重度(一重項のときは1、三重項のときは3)
を意味しており、不適切な値とするとエラーの原因になります。
なお、三重項状態を計算するときはルートセクションの計算手法を非制限法にする必要があるので注意してください。
(そのうち記事にするかもしれませんが、B3LYPを汎関数とするのであればUB3LYPのキーワードを用います。)
スピンが分子全体で揃っていれば一重項です。
様々な化学種の多重度と電荷については、PC CHEM BASICS.COMさんの記事: スピン多重度と電荷を設定しようを参考にしてみてください。
その下から続く以下の部分はxyzファイルの3行目以降をそのままコピーして貼り付ければよいです。
C -0.00590 0.00560 -0.06240 C 0.03290 -0.66660 -1.27820 C 0.06640 -2.05580 -1.30340 ・・・
3次元構造を記載したら最後に空行を入れますが、このとき注意が必要です。
①付加情報セクションの記載が不要なインプットファイルの場合
→空行を2行入れる
②付加情報セクションの記載が必要なインプットファイルの場合
→空行を1行入れる
付加情報セクション
付加情報セクションは、分子指定セクションの後に追加で記述されるセクションです。
特定の制約や追加のオプションを設定するために使用されます。
よく使うのは、座標の一部を固定した構造最適化や、relaxed scanを行う際の結合長、角度角、二面角の指定です。
このような最適化では、ルートセクションにおいてOpt=(modredundant)としたうえで、以下のような記載をします。
5 4 3 2 F !この例では5番目、4番目、3番目、2番目からなる二面角を固定化(F: Freeze)することを意味します。
詳細はrelaxed scanのインプットファイルを紹介する際に説明しようと思います。
計算の実行
以上でインプットファイル(benzene_single_point_e)が作成できました。
計算の実行はこのインプットファイルを用いて
g16 < benzene_single_point_e.com > benzene_single_point_e.log
のようなコマンドで行えます。
benzene_single_point_e.logは計算のログファイル名の指定なので、好みに応じて変更してかまいません。
スーパーコンピューターでは実行方法が変わるため、管理者に相談するようにしてください。
終わりに
今回はxyzファイルからGaussianのインプットファイルを作成する方法を説明しました。
計算オプションは多種多様で、今回だけではすべてを紹介できませんが気になる方はGaussianの公式サイトを御覧ください。
有機化学でよく使用するインプットファイルについてはChem-Stationさんの記事: Gaussian Input File データベースでもいくつか紹介されています。
このブログでも、役に立つインプットファイルについては別の記事で紹介するつもりです。
コメント・質問等はX (@chemmodelcomp)にてご連絡ください。