gnuplot

【インストール】
まず、gnuplotのホームページ(http://www.gnuplot.info/)からダウンロード。
windows版は、gp442win32.zip (2010.01現在)。 意味は、ver.4.4.2で、windows 32bit版。
zipファイルを解凍。インストール作業は不要。
解凍されたフォルダ内のwgnuplot.exeが実行ファイル。(注) gnuplot.exeではない。

【初期設定】
wgnuplot.exeを実行すると、gnuplotのウィンドウが開く。が、文字が読めない。
ウィンドウ上で右クリックして、Choose Fontタブをクリック。読みやすいフォントを設定。
再び、右クリックして、Update C:\... をクリックすると、設定したフォントが保存される。

【使い方】

• データファイルを読み込むには…
  データファイルの場所へ移動。
  　例）cd 'C:\Users\ishii\Desktop'
• プロット方法
  ２次元プロットのコマンド　plot
  　例）plot "data.txt"
  ３次元プロットのコマンド　splot
  　例）splot "data.txt"
• plotコマンドのオプション
  1) プロットの範囲を指定したい場合、
  　例として、x,y方向はデフォルトのまま。z方向だけ、0から0.05の範囲でプロットするには、
  splot [] [] [0:0.05] "input.dat"
  2) 線でプロットしたい場合、
  splot "input.dat" with lines
  3) 視点を変更したい場合、
  set view 60,20
  splot "input.dat"
  4) プロットした図をファイル出力したい場合。
  デフォルトのterminalは、パソコンのディスプレイ(wxt)なので、
  出力先のterminalを、postscript eps に変更する必要あり。
  そのために、次のコマンド入力。
  set terminal postscript eps
  続いて、出力ファイル名を入力。
  set output "sin.eps"
  最後に、出力させたいグラフを入力。
  plot sin(x)
  すると、sin曲線の書かれたsin.epsファイルが出来上がる。

Intel Visual Fortran

初期設定の問題点
デフォルトでは、スタック領域のサイズが小さすぎて、プログラムを実行させると以下のエラーが出る。
severe(170): Program Exception - stack overflow

スタック領域のサイズの変更方法
Intel Visual Studio の右にある「Solution explore」で右クリック。
そこから -> properties　-> linker　-> systemを選択。

stack reserve size 0
stack commit size 0

の数値を0から、100000000に変更すると、スタック領域を100MB確保することになる。
ちなみに、commit size は、ソフト起動時に予約する容量(らしい)。
ソフトの起動は遅くなるが、スタックの確保に時間はかからない。

GAMESSのホームページ（Mark Gordon's Quantum Theory Group）

入力データ

• 入力ファイルの例　（エチレン分子）
  

  指定するパラメータは、$CONTROL、$SYSTEMなどのGROUPに分けて記述。その書式は、必ず、$(GROUP名) ～ $END の形。
  “!”で始まる行は、コメント文。

  
  !Ethylene
  !B3LYP-D1/6-31G(d)
  !
    $CONTRL SCFTYP=RHF RUNTYP=ENERGY MULT=1 MAXIT=200 UNITS=ANGS NPRINT=3 $END
    $SYSTEM MWORDS=400 $END
    $DFT DFTTYP=B3LYP $END
    $BASIS GBASIS=N31 NGAUSS=6 NDFUNC=1 $END
    $GUESS GUESS=HUCKEL $END
    $SCF DIRSCF=.TRUE. $END
    $DATA
  ETHTLENE/B3LYP-D1/6-31+G(d)
  C1
  C 6 0.000000 0.000000 0.655000
  C 6 0.000000 0.000000 -0.655000
  H 1 0.000000 0.912325 1.214073
  H 1 0.000000 -0.912325 1.214073
  H 1 0.000000 -0.912325 -1.214073
  H 1 0.000000 0.912325 -1.214073
    $END
  
  
• 入力ファイルの意味
  
  

  $CONTRL	ICHARGE	系全体の電荷 =0, 中性 (default) =-1, －１価の陰イオン？ =+1, ＋１価の陽イオン？
  	MULT	系全体のスピン多重度。多重度は、不対電子のスピンをSとすると、2S+1で定義。 =1, Singlet (default)、不対電子無し。S=0。 =2, Doublet、不対電子１個。S=1/2。 =3, Triplet、不対電子２個。S=1。
  	DFTTYP	汎関数の種類 =None, Ab initio 計算 (分子軌道法を真面目に解く。) (default) =B3LYP, 密度汎関数法 (ハイブリッド汎関数B3LYPを使用。）
  	SCFTYP	計算の方法 =RHF, 制限 (閉殻系；一つの軌道に電子を反平行に２つずつ詰める。) =UHF, 非制限 (閉殻系 or 励起状態；一つの軌道に電子を１つずつ詰める。αスピンとβスピンの軌道に差が出る。）
  	RUNTYP	計算目的の指定 =ENERGY, シングルポイントエネルギー計算 (default) =OPTIMIZE, 構造最適化 =HESSIAN, 振動計算(調和振動解析)　($FORCE と $CPHF を作成) =MD, 分子動力学計算 ($MD を作成)
  	COORD	$DATAグループに格納する座標のフォーマット指定。 =UNIQUE, Cartesian座標 =ZMT, Gaussian形式のZ-matrix =ZMTMPC, MOPAC形式のZ-matrix
  	EXETYP	計算実行の種類を指定。 =RUN, 計算実行(default) =CHECK, Input と 必要なメモリーのチェック =DEBUG, デバッグ（出力可能な全てのデータを出力） =ROUTINE, 指定したroutineによるデータを全て出力
  	MAXIT	SCF計算の繰り返しの最大回数。 =30, (default)
  	UNITS	長さの単位。 =ANGS, Angstrom (default) =BOHR, Bohr atomic unit
  	NPRINT	出力データの条件。 =-7, extra printing from Boys localization. =-6, debug for geometry searches =-5, minimal output =-4, print 2e-contribution to gradient. =-3, print 1e-contribution to gradient. =-2, normal printing, no punch file =1, extra printing for basis,symmetry,ZMAT =2, extra printing for MO guess routines =3, print out property and 1e- integrals =4, print out 2e- integrals =5, print out SCF data for each cycle. Fock and density matrices, current MOs =6, same as 7, but wider 132 columns output. This option isn't perfect. =7, (default) normal printing and punching =8, more printout than 7. The extra output is (AO) Mulliken and overlap population analysis, eigenvalues, Lagrangians, ... =9, everything in 8 plus Lowdin population analysis, final density matrix.
  	QMTTOL	誤差の閾値。 =1.0E-6, (default) multi-ring moleculeに沢山のdiffuse関数を使う場合、SCF計算の収束を得るために、QMTTOL=3.0E-6に変更した方が良い(？)
  $SYSTEM	TIMLIM	計算時間の上限(単位：分) =525600, 分かりやすく言うと１年。(default)
  	MWORDS	確保するメモリサイズ(単位：MWord) =1, (default) ※1words=64bit
  $STATPT	NSTEP	構造最適化の際の探索回数の上限を指定。
  	OPTTOL	収束判定の際の閾値を指定。(単位 Hartree/Bohr) =0.0001, (default)
  $DFT	DC	Grimme's empirical dispersion correction (R^-6項)をするか否か。即ちDFT-D計算をするか否か。 =.FALSE., DFT-D計算しない。(default, 但しDFTTYP=B97-D, wB97X-D以外の場合 =.TRUE., DFT-D計算する。
  	IDCVER	DFT-Dのversionを指定。 =1, 1st in 2004; DFT-D1 =2, 2nd in 2006; DFT-D2 =3, 3rd in 2010; DFT-D3
  $BASIS	GBASIS	基底関数の基本セットの選択。 =STO, N21, N31, N311, G3L, G3LX, MINI, MIDI, DZV, DH, TZV, MC, CCn, ACCn, ... 他にもたくさん
  	NGAUSS	Gaussian関数の数を指定。(GBASIS=STO,N21,N31,N311のとき指定) =(整数値)
  	NDFUNC	加えるd-分極関数の数を指定。(重い原子？) 普通 3以下 =0, (default)
  	NFFUNC	加えるf-分極関数の数を指定。(Li - Cl) 普通、0 or 1 =0, (default)
  	FPFUNC	加えるp-分極関数の数を指定。(H, He) 普通 3以下 =0, (default)
  	DIFFSP	sp-diffuse関数を加えるか否か。flag to add diffuse sp (L) shell to heavy atoms. Heavy means Li-F, Na-Cl, Ga-Br, In-I, Tl-At. =.FALSE., 加えない。(default) =.TRUE., 加える。
  	DIFFS	s-diffuse関数を加えるか否か。flag to add diffuse s shell to hydrogens. =.FALSE., 加えない。(default) =.TRUE., 加える。
  $GUESS	GUESS	初期軌道の作成方法の指定。 =HUCKEL, 拡張Huckel法を使用。
  $SCF	DIRSCF	Direct SCF法を使うか否かの指定。 =.TRUE., 使う。
  $DATA	ETHYTLENE/B3LYP-D1/6-31+G(d) C1 C 6 0.000 0.000 0.655 C 6 0.000 0.000 -0.655 H 1 0.000 0.912 1.214 H 1 0.000 -0.912 1.214 H 1 0.000 -0.912 -1.214 H 1 0.000 0.912 -1.214	←　コメント文 ←　点群の指定　(通常 C1) ←　元素記号　原子番号　X Y Z (in Angstrom)　※Cartesian座標の場合



• 良く使われる B3LYP/6-31G(d) を指定するには
  $CONTRL DFTTYP=B3LYP
  $BASIS NGAUSS=6 GBASIS=N31 NDFUNC=1
  

実行の仕方

上記の入力ファイルを ethylene.inp、出力ファイル先を ethylene.log と指定する場合、

rungms ethylene.inp >& ethylene.log

で計算実行。

注）クラスター計算機で複数の計算を同時に行う場合、同じ名前の入力ファイルを使うと、 同一の計算nodeに投入された計算は、エラーを起こして停止するので注意。

出力データ

出力ファイルの最後の部分に、“exited gracefully.”と表示されていれば、とりあえず計算は正常終了。

• 出力ファイルの例

  　　　　　　　　　　　　　（途中まで省略）
  
  以下がTotal energy などの結果出力部分。
  ---------------------------------------------------------------------
  PROPERTIES FOR THE B3LYP DFT FUNCTIONAL (RHF TYPE) DENSITY MATRIX
  ---------------------------------------------------------------------
  
  -----------------
  ENERGY COMPONENTS
  -----------------
  
  WAVEFUNCTION NORMALIZATION = 1.0000000000
  
  ONE ELECTRON ENERGY = -170.0277680420
  TWO ELECTRON ENERGY = 58.2144599222
  NUCLEAR REPULSION ENERGY = 33.2774626572
  ------------------
  TOTAL ENERGY = -78.5358454626
  
  ELECTRON-ELECTRON POTENTIAL ENERGY = 58.2144599222
  NUCLEUS-ELECTRON POTENTIAL ENERGY = -247.7610113879
  NUCLEUS-NUCLEUS POTENTIAL ENERGY = 33.2774626572
  ------------------
  TOTAL POTENTIAL ENERGY = -156.2690888085
  TOTAL KINETIC ENERGY = 77.7332433459
  VIRIAL RATIO (V/T) = 2.0103250821
  
  　　　　　　　（途中省略）
  
  ----------------------------------------
  CPU timing information for all processes
  ========================================
  0: 46.676 + 0.101 = 46.778
  ----------------------------------------
  ddikick.x: exited gracefully.
  unset echo
  ----- accounting info -----
  Files used on the master node lion were:
  -rw-rw-r-- 1 ishii ishii 110467 Sep 25 15:45 /scr/ethtlene.dat
  -rw-rw-r-- 1 ishii ishii 586 Sep 25 15:45 /scr/ethtlene.F05
  -rw-rw-r-- 1 ishii ishii 1047040 Sep 25 15:45 /scr/ethtlene.F10
  -rw-rw-r-- 1 ishii ishii 1796336 Sep 25 15:45 /scr/ethtlene.F22
  Tue Sep 25 15:45:50 JST 2012
  0.115u 0.111s 0:50.37 0.4% 0+0k 0+0io 0pf+0w
  
  　　　　　　 （出力ファイルの最後）
  
  

• 出力ファイル中の単位
  エネルギー　Hartree （ 1 Hartree = 27.21138469 eV）

【要注意】GAMESSの実行時のエラー

(標準出力でのエラーメッセージ)
DDI Process 0: semget return an error.

入力ファイルにミスが無くても出るエラー。Semaphore Arraysを使い切ると出てくる。
GAMESSが異常終了すると、Semaphore arraysが解放されずに残ってしまうらしい。(参照URL：PCと遊ぶ日々の記録)

使用中(解放されてない) Semaphore arrays を確認するコマンド。
　ipcs -s
出力結果がずらずらと沢山出てきたら、怪しい。(以下は出力例)
　 ------ Semaphore Arrays --------
　key semid owner perms nsems
　0x00000000 10256412 root 666 1
　0x00000000 10253265 hishii 666 25

GAMESSを動かすには、２つの方法あり。

方法１：エラーを起こすノードを再起動。

方法２：再起動したくない場合、
次のコマンドで semid を指定して semaphore arrays を開放。
　ipcrm -s 10253265
沢山ある場合、面倒なので、次のスクリプトを利用すると楽。
　for i in `ipcs -s | sed -n '4,$p' | cut -d" " -f2 ` ;do ipcrm -s $i ;done
　計算クラスターの環境では、エラーを起こすすべてのノードでこの作業を行う必要あり。

GAMESSを利用時の注意

論文に下記の参考文献を引用しなければならない。

"General Atomic and Molecular Electronic Structure System"
M.W.Schmidt, K.K.Baldridge, J.A.Boatz, S.T.Elbert, M.S.Gordon, J.H.Jensen, S.Koseki, N.Matsunaga, K.A.Nguyen, S.Su, T.L.Windus, M.Dupuis, J.A.Montgomery,
J. Comput. Chem., 14, 1347-1363(1993).

• Maximaのホームページ（http://maxima.sourceforge.net/)
• Wikipedia (http://ja.wikipedia.org/wiki/Maxima)

グラフの描画

１つのコマンドの入力が終了したら、その都度、Shift+Enterでコマンド実行。

• 文字に数値代入、関数の定義。
  • 文字への数値を代入。コロン「:」を使う。（文字 a に 3 を代入）
    a:3
  • 関数も同様。（関数 f として、cos(x), 関数 g として、sin(y) を定義。）
    f:cos(x)
    g:sin(y)
    
  • 上で定義した文字 a や関数 f, g を使って、新たな関数も定義可能。（関数 h として、a * f * g を定義。）
    h:a*f*g
    即ち、定義した関数は、 h = 3 * cos(x) * sin (y) となっている。
• ２次元プロット
  • 上で定義した関数 f を、x=-πからx=+πの範囲で描画。
    plot2d(f,[x,-%pi,%pi])
  • 上で定義した関数 fとg を同時に、x=-πからx=+πの範囲で描画。
    plot2d([f,g],[x,-%pi,%pi])
• ３次元プロット
  • 上で定義した関数 h を、x=-πからx=+πの範囲、y=-πからx=+πの範囲で描画。
    plot3d(h,[x,-%pi,%pi],[y,-%pi,%pi])
  • グリッド数の指定。[オプション]
    plot3d(h,[x,-%pi,%pi],[y,-%pi,%pi],[grid,50,50])
  • 等高線の表示。[オプション]
    plot3d(h,[x,-%pi,%pi],[y,-%pi,%pi],[gnuplot_preamble,"set pm3d at bs"])