CONFLEX Tutorials

量子化学計算プログラムGaussianを用いた構造最適化と配座探索の実行

CONFLEXで行う構造最適化および配座探索は、通常CONFLEX本体に搭載されている古典的な分子力場を用いて実行しています。エネルギー極小構造を高速かつ効率的に網羅する上で分子力場は有用ですが、対象とする系に必要な原子タイプやパラメーターが無い場合は、これらを補わなければ極小構造を求めることができません。また、力場では記述するのが難しい電子状態での極小構造を求めたい場合、力場を使って算出した配座について量子化学計算を行った場合と、最初から量子化学計算を行って配座探索を実行した場合とで、得られる配座が異なることも考えられます。

CONFLEX8以降では、分子力場を用いる代わりに量子化学計算プログラムGaussianを外部プログラムとして呼び出し、最適化計算や配座探索を実行することができます。ここでは、Gaussianを呼び出して配座探索計算を行う場合の設定方法と、その結果の例を示します。

【環境設定】

外部プログラムとしてGaussian(16または09)を呼び出すためには、CONFLEXがインストールされているマシンにGaussianもインストールされている必要があります。またGaussianを実行する際には、お使いのアカウントでGaussian実行用の環境設定があらかじめ為されている必要があります。

Mac, Linuxの場合

環境変数”g16root”または”g09root”をあらかじめ設定してください。尚、環境変数「GAUSS_SCRDIR」を指定していない場合は、中間ファイルがカレント・ディレクトリに出力されます。

Windowsの場合

環境変数「GAUSS_EXEDIR」に、Gaussianがインストールされているフォルダー名(C:¥G16W、等)を設定してください。尚、キーワード「EXT_GAU_SCRDIR=」に中間ファイルを出力するフォルダーを指定しない場合は、「%GAUSS_EXEDIR%¥Scratch」が出力ファイルおよび中間ファイルを出力するフォルダーとして設定されます。

また、環境設定を含めた設定ファイルを別途用意することで、Gaussianの環境設定を計算実行時のみ有効にする方法もあります。

【Propylene glycolのB3LYP/6-31G(d)レベルによる配座探索計算】

ここではPropylene glycolを利用して、B3LYP/6-31G(d)レベルで配座探索を行う方法を説明します。 Propylene glycolの入力構造

Propylene glycolの座標データ(PropyleneGlycol.mol)

PropyleneGlycol.mol


 13 12  0  0  0  0  0  0  0  0  0    0
   -4.0461   -0.3459    0.0000 C   0  0  0  0  0  0  0  0  0  0  0  0
   -3.6894    0.1585    0.8737 H   0  0  0  0  0  0  0  0  0  0  0  0
   -3.6894    0.1585   -0.8737 H   0  0  0  0  0  0  0  0  0  0  0  0
   -5.1161   -0.3459    0.0000 H   0  0  0  0  0  0  0  0  0  0  0  0
   -3.5328   -1.7978    0.0000 C   0  0  0  0  0  0  0  0  0  0  0  0
   -3.8895   -2.3022    0.8737 H   0  0  0  0  0  0  0  0  0  0  0  0
   -1.9928   -1.7979   -0.0000 C   0  0  0  0  0  0  0  0  0  0  0  0
   -1.6361   -1.2935   -0.8737 H   0  0  0  0  0  0  0  0  0  0  0  0
   -1.6361   -1.2934    0.8736 H   0  0  0  0  0  0  0  0  0  0  0  0
   -4.0095   -2.4719   -1.1676 O   0  0  0  0  0  0  0  0  0  0  0  0
   -4.5428   -3.2264   -0.9068 H   0  0  0  0  0  0  0  0  0  0  0  0
   -1.5162   -3.1461    0.0000 O   0  0  0  0  0  0  0  0  0  0  0  0
   -0.5757   -3.1535    0.1927 H   0  0  0  0  0  0  0  0  0  0  0  0
  1  2  1  0  0  0  0
  1  3  1  0  0  0  0
  1  4  1  0  0  0  0
  1  5  1  0  0  0  0
  5  6  1  0  0  0  0
  5  7  1  0  0  0  0
  5 10  1  0  0  0  0
  7  8  1  0  0  0  0
  7  9  1  0  0  0  0
  7 12  1  0  0  0  0
 10 11  1  0  0  0  0
 12 13  1  0  0  0  0
M  END

[Interfaceから実行する場合]

PropyleneGlycol.molファイルをCONFLEX Interfaceを用いて開きます。 Interface Propylene Glycol

Calculationメニューから「CONFLEX」を選択し、開いた計算設定ダイアログのDetail Settingsをクリックします。詳細設定ダイアログが開きます。 Basic Settings

まず、詳細設定ダイアログの「General Settings」ダイアログにある「Calculation Type:」のプルダウンメニューから「Conformation Search」を選択します。さらに、「Estimating Method:」のプルダウンメニューから「External Program」を選択します。これにより、構造最適化などの計算が外部プログラム(Gaussian)を利用して行われるようになります。 General Settings

次に、「External Program」ダイアログで、Gaussianの計算設定を行います。
「Job Title:」はコメントであり、ここでは「Propylene Glycol B3LYP/6-31G(d) optimization」とします。「%Mem=」は「128」とし、単位のプルダウンメニューから「MW」を選択します。「%NProcShared=」は「2」とします。「Route Section:」は「B3LYP/6-31G(d) Opt」とし、B3LYP/6-31G(d)での構造最適化を行うように設定します。 Gaussian Settings

Gaussianの計算設定は、ファイルを用いて行うこともできます。【ファイルを用いたGaussianの計算設定】をご参照ください。

次に、「Conformation Search」ダイアログで、「Search Limit:」の値を「5.0」にし、配座探索空間の範囲を決めます。 Conformation Search 設定が終わりましたら、詳細設定ダイアログのEdit & Submitをクリックします。

開いたダイアログに「EDIF_HARD=1.0D-3」を追加します。これは、新たに得られた配座を保存するかどうかを判定する際のエネルギー差の閾値を0.001 kcal/molとする指定です。 Edit and Submit 設定が終わりましたらSubmitをクリックします。B3LYP/6-31G(d)レベルで配座探索が行われます。

[コマンドラインから実行する場合]

計算設定は、PropyleneGlycol.iniファイルにキーワードを記述することで行います。

PropyleneGlycol.iniファイル

CONFLEX 
SEL=5 
EXTERNAL_PROGRAM=GAUSSIAN 
EXT_GAU_LINK0=(mem=128MW) 
EXT_GAU_LINK0=(nprocshared=2) 
EXT_GAU_ROUTE=("B3LYP/6-31G(d) opt”) 
EXT_GAU_COMMENT=("Propylene Glycol B3LYP/6-31G(d) optimization”) 
EXT_GAU_CHARGE=0 
EXT_GAU_SPIN=1 
EXT_PREOPT=OFF 
EDIF_HARD=1.0D-3

「CONFLEX」キーワードは、配座探索計算を行うこと意味します。
「SEL=5」キーワードは、配座探索のエネルギー上限値を5.0 kcal/molと設定します。また、「EXTERNAL_PROGRAM=」は使用する外部プログラム名、「EXT_GAU_LINK0=」はLink0コマンド(複数指定可)、「EXT_GAU_ROUTE=」はルートセクション、「EXT_GAU_COMMENT=」はコメント、「EXT_GAU_CHARGE=」は系全体の電荷、「EXT_GAU_SPIN=」は系のスピン多重度をそれぞれ表しています。尚、「EXT_GAU_ROUTE=」と「EXT_GAU_COMMENT=」では他と異なり「(“」と「”)」で指定しますので、ご注意ください。「EXT_PREOPT=OFF」はGaussianを呼び出す前に分子力場による最適化を行わないようにする指定、「EDIF_HARD=1.0D-3」は新たに得られた配座を保存するかどうかを判定する際のエネルギー差の閾値を0.001 kcal/molとする指定です。

Gaussianの計算設定は、ファイルを用いて行うこともできます。【ファイルを用いたGaussianの計算設定】をご参照ください。

PropyleneGlycol.molとPropyleneGlycol.iniをフォルダに格納し、下記コマンドを実行してください。計算が始まります。

C:\CONFLEX\bin\flex9a_win_x64.exe   -par   C:\CONFLEX\par   PropyleneGlycolenter

上記は、Windowsの場合です。他のOSにおける実行コマンドについては、本書の「実行方法」を参照してください。

計算結果

PropyleneGlycol.ls1の出力内容とGaussianの入出力ファイルを示します。Gaussianの入出力ファイル名はls1ファイルの「CONF. ID」に対応しており、最安定構造はPropyleneGlycol_00006.logおよびPropyleyeGlycol_00006.fchkに出力されています。またls1ファイルの「ORIGINAL STERIC E」は、Gaussianの計算で得られた全エネルギーをkcal/molに変換したもの(最安定構造の全エネルギーは-269.56133526 a.u.なので、-269.56133526×627.5095=-169152.2987)です。尚この計算では、Gaussian 16 Rev. A.03を使用しています。

PropyleneGlycol.ls1:

 ====================================================================================
            CONF.           ORIGINAL               DISTRI-                    NO.
    NO.     ID             STERIC E   DELTA E     BUTION      INIT.  REOPT.  NEG.
------------------------------------------------------------------------------------
      1   00000006      -169152.2987    0.0000     21.2884       T*     F+     0
      2   00000017      -169152.2709    0.0278     20.3129       T*     F+     0
      3   00000019      -169152.1497    0.1490     16.5548       T*     F+     0
      4   00000011      -169151.9408    0.3579     11.6356       T*     F+     0
      5   00000001      -169151.9312    0.3675     11.4487       T*     F+     0
      6   00000008      -169151.6798    0.6189      7.4899       T*     F+     0
      7   00000012      -169151.4409    0.8578      5.0044       T*     F+     0
      8   00000018      -169151.1376    1.1611      2.9993       T*     F+     0
      9   00000007      -169150.9622    1.3365      2.2310       T*     F+     0
     10   00000016      -169149.5319    2.7668      0.1995       T*     F+     0
     11   00000020      -169149.3686    2.9301      0.1515       T*     F+     0
     12   00000009      -169149.2519    3.0468      0.1244       T*     F+     0
     13   00000022      -169149.2130    3.0857      0.1165       T*     F+     0
     14   00000021      -169149.0687    3.2300      0.0913       T*     F+     0
     15   00000023      -169149.0383    3.2604      0.0868       T*     F+     0
     16   00000024      -169149.0232    3.2755      0.0846       T*     F+     0
     17   00000013      -169148.9887    3.3101      0.0798       T*     F+     0
     18   00000003      -169148.8937    3.4050      0.0680       T*     F+     0
     19   00000010      -169147.6211    4.6776      0.0079       T*     F+     0
     20   00000014      -169147.4593    4.8394      0.0060       T*     F+     0
     21   00000004      -169147.3564    4.9423      0.0051       T*     F+     0
     22   00000015      -169147.3371    4.9616      0.0049       T*     F+     0
     23   00000005      -169147.3181    4.9806      0.0048       T*     F+     0
     24   00000002      -169147.2083    5.0904      0.0040       T*     F+     0
------------------------------------------------------------------------------------
MINIMUM ENERGY: -169152.2987  KCAL/MOL
AVERAGE ENERGY: -169151.9983  KCAL/MOL (BASED ON THE BOLTZMANN LAW)
====================================================================================

Gaussian入出力ファイルの一覧:

PropyleneGlycol_00001.fchk  PropyleneGlycol_00008.fchk  PropyleneGlycol_00015.fchk
PropyleneGlycol_00001.gjf   PropyleneGlycol_00008.gjf   PropyleneGlycol_00015.gjf
PropyleneGlycol_00001.log   PropyleneGlycol_00008.log   PropyleneGlycol_00015.log
PropyleneGlycol_00002.fchk  PropyleneGlycol_00009.fchk  PropyleneGlycol_00016.fchk
PropyleneGlycol_00002.gjf   PropyleneGlycol_00009.gjf   PropyleneGlycol_00016.gjf
PropyleneGlycol_00002.log   PropyleneGlycol_00009.log   PropyleneGlycol_00016.log
PropyleneGlycol_00003.fchk  PropyleneGlycol_00010.fchk  PropyleneGlycol_00017.fchk
 ・・・・・

Propylene glycolをMMFF94sで配座探索を行う(Search Limit: 7(SEL=7))と、21個の配座異性体が得られます。これらをB3LYP/6-31G(d)レベルで再度最適化すると、MMFF94sで5番目に安定な配座が最も安定になります(下図)。 External Fig2 Propylene glycolのB3LYP/6-31G(d)レベルでの最安定構造

また、MMFF94sで10番目に安定な構造をB3LYP/6-31G(d)で最適化すると他の構造と同一になる(下図)ため、最終的にB3LYP/6-31G(d)で得られた配座は20個となりました。 External Fig3 MMFF94sで10番目に安定な構造(左)をB3LYP/6-31G(d)で再度最適化すると、異なる配座異性体(右)として収束する

先にMMFF94s力場を用いて探索し、得られた配座異性体をGaussianで再度最適化した結果と、Gaussianのみを利用して配座探索を行った結果を比較しますと、次の図に示す4つの構造が新たに見つかったことがわかりました。 External Fig4 B3LYP/6-31G(d)による配座探索で新たに得られた構造

【ファイルを用いたGaussianの計算設定】

Link0コマンドやルートセクションなどを個別のキーワードで設定する以外に、ファイルで設定する方法があります。ここでは、ファイルを用いたGaussianの計算設定について説明します。

*例1

Gaussian実行用の入力ファイルの内容を下記とします。

%mem=32MW
%nprocshared=2
#HF/3-21G opt

test of HF/3-21G calculation
0   1
molecular specification
...

分子構造データより上の設定を、gaussian.hdrファイルとして保存します。

gaussian.hdrファイル:

%mem=32MW
%nprocshared=2
#HF/3-21G opt

test of HF/3-21G calculation

0   1

[Interfaceから実行する場合]

詳細設定ダイアログの「External Program」ダイアログで、「Header File:」のチェックボックスにチェックを入れます。 External Header File

次に、「Select...」をクリックして、「gaussian.hdr」ファイルを選択してください。gaussian.hdrの内容を元に、Gaussianの計算設定がなされます。 External Header File set

[コマンドラインから実行する場合]

PropyleneGlycol.iniファイルは、下記のように変更します。

PropyleneGlycol.iniファイル

CONFLEX
SEL=5
EXTERNAL_PROGRAM=GAUSSIAN
EXT_GAU_HEADER_FILE=(gaussian.hdr)
EXT_PREOPT=OFF
EDIF_HARD=1.0D-3

PropyleneGlycol.molとPropyleneGlycol.ini、gaussian.hdrをフォルダに格納し、コマンドを実行します。

*例2

キーワード「Gen」を指定して基底関数を設定する時のように分子構造データの後に設定を記述する場合も、別ファイルに記述して読み込ませることが可能です。Gaussian実行用の入力ファイルの内容を下記とします。

%mem=20MW
%nprocshared=2
# B3LYP/gen Opt

TEST OF B3LYP/6-31G(d) CALCULATION

0   1
molecular specification
...

H C O 0
6-31G(d)
****

分子構造データより上の設定を、gaussian.hdrファイルとして保存します。

gaussian.hdrファイル:

%mem=32MW
%nprocshared=2
#HF/3-21G opt

test of HF/3-21G calculation

0   1

分子構造データより下の設定を、gaussian.ftrファイルとして保存します。

gaussian.ftrファイル:

H C O 0
6-31G(d)
****

[Interfaceから実行する場合]

詳細設定ダイアログの「External Program」ダイアログで、「Header File:」のチェックボックスにチェックを入れます。 External Program Dialog 次に、「Select...」をクリックして、「gaussian.hdr」ファイルを選択してください。 Select Gaussian Header

また、「Footer File:」のチェックボックスにチェックを入れます。 External Footer 次に、「Select...」をクリックして、「gaussian.ftr」ファイルを選択してください。gaussian.hdrとgaussian.ftrの内容を元に、Gaussianの計算設定がなされます。 Select Gaussian Footer

[コマンドラインから実行する場合]

PropyleneGlycol.iniファイルは、下記のように変更します。

PropyleneGlycol.iniファイル

CONFLEX
SEL=5
EXTERNAL_PROGRAM=GAUSSIAN
EXT_GAU_HEADER_FILE=(gaussian.hdr)
EXT_GAU_FOOTER_FILE=(gaussian.ftr)
EXT_PREOPT=OFF
EDIF_HARD=1.0D-3

PropyleneGlycol.molとPropyleneGlycol.ini、gaussian.hdr、gaussian.ftrをフォルダに格納し、コマンドを実行します。

【ジョブスケジューラを使ったGaussianの実行(Linuxのみ)】

CONFLEXから実行されるGaussianジョブを、ジョブスケジューラー経由で実行する方法を説明します。ここでは、ジョブスケジューラーはUniva Grid Engineとします。

[コマンドラインから実行する場合]

PropyleneGlycol.iniファイルは、下記のように変更します。

PropyleneGlycol.iniファイル

CONFLEX
SEL=5
EXTERNAL_PROGRAM=GAUSSIAN
EXT_JOB_COM=("qsub -sync yes”)
EXT_JOB_FILE=qsub.sh
EXT_GAU_LINK0=(mem=128MW)
EXT_GAU_LINK0=(nprocshared=2)
EXT_GAU_ROUTE=("B3LYP/6-31G(d) opt”)
EXT_GAU_COMMENT=("Propylene Glycol B3LYP/6-31G(d) optimization”) 
EXT_GAU_CHARGE=0
EXT_GAU_SPIN=1
EXT_PREOPT=OFF
EDIF_HARD=1.0D-3

「EXT_JOB_COM=」は、ジョブスケジューラーの実行コマンド(qsub -sync yes)を設定しています。ジョブスケジューラーを用いる場合は必ずインタラクティブに実行するオプション(Univa Grid Engineでは「-sync yes」)を加えてください。また、「EXT_JOB_FILE=」は、ジョブスケジューラーに対する設定などが記述されたファイルを指定します。

qsub.shファイル:

#!/bin/bash
#$ -S /bin/bash
#$ -N g16_from_conflex
#$ -cwd
#$ -V
#$ -j y
#$ -q all.q
#$ -pe sme 2

export g16root=/usr/local

これにより、Gaussianのジョブは「qsub -sync yes qsub.sh」で実行され、ジョブスケジューラー経由でのジョブ実行となります。なお、qsub.shにはGaussianの実行コマンドがありませんが、CONFLEXが自動的に追加します。qsub.shは、PropyleneGlycol.molとPropyleneGlycol.iniと同じフォルダに格納してください。

【一時的な環境設定の有効化】

Gaussianの環境設定を計算実行時のみ有効にする方法を説明します。

Windowsの場合

[Interfaceからの実行]

詳細設定ダイアログのEdit & Submitをクリックし、開いたダイアログに「EXT_JOB_FILE=setup.txt」を追加します。 Edit Submit

[コマンドラインから実行する場合]

PropyleneGlycol.iniファイルは、下記のように変更します。

PropyleneGlycol.iniファイル

CONFLEX
SEL=5
EXTERNAL_PROGRAM=GAUSSIAN
EXT_JOB_FILE=setup.txt
EXT_GAU_LINK0=(mem=128MW) 
EXT_GAU_LINK0=(nprocshared=2) 
EXT_GAU_ROUTE=("B3LYP/6-31G(d) opt”) 
EXT_GAU_COMMENT=("Propylene Glycol B3LYP/6-31G(d) optimization”) 
EXT_GAU_CHARGE=0 
EXT_GAU_SPIN=1 
EXT_PREOPT=OFF
EDIF_HARD=1.0D-3

ここで「EXT_JOB_FILE=」は設定ファイル名を表します。設定ファイルsetup.txtの内容を以下のようにすると、CONFLEXからGaussianを呼び出す毎に、このファイルの内容とGaussian実行コマンドを含むbatファイルを作成し、計算を実行します。

setup.txtファイル:

set GAUSS_EXEDIR=C:\G16W

また、Edit & Submitをクリックし開いたダイアログやiniファイルに「EXT_GAU_SCRDIR=D:\Scratch」を加えると、中間ファイルの保存先をD:¥Scratch以下に変更できます。加えない場合は、「%GAUSS_EXEDIR%¥Scratch」ファルダ−に保存されます。

Mac, Linuxの場合

[Interfaceから実行する場合]

詳細設定ダイアログのEdit & Submitをクリックし、開いたダイアログに「EXT_JOB_COM=(“sh")」と「EXT_JOB_FILE=setup.sh」を追加します。 Edit and Submit

[コマンドラインから実行する場合]

PropyleneGlycol.iniファイルは、下記のように変更します。

PropyleneGlycol.iniファイル

CONFLEX
SEL=5
EXTERNAL_PROGRAM=GAUSSIAN
EXT_JOB_COM=("sh”)
EXT_JOB_FILE=setup.sh
EXT_GAU_LINK0=(mem=128MW) 
EXT_GAU_LINK0=(nprocshared=2) 
EXT_GAU_ROUTE=("B3LYP/6-31G(d) opt”) 
EXT_GAU_COMMENT=("Propylene Glycol B3LYP/6-31G(d) optimization”) 
EXT_GAU_CHARGE=0 
EXT_GAU_SPIN=1
EXT_PREOPT=OFF
EDIF_HARD=1.0D-3

ここで「EXT_JOB_COM=」は実行コマンド、「EXT_JOB_FILE=」は設定ファイル名をそれぞれ表します。設定ファイルsetup.shの内容を以下のようにすると、CONFLEXからGaussianを呼び出す毎に、このファイルの内容とGaussian実行コマンドを含むスクリプトファイルを作成し、shコマンドで計算を実行します。

setup.shファイル:

export g16root=/usr/local
export GAUSS_SCRDIR=/tmp