pymatgenを用いてVASP POSCARファイルからdensityを計算する

　<b>pymatgenを用いてお気軽に密度を求める</b>コードを作りました。秘匿すべきノウハウではないと判断し、ここに書き記すこととします。

　結果だけであれば<b>「1. まずは帰結」</b>をご覧ください。

　通常外部ライブラリを活用する場合、私は普段サンプルコードを探します。<b>APIリファレンスというものの見方が分からない状態が長かった</b>のですが、最近そのあたりの勘所を得るに至っています。今回は<b>pymatgenのAPIリファレンスを見ながら自力解決を行ったプロセス</b>というものをお見せできればと思います。今回は密度算出をテーマとしていますが、この思考回路があればそのときどきの問題に応じ自力解決を進めやすくなると考えています。ご興味持たれましたら<b>「3. API紐解きの思考回路」</b>をご笑覧ください。


# 1. まずは帰結

```python
from pymatgen.io.vasp.inputs import Poscar

### read POSCAR file (edit file name 'xxxxx.vasp' as you like)
poscar = Poscar.from_file('MgO_mp-1265_conventional_standard.vasp', 
                          check_for_POTCAR=False, read_velocities=False)

### convert into structure object
str = poscar.structure

### calculate density
print(str.density)
```

```:出力結果
3.4714308905526954
```

以上です。MgOに関するPOSCARファイルから密度（g/cm3）を得ました。

# 2. なぜこれをやろうと考えたか

　昨今、機械学習やシミュレーションやプログラミングで実現されることは山ほどあると思います。「やったらいいこと」は山のようにあふれていて、「（ほかのタスクを押しのけて）やらなければならないこと」を<b>能動的に選ぶ（whyにこだわる）</b>ことが重要な時代だと強く感じます。今回のワークのwhyは下記の通りです。

* 材料編集の影響を定性理解するため、母構造を編集（置換・欠損等）し第一原理計算を実行することがある。
* 大々的に編集を行った材料に対し第一原理計算や第一原理分子動力学法（AIMD）で尤もらしい平衡体積を得ることは容易でない。例えば、ベース構造から原子をいくつか抜くと密度が下がるが、そのような構造から最適化を実施しても密度が適切なレベルまで上がるとは限らない。
* そうしたとき、（実験などで見積もった）ターゲット密度を定め、構造最適化・AIMDの初期状態としてある程度格子定数等を調整した状態を用いることがある。
* 編集した構造（第一原理計算の入力構造）の密度がターゲットの密度とある程度一致するかを確認することがこのワークのモチベーション。
* VESTAやOVITOでは密度を算出する機能が見いだせなかった。
* この手のユーティリティコードを勤務先にはいろいろと保有していますが、休日にはアクセスできないため今回のコードを<b>趣味の範囲</b>で作成した（＝業務成果物ではない）。


# 3. API紐解きの思考回路

### 1. 仮説パイプラインの想定
　最初に以下の仮説を立てました。

　1. POSCARを読むためのクラスおよびメソッドがpymatgenにあるだろう。
　2. pymatgen.core.structureオブジェクトに変換できるだろう。
　3. pymatgen.core.structure moduleには密度を求めるmethodがあるだろう。

　ここで振り返ってみて気づいたことは、こうした仮説を思いつくには<b>ある程度pymatgenを使う場数を踏んでいる</b>ことが必要だったということでした。要は<b>慣れ</b>です（残念ながら聖杯は提供できません）。いつも仮説のとおり進むわけではありませんが、今回は上記のステップで実際に完了します。

### 2. "pymatgen POSCAR"でGoogle検索（ステップ1）
　pymatgen.io.vasp.inputs moduleのページに当たります。
　https://pymatgen.org/pymatgen.io.vasp.inputs.html

　このページにおいて「read」などのwordでページ内検索をかけるとfrom_fileメソッドが見つかり、
```python
poscar = Poscar.from_file('MgO_mp-1265_conventional_standard.vasp', 
                          check_for_POTCAR=False, read_velocities=False)
```
が導かれます。

### 3. オブジェクトのタイプやメソッドを調べる（Poscarオブジェクト；ステップ2）
```python
print(type(poscar))
print(dir(poscar))
```
```:出力結果
pymatgen.io.vasp.inputs.Poscar
['REDIRECT',
 '__class__',
 '__delattr__',
 '__dict__',
  (中略)
 '__weakref__',
 'as_dict',
 'comment',
 'from_dict',
 'from_file',
 'from_string',
 'get_string',
 'natoms',
 'predictor_corrector',
 'predictor_corrector_preamble',
 'selective_dynamics',
 'set_temperature',
 'site_symbols',
 'structure',
 'temperature',
 'to_json',
 'true_names',
 'unsafe_hash',
 'validate_monty',
 'velocities',
 'write_file']
```

　<b>type()</b>はオブジェクトのタイプを返しますが、ライブラリ中の階層構造も表していることが重要です。<b>dir()</b>は当該オブジェクトで使用可能なメソッドの一覧を示します。「名は体を表す」という言葉があるように、
* as_dict: 内容を辞書で返す
* to_json: 内容をJSON形式で返す
* write_file: 内容をファイルに書き出す
 
といったことはおおよその見当がつきます。pymatgenではStructureオブジェクトが基本中の基本で、ここでもstructureメソッドを<b>「試し打ち」</b>する
```python
print(poscar.structure)
```
```:出力結果
（略）
```
```python
print(type(poscar.structure))
```
```:出力結果
pymatgen.core.structure.Structure
```
ことで、実際にpymatgenのStructureオブジェクトが得られていることが分かります。

　さて、いつもdir()で示されるメソッドの内容が想像できるわけではありません。そのときには先に示したAPIリファレンスページ（<href>https://pymatgen.org/pymatgen.io.vasp.inputs.html</href>）で当該クラスのメソッド（下図のグレーのボックス（structure, comment, true_names, ...）と付加的な説明）をパラパラとみていくだけで狙ったものに引っ掛かることがあります。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/267549/907d769a-e35d-95ad-69ab-bbb799693d78.png)
　このように、①次に行うこと（Structureへの変換）を仮説として持ち、②type()で今注目するオブジェクトの型を把握し、③dir()やAPIリファレンスにより当該クラスのメソッドをパラパラ読む、ということをすれば、目的を達成する方法にヒットする確率が高まります（「そのような仮説を満たす機能はない」ということが解である場合ももちろんあります）。

### 4. オブジェクトのタイプやメソッドを調べる（Structureオブジェクト；ステップ3） 
```python
str = poscar.structure
print(dir(str))
```
```:出力結果
['DISTANCE_TOLERANCE',
 'REDIRECT',
 '__abstractmethods__',
 '__class__',
（中略）
 'count',
 'density',
 'distance_matrix',
（中略）
 'unsafe_hash',
 'validate_monty',
 'volume']
```

　ここでも、仮説（Structureオブジェクトには密度を求めるメソッドがあるのでは？）がスタート地点です。dir()で出てくるメソッドは数が多いですが、密度（density）をゴールとして念頭においておけばdensityメソッドが最初にトライすべきものとしてアンテナにひっかかるのではないでしょうか。ちなみに、Structureオブジェクトについては<href>https://pymatgen.org/pymatgen.core.structure.html</href>に詳細が示されています。

### API紐解きに関するまとめ
　以上のように、「こういうパイプラインを組めそう」「こういう機能（メソッド）がありそう」という仮説をたて、オブジェクトの型把握やメソッドの一覧化・検索を行う、というアプローチが今のところ外部ライブラリを使いこなすための近道だと考えています。もっとよい方法があればぜひご教示ください。