#+TITLE: dump2sq
#+AUTHOR: 大窪 貴洋
#+LANGUAGE: en
#+EMAIL: ohkubo.takahiro@faculty.chiba-u.jp
#+OPTIONS: toc:2 num:nil author:t creator:nil LaTeX:nil \n:nil ^:{}
#+SETUPFILE: https://amorphous.tf.chiba-u.jp/org-html-themes/org/theme-readtheorg.setup
#+MACRO: color @@html:$2@@
#+MACRO: kakomi @@html:
$1
@@
\\
# * Author
# [[https://amorphous.tf.chiba-u.jp/][大窪 貴洋]] (千葉院工)
* はじめに
複数の時間ステップからなる時間ステップの構造データから、動径分布関数
$g(r)$ 、積算配位数 $CN(r)$ 、構造因子 $S(Q)$ を計算するプログラムです。
Triclinicなセルでも計算できます。原子ペア毎の$S(Q)$ を計算した後、原子
散乱因子または中性子散乱長を使って、X-rayと中性子回折実験的で観測され
る $S(Q)$ を計算します。$g(r)$ は、セルサイズ $L/2$ を超えても計算でき
ます。$L/2\sqrt{3}$ まで確保しても良さそうです。産総研の土田さんにアド
バイスいただきました。
https://www.cmu.edu/biolphys/deserno/pdf/gr_periodic.pdf
1ステップ毎に$g(r)$ と$S(Q)$ を計算して最後に時間平均しています。よっ
て、NPT計算でセルサイズが変わったとしても問題なく計算できます。計算コ
ストのほとんどは、原子間距離の計算です。この計算はopenmpで並列化してい
ます。
* 入力ファイル
拡張子で入力ファイルのformatを判断しています。
** LAMMPS
dumpコマンドで生成したtrajectryファイルです。1ステップ目で原子のな
らび決めて、2ステップ目以降は、同じ原子の並びとして計算しています。
時間ステップで原子の並びが変わると、正しい結果にならないので注意し
てください。データの並びを決め打ちしています。そのため、以下のよう
にdumpとdump_modifyで元素記号の出力し、idでsortしておく必要がありま
す。
#+begin_example
dump 1 all custom 10 hoge.lammpstrj id type element x y z
dump_modify 1 element Si O sort id
#+end_example
** XYZ
xyz形式です。1行目に原子数、2行目にStep数と、格子定数を記述しておき
ます。alpha, beta, gammaがない場合は、90, 90, 90とします。3行目以降
にデータがならびます。2ステップ目以降も同じ形式でデータを並べます。
#+begin_example
420
Step 225837 18.4116484 18.4116484 18.4116484 90 90 90
Si 10.560073 17.438902 6.614775
O 17.812448 0.110128 0.3462461
.
.
.
#+end_example
* 出力ファイル
4つのファイルに拡張子を追加して出力します。basenameは、指定しなけれ
ば、インプットファイルと同じになります。
** 積算配位数
拡張子はcnです。原子ペア毎の $CN(r)$ が保存されています。
** 動径分布関数
拡張子はgrです。原子ペア毎の $g(r)$ と全原子 $G(r)$ が保存されています。
** 構造因子
拡張子はsqです。原子ペア毎の $S(Q)$ とX線と中性子で観測される $S(Q)$ が保存されています。
** 原子散乱因子と中性子散乱長
拡張子はcoeffです。中性子散乱長と原子散乱因子を出力します。原子散乱
因子は、[[https://www.szfki.hu/~nphys/rmc++/downloads.html][xcoeff.f]]から抽出しています。中性子散乱長は、[[https://www.ncnr.nist.gov/resources/n-lengths/][NIST]]のデータを
参照しています。
* インストール
プログラムはCで作成しています。Cのコンパイラを準備してmakeしてください。
ダウンロードして単にmakeしてください。
#+begin_example
make
#+end_example
icc用のMakefileも準備しています。
#+begin_example
make -f Makefile.icc
#+end_example
これで、dump2sqが生成されます。
* 使い方
並列計算のために、環境変数をセットしてスレッド数を指定しておく。
#+begin_example
export OMP_NUM_THREADS=4
#+end_example
example以下のglass.lammpstrjを引数で渡して実行してみます。
#+begin_example
dump2sq example/glass.lammpstrj
#+end_example
デフォルトで、$g(r)$ を計算する$r$ の最大値は、もっとも小さいセル長
さの半分としています。プログラムが実行されると、次のように出力されて、
4つのファイルが出力されます。
#+begin_example
basename: example/glass
Maximum cell length = 17.284100
Rmax = 8.642050
Start Gr and Sq calculation...
Calculation 0 steps...
Calculation 10 steps...
Calculation 20 steps...
Calculation 30 steps...
Calculation 40 steps...
Calculation 50 steps...
Calculation 60 steps...
Calculation 70 steps...
Calculation 80 steps...
Calculation 90 steps...
Total: 99 steps
example/glass.coeff was generated.
example/glass.cn was generated.
example/glass.gr was generated.
example/glass.sq was generated.
#+end_example
$g(r)$ の$r$ ステップや最大値、$Q$ の範囲をコマンドラインで指定できます。詳細
は、ヘルプを参照してください。ヘルプで表示される[ ]の値はデフォルト
値です。GRmaxは、セルサイズを超えても計算できますが、物理的な意味は
ないでしょう。
#+begin_example
./dump2sq -h
usage: dump2sq [options] inputfile
inputfile format
LAMMPS: *.lammpstrj
FEMTECK: *.10
options:
-h --help show this help message
--SQmin=%f Q minimum [0.2]
--SQmax=%f Q maximum [30]
--SQdel=%d Q delta [0.02]
--GRmin=%f r minimum [0.025]
-r --GRmax=%f r maximum [L(minium cell length)*0.5]
-s --GRdel=%d r delta [0.05]
-o --output=%s Output basename [basename]
#+end_example