4.1. 多元系イオン結晶 \({\rm Mg}{\rm Al}_2 {\rm O}_4\) の, Mg, Al 原子の反転度計算

このチュートリアルの入力ファイルは examples/standard/spinel にあります。 以下、QuantumESPRESSOおよびVASPでの例について順に記載します。

4.1.1. QuantumESPRESSO を用いた例

abICS の入力ファイルの準備

abICS の入力ファイルは input_qe.toml です。

[replica]
nreplicas = 2
nprocs_per_replica = 1

kTstart = 1000.0
kTend = 1200.0

nsteps = 5  # Number of steps for sampling
RXtrial_frequency = 1
sample_frequency = 1
print_frequency = 1

レプリカ交換モンテカルロ法に関する情報です。 この例では簡単のため、レプリカ数を2つに設定し、全ステップ数は 5 にしてあります。

[solver]
type = 'qe'
path= './pw.x'
base_input_dir = './baseinput'
perturb = 0.0
run_scheme = 'mpi_spawn'

第一原理エネルギー計算を行うソルバーの情報です。 この例では QuantumESPRESSO (QE) を用いるために type'qe' にしています。 QE のエネルギー計算ソルバー pw.x が、実行箇所にあることを path で指定しています(シンボリックリンクでも構いません)。 ソルバーごとに固有の入力パラメータファイルが収められたディレクトリとして、 base_input_dir を用いて ./baseinput を指定しています。 perturb は構造最適化のために原子位置を乱数でずらすパラメータですが、今回は構造最適化を行わないために 0 にしてあります(このパラメータは緩和過程で、高い対称性に起因した鞍点にトラップされてしまうのを避ける際に有用です)。 abICS は pw.xMPI_Comm_spawn で起動するため、 run_scheme として mpi_spawn を与えます。

[config]
unitcell = [[8.1135997772, 0.0000000000, 0.0000000000],
            [0.0000000000, 8.1135997772, 0.0000000000],
            [0.0000000000, 0.0000000000, 8.1135997772]]
supercell = [1,1,1]

[config] は原子位置を指定するセクションです。 unitcell および supercell でシミュレーションする空間を(オングストローム単位で)指定します。

[[config.base_structure]]
type = "O"
coords = [
     [0.237399980, 0.237399980, 0.237399980],
     [0.762599945, 0.762599945, 0.762599945],
     [0.512599945, 0.012600004, 0.737399936],
     [0.487399966, 0.987399936, 0.262599975],
     [0.012600004, 0.737399936, 0.512599945],
     [0.987399936, 0.262599975, 0.487399966],
     [0.737399936, 0.512599945, 0.012600004],
     [0.262599975, 0.487399966, 0.987399936],
     [0.987399936, 0.487399966, 0.262599975],
     [0.012600004, 0.512599945, 0.737399936],
     [0.487399966, 0.262599975, 0.987399936],
     [0.512599945, 0.737399936, 0.012600004],
     [0.262599975, 0.987399936, 0.487399966],
     [0.737399936, 0.012600004, 0.512599945],
     [0.237399980, 0.737399936, 0.737399936],
     [0.762599945, 0.262599975, 0.262599975],
     [0.512599945, 0.512599945, 0.237399980],
     [0.487399966, 0.487399966, 0.762599945],
     [0.012600004, 0.237399980, 0.012600004],
     [0.987399936, 0.762599945, 0.987399936],
     [0.987399936, 0.987399936, 0.762599945],
     [0.012600004, 0.012600004, 0.237399980],
     [0.487399966, 0.762599945, 0.487399966],
     [0.512599945, 0.237399980, 0.512599945],
     [0.737399936, 0.237399980, 0.737399936],
     [0.262599975, 0.762599945, 0.262599975],
     [0.237399980, 0.512599945, 0.512599945],
     [0.762599945, 0.487399966, 0.487399966],
     [0.762599945, 0.987399936, 0.987399936],
     [0.237399980, 0.012600004, 0.012600004],
     [0.737399936, 0.737399936, 0.237399980],
     [0.262599975, 0.262599975, 0.762599945],
     ]

[[config.base_structure]] は、モンテカルロ計算で動かさない原子を指定するセクションです。 type で原子名を指定し、 coords でその原子が位置する分率座標を指定します。 coords ではリストのリストだけではなく、空白区切り3列で、原子数だけ行があるような数値行列を表す文字列を与えることもできます。

[[config.defect_structure]]
coords = [
     [0.000000000, 0.000000000, 0.000000000],
     [0.749999940, 0.249999985, 0.499999970],
     [0.249999985, 0.749999940, 0.499999970],
     [0.249999985, 0.499999970, 0.749999940],
     [0.749999940, 0.499999970, 0.249999985],
     [0.499999970, 0.749999940, 0.249999985],
     [0.499999970, 0.249999985, 0.749999940],
     [0.000000000, 0.499999970, 0.499999970],
     [0.749999940, 0.749999940, 0.000000000],
     [0.249999985, 0.249999985, 0.000000000],
     [0.249999985, 0.000000000, 0.249999985],
     [0.749999940, 0.000000000, 0.749999940],
     [0.499999970, 0.000000000, 0.499999970],
     [0.000000000, 0.749999940, 0.749999940],
     [0.000000000, 0.249999985, 0.249999985],
     [0.499999970, 0.499999970, 0.000000000],
     [0.374999970, 0.374999970, 0.374999970],
     [0.624999940, 0.624999940, 0.624999940],
     [0.374999970, 0.874999940, 0.874999940],
     [0.624999940, 0.124999993, 0.124999993],
     [0.874999940, 0.874999940, 0.374999970],
     [0.124999993, 0.124999993, 0.624999940],
     [0.874999940, 0.374999970, 0.874999940],
     [0.124999993, 0.624999940, 0.124999993],
     ]

[[config.defect_structure]] はモンテカルロ計算で動かす原子の位置を指定するセクションです。

[[config.defect_structure.groups]]
name = 'Al'
# species = ['Al']    # default
# coords = [[[0,0,0]]]  # default
num = 16
[[config.defect_structure.groups]]
name = 'Mg'
# species = ['Mg']    # default
# coords = [[[0,0,0]]]  # default
num = 8

[[config.defect_structure.groups]] はモンテカルロ計算で動かす原子の名前と数を指定するセクションです。 Al 原子 16 個と、 Mg 原子 8 個を [[config.defect_structure]] で指定した座標中に動かします。

QE の入力ファイルの準備

abICS では、内部で保持した原子座標をもとに QE の入力ファイルを生成しますが、原子座標以外の情報、例えば擬ポテンシャルの情報は、ユーザが与える必要があります。 QE の場合は、この情報は base_input_dir で指定したディレクトリの中に scf.in という名前で、 pw.x で使うファイルをそのまま置きます。 このファイルをもとに、ユニットセルや原子座標などを書き換えた入力ファイルが自動生成されます。

  • 注意点

    • 擬ポテンシャルディレクトリ pseudo_dir は絶対パスで渡す必要があります。 行頭の ~ はユーザのホームディレクトリに展開されます。

    • calculation はSCF 計算 scf と構造最適化 relax のみが利用可能です。

    • disk_ionone にすると計算が進まなくなります。 low にしてください。

実行・解析

MPI 実行におけるプロセス数やオプションの設定については 基本的な使用方法 を参照してください。

mpiexec -np 2 abics input_qe.toml

abics 計算終了後、反転度(もともと Mg 原子があった位置にAl 原子が占めている割合)を計算します。 そのために examples/standard/spinel ディレクトリにある separateT.pycalc_DOI.py スクリプトを用います。

まず、 separateT.py を使い、レプリカごとに分かれていた構造情報を温度ごとにまとめます:

mpiexec -np 2 python3 ./separateT.py input_qe.toml

結果が Tseparate ディレクトリに格納されます。 このディレクトリが既にあった場合は、 Tseparate.DATE というディレクトリに退避されます(DATE は現在時刻)。

つぎに、 calc_DOI.py スクリプトで、反転度を計算します:

mpiexec -np 2 python3 ./calc_DOI.py input_qe.toml

結果が Tseparate/DOI_T.dat に出力されます。

../_images/doi_2.png

レプリカ数を増やすと様々な温度で計算ができます。 また、ステップ数を増やすことでよりきれいな結果を得られます。 モンテカルロ平均で最初に捨てるサンプル数は calc_DOI.py スクリプトにある throwout の値で指定できます。 実際の計算ではさらに、DFT 計算の前に構造最適化をする必要があります。 QE の場合には baseinput/scf.inCALCULATION パラメータを relax にしてください。 例えば、17 レプリカ、1000 ステップで、 perturb = 0.1 のもとで計算をすると

../_images/doi_qe.png

といった結果が得られます。 温度を上げるほど反転度が上がっていくのがわかります。

4.1.2. VASP を用いた例

abICS の入力ファイルの準備

abICS の入力ファイルは input_vasp.toml です。 以下、 examples/standard/spinel にある input_vasp.toml を例に説明します。 [solver] セクションのみ、QuantumESPRESSOと異なり、以下のように指定されています。

[solver]
type = 'vasp'
path = './vasp'
base_input_dir = './baseinput'
perturb = 0.0
run_scheme = 'mpi_spawn_ready'

VASP を用いるために type'vasp' にしています。 また、VASP のエネルギー計算ソルバー vasp が、実行箇所にあることを path で指定しています(シンボリックリンクでも構いません)。 ソルバーごとに固有の入力パラメータファイルが収められたディレクトリとして、 base_input_dir を用いて ./baseinput を指定しています。 perturb は構造最適化のために原子位置を乱数でずらすパラメータですが、今回は構造最適化を行わないために 0 にしてあります。 abICS は vaspMPI_Comm_spawn で起動するため、 run_scheme として mpi_spawn_ready を与えます(VASPをソルバーとして利用する際には、MPI_COMM_SPAWNを利用するためのパッチをあてる必要があります。利用されたい場合には、お問い合わせ のその他に記載された連絡先までご連絡ください)。

VASP の入力ファイルの準備

abICS では、内部で保持した原子座標をもとに VASP の入力ファイルを生成します。それ以外の情報は、ユーザが与える必要があります。 この情報は base_input_dir で指定したディレクトリの中に INCARPOSCARKPOINTS があります。 上記のファイル以外にも、 POTCAR がVASPの実行には必要となりますが、ライセンスの関係上、サンプルファイルにはおいてありません。 実行前には O, Al, Mg の擬ポテンシャルファイルから POTCAR ファイルを作成してください。 これらのファイルをもとに、ユニットセルや原子座標などを書き換えた入力ファイルが自動生成されます。

  • 注意点

    • POSCAR の座標情報はabICSの入力情報で上書きされますが、内部処理を行うために記載する必要があるのでご注意ください。

    • POTCAR は原子のアルファベット順に記載してください。

実行・解析

実行・解析の手順はQEの場合と同様です. QEの例と同じように、17 レプリカ、1000 ステップで計算をすると

../_images/doi_vasp.png

といった結果が得られ、同様の結果が得られることがわかります。