はじめに

CrySPY: Python製の結晶構造探索ツール．https://tomoki-yamashita.github.io/CrySPY_doc/ja/index.html
ASE: 第一原理計算や古典力場計算の様ざなコードのインターフェースがあるPythonライブラリ．材料・物性の世界ではよく使われる．

今回はとりあえずのお試しで，簡単な古典力場計算（ASEに付属）を用いてCu 8原子の構造探索をランダムサーチで行う．計算はすぐ終わるものなので，手持ちのPCで実行できる．（windowsならWSLが必要）

Python3も必要．

インストール

Linux, Mac(x86 or rosetta)

pip3 install csp-cryspy

pipが使えればこのコマンドでCrySPYもASEも必要なものはインストールされる．
cryspyという実行スクリプトがpathが通ったところにインストールされる．

管理者権限がない場合などで，~/.local以下にインストールした場合は~/.local/binにpathを通しておかないとダメ．

Mac (arm64)

Apple siliconの環境では，pyshtoolsというライブラリのインストールが少し難しい（CrySPYが内部で使用しているPyXtalのインストールに必要）
こちらを参考：https://tomoki-yamashita.github.io/CrySPY_doc/ja/installation/cryspy/cryspy_1.0/index.html#arm64-on-macos-without-rosseta-2

チュートリアル

exampleファイル

https://tomoki-yamashita.github.io/CrySPY_doc/examples/ase_Cu8_RS.tar.gz
上記をクリックするか，下記コマンドでダウンロードする

wget https://tomoki-yamashita.github.io/CrySPY_doc/examples/ase_Cu8_RS.tar.gz

展開して，ディレクトリに移動する

tar zxf ase_Cu8_RS.tar.gz
cd ase_Cu8_RS

インプットファイル

cryspy.in

cryspy.inがCrySPYのインプットファイル．ランダムに5構造生成して，ASEを使って構造最適化，エネルギーを評価するような設定になっている．

[basic]
algo = RS
calc_code = ASE
tot_struc = 5
nstage = 1
njob = 2
jobcmd = zsh
jobfile = job_cryspy

[structure]
natot = 8
atype = Cu
nat = 8

[ASE]
ase_python = ase_in.py

[option]

calc_in/

calc_in/ の中にASEによる構造最適化のインプットとジョブファイルがある．
CrySPYでは構造最適化のインプットファイルにはxxx_1, xxx_2のようにステージに応じて番号をつけることになってる．（今回は1ステージなのでxxx_1のみ）

ase_in.py_1

from ase.constraints import ExpCellFilter, StrainFilter
from ase.calculators.emt import EMT
from ase.calculators.lj import LennardJones
from ase.optimize.sciopt import SciPyFminCG
from ase.optimize import BFGS
from ase.spacegroup.symmetrize import FixSymmetry
import numpy as np
from ase.io import read, write

# ---------- input structure
# CrySPY outputs 'POSCAR' as an input file in work/xxxxxx directory
atoms = read('POSCAR', format='vasp')

# ---------- setting and run
atoms.calc = EMT()
atoms.set_constraint([FixSymmetry(atoms)])
atoms = ExpCellFilter(atoms, hydrostatic_strain=False)
opt = BFGS(atoms)
#opt=SciPyFminCG(atoms)
opt.run()

# ---------- opt. structure and energy
# [rule in ASE interface]
# output file for energy: 'log.tote' in eV/cell
#                         CrySPY reads the last line of 'log.tote'
# output file for structure: 'CONTCAR' in vasp format
e = atoms.atoms.get_total_energy()
with open('log.tote', mode='w') as f:
    f.write(str(e))

write('CONTCAR', atoms.atoms, format='vasp')

上記はASEに付属のEMT（精度は良くない）で構造最適化する例．CrySPYではエネルギーをeV/cellの単位でlog.toteというファイルに，最適化後の構造情報をVASPフォーマットでCONTCARに出力する約束になっている．

calculatorは簡単に変えられるので，OpenKIMや機械学習ポテンシャルなんかも対応できるはず．

job_cryspyの中身

#!/bin/sh

# ---------- ASE
python3 ase_in.py

# ---------- CrySPY
sed -i -e '3 s/^.*$/done/' stat_job

ase_in.py（内部ではxxx_1の番号部分は自動的に削除されてコピーされる）をpythonで実行しているだけのスクリプト．CrySPYでは最後の1行を書く約束(ジョブが終わったかどうかの判断に必要)．

cryspy.inでjobcmd = zsh，jobfile = job_cryspyとした場合，CrySPYは内部でzsh job_cryspyを別プロセスとして実行している．

CrySPYの実行

cryspy.inと同じディレクトリで，cryspyと打ってEnter

cryspy

すると初回はこんな感じで5つ構造生成されるはず．

[2023-07-11 15:53:46,821][cryspy_init][INFO] 


Start CrySPY 1.2.0


[2023-07-11 15:53:46,821][cryspy_init][INFO] # ---------- Read input file, cryspy.in
[2023-07-11 15:53:46,822][read_input][INFO] Save input data in cryspy.stat
[2023-07-11 15:53:46,823][cryspy_init][INFO] # ---------- Initial structure generation
[2023-07-11 15:53:46,823][cryspy_init][INFO] Number of MPI processes: 1
[2023-07-11 15:53:46,823][gen_init_struc][INFO] # ------ mindist
[2023-07-11 15:53:46,826][struc_util][INFO] Cu - Cu: 1.32
[2023-07-11 15:53:46,826][gen_init_struc][INFO] # ------ generate structures
[2023-07-11 15:53:46,840][gen_pyxtal][WARNING] Compoisition [8] not compatible with symmetry 170: spg = 170 retry.
[2023-07-11 15:53:46,868][gen_pyxtal][INFO] Structure ID      0 was generated. Space group: 218 --> 223 Pm-3n
[2023-07-11 15:53:46,910][gen_pyxtal][INFO] Structure ID      1 was generated. Space group:   6 -->   6 Pm
[2023-07-11 15:53:46,925][gen_pyxtal][INFO] Structure ID      2 was generated. Space group:  17 -->  17 P222_1
[2023-07-11 15:53:46,969][gen_pyxtal][INFO] Structure ID      3 was generated. Space group:  83 -->  83 P4/m
[2023-07-11 15:53:47,179][gen_pyxtal][INFO] Structure ID      4 was generated. Space group:  39 -->  39 Aem2
[2023-07-11 15:53:47,255][cryspy_init][INFO] Elapsed time for structure generation: 0:00:00.432064

./data/init_POSCARSに生成された構造の情報がテキストで出力されるのでこれはVESTAで読み込んで可視化することができる．例えばID 0はこのような構造．

もう一度cryspyを実行すると，生成した構造の最適化に移る．なお，CrySPYはcryspy.statというファイルが存在すれば途中再開，なければゼロからスタートする．

今回はcryspy.inでnjob = 2という設定にしているので2つずつ構造最適化のジョブを投げる．このジョブは別プロセスとして裏で実行されるので，CrySPYはジョブを投げたら直ちに一度終了する．

cryspy

[2023-07-11 16:02:25,591][cryspy_restart][INFO] 


Restart CrySPY 1.2.0


[2023-07-11 16:02:25,608][ctrl_job][INFO] # ---------- job status
[2023-07-11 16:02:25,645][ctrl_job][INFO] ID      0: submit job, Stage 1
[2023-07-11 16:02:25,654][ctrl_job][INFO] ID      1: submit job, Stage 1

先ほど投げたジョブが終わった頃（今回は軽い計算なので一瞬で終わる）にcryspyをもう一度実行すると結果を回収して次のジョブを投げる

cryspy

[2023-07-11 16:04:05,735][cryspy_restart][INFO] 


Restart CrySPY 1.2.0


[2023-07-11 16:04:05,740][ctrl_job][INFO] # ---------- job status
[2023-07-11 16:04:05,740][ctrl_job][INFO] ID      0: Stage 1 Done!
[2023-07-11 16:04:05,775][ctrl_job][INFO]     collect results: E = 0.09291249537112356 eV/atom
[2023-07-11 16:04:05,805][ctrl_job][INFO] ID      1: Stage 1 Done!
[2023-07-11 16:04:05,806][ctrl_job][INFO]     collect results: E = 0.09295598446685593 eV/atom
[2023-07-11 16:04:05,819][cryspy][INFO] 

recheck 1

[2023-07-11 16:04:05,819][ctrl_job][INFO] # ---------- job status
[2023-07-11 16:04:05,827][ctrl_job][INFO] ID      2: submit job, Stage 1
[2023-07-11 16:04:05,835][ctrl_job][INFO] ID      3: submit job, Stage 1

全部で構造は5つなので最後（ID 4）もう1回実行する

cryspy

[2023-07-11 16:05:15,801][cryspy_restart][INFO] 


Restart CrySPY 1.2.0


[2023-07-11 16:05:15,806][ctrl_job][INFO] # ---------- job status
[2023-07-11 16:05:15,806][ctrl_job][INFO] ID      4: Stage 1 Done!
[2023-07-11 16:05:15,841][ctrl_job][INFO]     collect results: E = 0.12635875495049875 eV/atom
[2023-07-11 16:05:15,859][cryspy][INFO] 
Done all structures!

結果は./data/cryspy_rslt_energy_ascなどのファイルに出力される

cat data/cryspy_rslt_energy_asc

   Spg_num Spg_sym  Spg_num_opt Spg_sym_opt  E_eV_atom  Magmom      Opt
2       17  P222_1           64        Cmce   0.014612     NaN  no_file
0      223   Pm-3n          223       Pm-3n   0.092912     NaN  no_file
1        6      Pm            6          Pm   0.092956     NaN  no_file
4       39    Aem2           39        Aem2   0.126359     NaN  no_file
3       83    P4/m           83        P4/m   0.695032     NaN  no_file%

最適化後の構造は./data/opt_CIFS.cifまたは./data/opt_POSCARSに入っていて，どちらもVESTAで可視化できる．今回1番エネルギーが低かったID 2の構造は下記の通り．

5構造ではFCCのような安定構造に辿りつけなかったが，構造数を増やせばそのうち出てくる．