VASPのサンプルやってみたシリーズ。今回はNudged Elastic Band（NEB）法。
Nudged Elastic Band法
結晶中や表面にいる原子が今いる位置から隣の位置に移動するする（拡散）ことを考える。原子は結晶中ではもっとも安定な場所（ポテンシャルの低い場所）に位置しており、隣の場所に移動するためには一度ポテンシャルの山を登ってからまた下らないといけない。下図のような凹みから凹みに移るイメージだ。
この時隣のサイトに移る頻度は、
(1) 
で表すことができ、ここで
は移動の活性化エネルギーである。この活性化エネルギーは隣に移動する際の山の高さに相当し、これを第一原理（平面波基底）で求める方法が、ナッジドエラスティックバンド（Nudged Elastic Band: NEB）法である。
VASPでNEB法で計算をするには、最初と最後の位置を表すPOSCARに加えて、経路の途中を何点か抜き出してその原子位置（POSCAR）を用意しておく必要がある。この経路の途中においた点をイメージと呼ぶ。このイメージ同士はバネで結ばれていて、各イメージで左右のバネ（下図の赤いバネ）にかかる力が同じになるように緩和計算が行われる。（多分）これだけだとエラスティックバンド法と呼ばれるもので、これだと最小のエネルギー経路（minimum energy path: MEP）を通らずに経路が近道を通ることがある（コーナーカッティング）という欠陥があるらしい。そこで経路に垂直な方向にもバネを付けて（下図の青いバネ）垂直方向にかかる力も0になるように緩和計算を行うのがNEB法である。（…という認識でいる。間違っていたら教えてほしい。）
とりあえずVASPのサンプルに沿ってやってみる。
このサンプルはプラチナ (Pt) の001表面にある原子が隣の場所に移動する場合を扱っている。（下図）
下は上から見た図。
Pt17がPt13を飛び越えて隣に行くのではなく、Pt17がPt13と位置を交換しながら隣に動くという感じだ。もちろんこれ以外にも色んな拡散経路が考えられるが、とりあえずサンプルではこれを例に計算している。
INCAR
INCARは構造緩和に必要なパラメータに、下記の2つのパラメータを付け加える
| SPRING: | バネ定数。デフォルトは-5。この値が0だと各イメージは経路にそって垂直にのみ緩和される。値がネガティブ（0）の場合は、NEB法になる。つまり値がプラスだとただのEB法。 |
| IMAGES: | 経路途中のイメージの数。最初と最後合わせてIMAGES+2個の状態を扱うことになる。 |
Sub directory
NEB法をする場合、IMAGEの数+2個のサブディレクトリを作っておく必要がある。例えばIMAGES = 2の時、00, 01, 02, 03と4つのディレクトリを作る。00は始点、03は終点である。各ディレクトリにはPOSCARを作って保存しておく。01, 02のPOSCARは自分で作ってもいいと思うが、VASPのサンプルページにあるinterpolatePOSCARというスクリプトを使えば作ってくれる。始点と終点のPOSCARとイメージの数を記載したファイルPOSCAR1_POSCAR2を作って、interpolatePOSCARの引数として与えてあげればよい。AWKはよくわからないが、パッと見、始点と終点の座標を等分して出力しているだけの単純なスクリプト。→その後自分でスクリプトを作りました。
INCAR, KPOINTS, POTCARはメインのディレクトリに置いておく。サブディレクトリに入れる必要はない。
実行
上の2つのことをすれば、メインのディレクトリでVASPを実行すれば計算してくれる。なんて簡単。
結果
計算はイメージに対して行うのみで、始点と終点は計算されない。始点と終点は動かさないんだから当然か。始点と終点はsurface/hollowを使って計算したものを用いた。
こうやって見ると、01と02の間が気になる。多分そこが一番エネルギーが高いところ。というわけで、このサンプルではイメージの数を奇数にすべきだったのだ。ちなみに00と01の差は0.85eVが得られた。
