# dump2sq

## Table of Contents

- はじめに
- 入力ファイル
  - LAMMPS
  - XYZ
- 出力ファイル
  - 積算配位数
  - 動径分布関数
  - 構造因子
  - 原子散乱因子と中性子散乱長
- インストール
- 使い方
- tools
  - coeffploter.py
  - gr2sq.py
  - cn2Iq.c
  - cn2Iq.py

## はじめに

複数の時間ステップからなる時間ステップの構造データから、動径分布関数 \(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計算でセルサイズが変わったとしても問題なく計算できます。計算コ ストのほとんどは、原子間距離の計算です。

## 入力ファイル

拡張子で入力ファイルのformatを判断しています。

### LAMMPS

dumpコマンドで生成したファイルで、拡張子は.trajectryです。1ス テップ目で原子のならび決めて、2ステップ目以降は、同じ原子の並びとし て計算しています。時間ステップで原子の並びが変わると、正しい結果に ならないので注意してください。データの並びを決め打ちしています。そ のため、以下のようにdumpとdump_modifyで元素記号の出力し、idでsortし ておく必要があります。

    dump              1 all custom 10 hoge.lammpstrj id type element x y z
    dump_modify       1 element Si O sort id

### XYZ

xyz形式で、拡張子は.xyzまたは.10です。す。1行目に原子数、2行目に Step数と、格子定数を記述しておきます。alpha, beta, gammaがない場合 は、90, 90, 90とします。3行目以降にデータがならびます。2ステップ目 以降も同じ形式でデータを並べます。

    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
    .
    .
    .

## 出力ファイル

4つのファイルに拡張子を追加して出力します。basenameは、指定しなけれ ば、インプットファイルと同じになります。

### 積算配位数

拡張子はcnです。原子ペア毎の \(CN(r)\) が保存されています。

### 動径分布関数

拡張子はgrです。原子ペア毎の \(g(r)\) と全原子 \(G(r)\) が保存されています。

### 構造因子

拡張子はsqです。原子ペア毎の \(S(Q)\) とX線と中性子で観測される \(S(Q)\) が保存されています。

### 原子散乱因子と中性子散乱長

拡張子はcoeffです。中性子散乱長と原子散乱因子を出力します。原子散乱 因子は、xcoeff.fから抽出しています。中性子散乱長は、NISTのデータを 参照しています。

## インストール

プログラムはCで作成しています。Cのコンパイラを準備してmakeしてください。 ダウンロードして単にmakeしてください。

dump2sq-1.0.tar.gz

    make

## 使い方

並列計算のために、環境変数をセットしてスレッド数を指定しておく。

    export OMP_NUM_THREADS=4

example以下のglass.lammpstrjを引数で渡して実行してみます。

    dump2sq example/glass.lammpstrj

デフォルトで、\(g(r)\) を計算する\(r\) の最大値は、もっとも小さいセル長 さの半分としています。プログラムが実行されると、次のように出力されて、 4つのファイルが出力されます。

    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.

\(g(r)\) の\(r\) ステップや最大値、\(Q\) の範囲をコマンドラインで指定できます。詳細 は、ヘルプを参照してください。ヘルプで表示される[ ]の値はデフォルト 値です。GRmaxは、セルサイズを超えても計算できますが、物理的な意味は ないでしょう。

    ./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]

## tools

### coeffploter.py

元素を指定して、原子散乱因子と中性子散乱長をプロットする。引数に元 素名を複数していできる。

### gr2sq.py

Qの刻みを変えて再計算するために、dump2sqで生成したgrファイルから S(Q)を計算しなおす。結晶のパターンを計算するために、\(2\theta\) で\(Q\) を指定して計算することもできる。\(\theta\) と\(Q\) の関係は以下のとおり。 defaultで波長は、CuのK alpha線の1.54183 Åとしている。 \[ Q = \frac{4\pi}{\lambda}\sin(\theta) \] 計算した$S(Q)$はpngとxlsxを生成する

### cn2Iq.c

dump2sqで生成された*.cnと*.coeffから$I(q)$を計算するCプログラム。 openmpでスレッド並列する。f*fは足していない。sincのwindow関数をかけ る。コンパイルは、

    gcc -fopenmp -O3 -Wall cn2Iq.c  -o cn2Iq

実行は、*.cnと*.coeffを渡す。I(q)は標準出力に出力されうｒ．

    ./cn2Iq fort.cn fort.cn.coeff

複数のcnファイルやxslxに保存したい場合は、cn2Iq.pyを使う。

### cn2Iq.py

subprocessでCのcn2Iqを呼び出して、I(Q)のpngとxlsxを生成する。Qの範囲 またはthetaの範囲を指定できる。

Author: 大窪 貴洋

Created: 2023-04-24 月 09:01

Validate