RMC計算用プログラム

目次

はじめに

SPring8やJ-PARCで実験を行って構造因子を取得したら、Reverse Monte Carlo 法(RMC)を使って分子構造を評価します。O氏が姫路で遊んでいる間と新幹線の 移動中にRMC計算用のプリポストプログラムを作成しました。学内からのみプ ログラムをダウンロードできます。

  • 複数データのフィット等、いくつかのプログラムをupdateしました。 (20-Apr-2020)
  • O氏に関連する文章を修正しました。 (27-Nov-2021)

RMCプログラム

RMCプログラムはRMC++を利用します。マニュアルを参照してインプットファイ ルを理解します。https://www.szfki.hu/~nphys/rmc++/opening.html

  • wsl上ですぐに実行できるようインテルコンパイルで作成したコンパイル済 のrmc_pot-1.4とpython一式を含むパッケージを作成しました。次のコマン ドでlinux用の実行ファイル一式を$HOME/dist以下に展開してインストールし ます。自前で準備したPythonとrmc_pot-1.4を使う場合、インストールの必 要はありません。

    curl https://amorphous.tf.chiba-u.jp/memo.files/rmc/linux_dist.tgz | tar xzf - -C ~/
    

    .bashrcにパス等の環境変数をセットしておきます。

    alias emacs='emacsclient -a "" -n -c'
    alias ls="ls -F --color"
    umask 022
    export PS1='\u@\h[\w] '
    export PATH=$HOME/bin/rmc:$PATH
    export MSROOT=$HOME/dist
    export PATH=$MSROOT/bin:$PATH
    export PATH=$MSROOT/py:$PATH
    export PATH=$MSROOT/miniconda3/bin:$PATH
    
  • 自作のpython解析プログラムをインストールします。$HOME/bin/rmc以下に ファイルが展開されます。

    curl https://amorphous.tf.chiba-u.jp/memo.files/rmc/rmcpy.tgz | tar xzf - -C ~/
    
  • 練習用に解析ファイル一式をダウンロードします。$HOME/example以下にデー タの入った00AとSiO2が展開されます。
curl https://amorphous.tf.chiba-u.jp/memo.files/rmc/rmc_example.tgz | tar xzf - -C ~/

測定データの確認

rmc_datplot.py

測定で得られたdatファイルをとりあえずプロットして確認するプログラムです。

rmc_datplot.py 00A_altered-2.dat 

rmc_datplot.png

図1: rmc_fitplot.pyの出力例

rmc_dat2interpolate.py

ビームが落ちた等で強度が妙なデータを前後のポイントで内挿して補完するプ ログラムです。補完するポイントと補完のための前後のポイント数を指定しま す。100-200のデータを表示して、111と112のポイントを106-110と113-117の ポイントで内挿して補完する例です。

rmc_dat2interpolate.py -x 100 200 -p 111 112 -s 5 04A_altered-2.dat -xl 100 120

補完したデータは、04A_altered-2_interp.datに保存されます。

04A_altered-2_interp.dat was created.
04A_altered-2_interp.png was created.

rmc_dat2interpolate.png

図2: rmc_dat2interpolate.pyの出力例

RMC用のインプットファイル作成

古典MDで構造をつくって、BL02B2で測定した散乱因子からインプットファイ ルを作成します。手順は以下のとおり。

  • 古典MDで作成した原子構造(lammpstrj)からrmc用のcfgファイルを作成 (rmc_dump2cfg.py)
  • 原子の数とセルサイズから原子を格子上にランダムに置いたrmc用のcfgファ イルを作成(rmc_makecfg.py)
  • 実験データをrmc用のフォーマットに変換 (rmc_cfg2sq.py)
  • rmc用のインプットファイルの作成(rmc_cfg2dat.py)

古典MDで構造が作れない場合は、rmc_makecfg.pyでランダムな初期構造kを作る。

rmc_convert2cfg.py

lammpstrjまたはxyz形式のデータをcfg形式に変換します。

rmc_dump2cfg.py 00A.lammpstrj

単位胞のサイズや原子種、原子座標の入った00A.cfgが生成されます。 xyzファイルの場合は、2行目にセルサイズを含む必要があります。

rmc_makecfg.py

引数に構成元素、格子サイズを指定してrmcのための初期構造を作成します。 セル中の原子位置は、格子状にサイトを生成してランダムに配置します。

rmc_makecfg.py  -l 34.350974 -c 678Si 360B 284Na 2038O  -o test.cfg

これで、下図のような1辺34.350974 Å のセルで678Si-360B-284Na-2038Oの 原子で構成されるtest.cfgが生成されます。

rmc_makecfg.png

rmc_cfg2sq.py

引数にcfgファイルと実験から求めたS(Q)のファイルを渡します。

rmc_cfg2sq.py -c 00A.cfg -s Sq_00A_pristine.txt

原子散乱因子を含むrmc用の00A_xrdraw.sqが生成されます。

rmc_cfg2dat.py

rmcのインプットファイルを作成します。引数にcfgファイルとrmc_cfg2sq.py で作成したS(Q)のファイルを渡します。

rmc_cfg2dat.py -c 00A.cfg -s 00A_xrdraw.sq

計算条件の00A.datが生成されます。このファイルを編集して束縛条件等の原 子を動かす条件を付け加えます。

計算の実行

作成したインプットファイルで正しくrmcが動作するか確認するために、1ステッ プのrmc計算を流します。

rmcp_onestep 00A

logファイルを見て、指定した通りの束縛条件になっているか確認します。指 定したとおりに実行できていたら、rmcを実行します。teeで標準出力に結果を 表示しつつ、標準出力を00A.logに保存します。00A.logは収束状況を確認する 時に使います。

rmcp 00A | tee 00A.log

ポスト処理

計算を流しながら、rmcが生成するファイルから収束状況やフィッテングの具 合、構造の妥当性をチェックします。計算結果を確認するためのポスト処理プ ログラムを作成しています。

rmc_fitplot.py

カレントディレクトリにあるfitファイルを読んでフィッテングの具合を確認 します。引数でファイルを指定しなければ、ディレクトリにあるfitファイル をプロットします。複数データをfitしている場合、全てのデータをプロット します。

rmc_fitplot.py

fitplot.png

図3: rmc_fitplot.pyの出力例

rmc_logplot.py

fitファイルを読んで収束状況確認します。引数を与えなければ、ディレクト リにあるlogファイルをプロットします。

rmc_logplot.py

logplot.png

図4: rmc_logplot.pyの出力例

rmc_pfqplot.py

部分S(Q)のpfqファイルを読んでを確認します。引数がなければ、カレント ディレクトリにあるpfqを読んで全ての部分S(Q)をプロットします。-fで pfqのファイルを指定でけきます。また特定のペアを指定してS(Q)をプロット できます。ペアは数字で指定します。例えば、1-1と1-3ペアのS(Q)をプロッ トしたい場合は次のようにします。

rmc_pfqplot.py -p 1-1 1-3

pfqplot.png

図5: rmc_pfqplot.pyの出力例

rmc_ppcfplot.py

部分G(r)のppcfファイルを読んでを確認します。引数がなければ、カレン トディレクトリにあるppcfを読んで全ての部分G(r)をプロットします。-f でppcfのファイルを指定でけきます。また特定のペアを指定してG(r)をプロッ トできます。ペアは数字で指定します。例えば、1-1と1-3ペアのG(r)をプロッ トしたい場合は次のようにします。

rmc_ppcfplot.py -p 1-1 1-3

ppcfplot.png

図6: rmc_ppcfplot.pyの出力例

rmc_cncheck.py

原子ペアとカットオフ距離から配位数分布をcfgファイルから求めます。引数 でcfgファイルの原子の並びを指定します。配位数を判定するための距離は、 rmc_check.pyの中で、以下のように辞書変数として直接指定しています。辞 書変数のキーが配位の種類、値が配位数をカウントする最小と最大距離で す。'O-Si,B,Al'の指定は、Oに対して、Si,B,Alの配位数をカウントするよう 指定しています。

CUTOFF['Si-O'] = [[1.0, 2.1]]
CUTOFF['B-O'] = [[1.0, 2.1]]
CUTOFF['Al-O'] = [[1.0, 2.3]]
CUTOFF['O-Si,B,Al'] = [[1.0, 2.1], [1.0, 2.1], [1.0, 2.3]]

プログラム中で指定したCUTOFF全てのペアについて配位数をカウントします。 配位数を計算するペアを追加したければ、rmc_check.pyの中にCUTOFFを追加 すれば良いでしょう。CUTOFFで指定した原子がcfg中に存在しない場合、エラー になります。解析毎にrmc_cncheck.pyを作成して解析すれば良いでしょう。

実行例です。

rmc_cncheck.py 00A.cfg -e Si B Na O

標準出力に結果が表示されます。

 Counting Si-O [[1.0, 2.1]]...
----------------------------------------
[cn]: N(%) Si for O
----------------------------------------
[3]:   26(3.8%)
[4]:  638(94.1%)
[5]:   14(2.1%)
Counting B-O [[1.0, 2.1]]...
----------------------------------------
[cn]: N(%) B for O
----------------------------------------
[2]:   11(3.1%)
[3]:  180(50.0%)
[4]:  168(46.7%)
[5]:    1(0.3%)
Counting O-Si,B [[1.0, 2.1], [1.0, 2.3]]...
----------------------------------------
[cn]: N(%) O for Si,B
----------------------------------------
[0]:    3(0.1%)
[1]:  154(7.6%)
[2]: 1823(89.5%)
[3]:   58(2.8%)

定義した辞書変数CUTOFF毎に配位数の数割合が表示されます。

rmc_cfg2convert.py

rmcのcfgファイルを各種フォーマットに変換します。cfgのデータに並びにした がった元素名を-eで指定する必要があります。-t で出力するフォーマットを 選べます。

rmc_cfg2convert.py hoge.cfg -e Si B O -t cif

指定できる出力フォーマットは、lammpstrj, cif, xyzです。

rmc_clean

中間ファイル一式を削除します。rmcpが生成するファイルをすべて削除する ので注意が必要です。

著者: 大窪 貴洋

Created: 2022-11-30 水 20:59

Validate