Pymatgenチュートリアル① 化合物の情報を取り出す

はじめに

この文章を書いた理由

ここでは、材料化学の研究者の視点から、「Pymatgen」というライブラリの使い方をまとめていく。このライブラリそのものの説明に関しては、例えば以下のサイトに詳しいが、単純に実験家（要するに使い手側）を中心として考えると、「計算科学ベースでの無機化合物データベースを扱うためのPythonライブラリ」と捉えればいいと理解している。
MateriApps -pymatgen

ある程度使えるようになれば、無機材料の研究を行う上で役に立つツールになりそうだ。しかし、インターネット上には、実験系の研究者視点での情報が少ない。別の言い方をすると、コーディングが達者な人が書いたと思しき情報が多く、経験の多くが実験ベースになっている人間（私）には、なかなかとっつきにくかった。

そういう状況ではあったが、どうにか最低限の実務に活かせそうな程度は扱えるようになった。これまでに得た知見は、他に実験ベースで研究しており、しかしデータベースをうまく扱って研究に活かしたいと考えている人の役に立つだろう、と考え、このような文章を公開することにした。とはいえ、ネット受けする文章は書けないので、自分以外に1人か2人でも読んでもらえたらいいな、と思っている。どちらかと言うと、書くことで自分の勉強になるだろう、という気持ちが強い。

想定している読者

上記の通り、「実験ベースで材料の研究をしている人」を想定して書く。ただし、PCまたはMacの環境の整え方を解説するつもりはない（というか、それには能力が足りない）ので、そこは自力でできる人を想定している。

使用するバージョン

pymatgen 2022.0.8
どうも1~2年くらい前に大幅に改訂があって、それ以前の情報が使いにくくなっている。（よくよく見たら関数などの呼び出し方が変わっただけだったみたいだが、Python初心者には辛かった。）
ちなみに、作業はVS code上で行っている。実行結果は、まんまコピペしているので、環境によっては多少見た目が違っているだろうが、気にしなくていいと思う。

実際の操作

化合物の扱い方①　元素の情報を取得する

早速、シリコン（Si）の情報を取り出してみる。

#Element：元素を扱うための関数
from pymatgen.core.periodic_table import Element

#普通の元素記号でデータを持ってくることができる
si = Element('Si')
print(type(si))

これを実行し、以下の結果を得る。

<enum 'Element'>

つまり、'Element'という関数に、文字列で元素記号を入れると、'Element'という、pymatgen独自の型で、変数siがセットされた。この'Element'型のsiから、いろいろと情報を取り出すことができる。入っている情報のうち、基本的なものを一覧するには、'.data'を付けるといい。つまり、

display(type(si.data))
display(si.data)

を実行すると、以下のように、長々と情報が出てくる。

<class 'dict'>
{'Atomic mass': 28.0855,
 'Atomic no': 14,
 'Atomic orbitals': {'1s': -65.184426,
  '2p': -3.514938,
  '2s': -5.075056,
  '3p': -0.153293,
  '3s': -0.398139},
 'Atomic radius': 1.1,
 'Atomic radius calculated': 1.11,
（中略）
 'Ionization energies': [8.15168,
  16.34585,
  33.493,
  45.14179,
  166.767,
  205.279,
  246.57,
  303.59,
  351.28,
  401.38,
  476.273,
  523.415,
  2437.65815,
  2673.17755],
 'Electron affinity': 1.38952128}

シレッと冒頭にtypeを入れたが、見ての通り、hoge.dataは辞書型になっている。ということは、keyを指定すれば、情報を取り出すことができる。例えば、

print(si.data['Atomic mass'])
#実行結果：28.0855

といった要領。pymatgenには、独自のデータの取り出し方も準備されているが、こっちの方が汎用性があるように思うので、現時点では、私はこちらを積極的に使っている。ただし、上記の方法では、取り出せない情報も、Element型には含まれている。これを取り出すには、本家のページを参考にして、以下のように取り出す。

si.full_electronic_structure
#実行結果：[(1, 's', 2), (2, 's', 2), (2, 'p', 6), (3, 's', 2), (3, 'p', 2)]

要するに、上記のページ内で、propertyと書いてある後のフレーズを.（ピリオド）でElement型の変数の後につけると、その情報を取り出すことができる。以下のものとかはよく使いそう。

#金属かどうかを判定、bool型の値を返す
#他にもis_metalloidとかis_noble_gasなどがある
si.is_metal
#実行結果：False

#原子半径を調べる
si._atomic_radius
#実行結果：1.1

化合物の扱い②　組成式の情報を取得する

以上のように、元素の情報を取り出すことができた。次は、化合物から情報を取り出す。

#Composition：組成式を扱うための関数
from pymatgen.core.composition import Composition

#Compositionの中に、文字列で普通に組成式を入れればいい
mat_0 = Composition("Ti3O5")
print(type(ti3o5))
#実行結果：<class 'pymatgen.core.composition.Composition'>

Elementで、.（ピリオド）区切りでいろいろな情報を取り出したように、いろいろな情報を取り出すことができる。

#組成に関する基本的な情報を得る
print(ti3o5.to_data_dict)
#実行結果:{'reduced_cell_composition': Comp: Ti3 O5,
# 'unit_cell_composition': defaultdict(float, {'Ti': 3.0, 'O': 5.0}),
# 'reduced_cell_formula': 'Ti3O5',
# 'elements': ['Ti', 'O'],
# 'nelements': 2}

#組成式を元に、各元素の価数を予想してくれる
display(mat_0.to_data_dict)
#実行結果：({'Ti': 3.3333333333333335, 'O': -2.0},)

#含まれている元素をリストで取得する
print(mat_0.elements)
#実行結果：[Element Ti, Element O]

#Composition型を辞書型に直す。keywordsが元素の羅列になる
print(mat_0.to_reduced_dict)
print(mat_0.to_reduced_dict.keys())
#実行結果：defaultdict(<class 'float'>, {'Ti': 3.0, 'O': 5.0})
#dict_keys(['Ti', 'O'])

#最も簡単な組成式に直してくれる。
print(mat_0.reduced_formula)
#実行結果：Ti3O5

最後の例は、Ti3O5だとわかりにくいので、補足的に次のコードも記しておく。

#約分できる組成の化合物を持ってくる
mat_1 = Composition('Ti6O10')
print(mat_1)
#実行結果：Ti6 O10

print(mat_1.reduced_formula)
#実行結果：Ti3O5

#逆に、1以下の添字で表された化合物は、以下のように扱える。
mat_2 = Composition('Ti0.6O1')
print(mat_2)
#実行結果：Ti0.6 O1

print(mat_2.get_integer_formula_and_factor())
#実行結果：('Ti3O5', 0.19999999999999996)

ついでに、少数点を用いた式で記された化合物を、組成式に直す方法も記した。こちらの実行結果に含まれる「0.19999999999999996」は、「0.19999999999999996で割って、組成式を得た」という意味合いの数字。なお、こちらの方が汎用的と感じるかもしれないが、get_integer_formula_and_factor()の出力結果の'Ti3O5'は、文字列であって、Composition型ではないことに注意したい。

以上、今回まとめた知識を活用すると、次のように化合物を扱うことができる。

#よくある小数点の添字をもった化合物を扱う
#半端な組成の化合物を持ってくる
mat_3_0 = Composition('Ti2.98O5.10')

#添字の数字だけを取り出す
sub_num = mat_3_0.to_reduced_dict.values()
#添字の数字の中の最大値で添字を全て割る = 規格化
max_sub = max(sub_num)
std_sub = [i/max_sub for i in sub_num]

#規格化した添字を持つ式に直す
#もうちょいうまい書き方があるが、ここはこれでいいことにする
mat_3_1 = Composition(str(mat_3_0.elements[0]) + str(std_sub[0])
+str(mat_3_0.elements[1]) + str(std_sub[1]))
print(mat_3_1)

#モデルにしたい化合物の組成式を持ってきて、上と同じことをする
mat_4_0 = Composition('Ti3O5')
sub_num = mat_4_0.to_reduced_dict.values()
max_sub = max(sub_num)
std_sub = [i/max_sub for i in sub_num]

mat_4_1 = Composition(str(mat_4_0.elements[0]) + str(std_sub[0])
+str(mat_4_0.elements[1]) + str(std_sub[1]))
print(mat_4_1)

#2つの式が、ある許容範囲以内で同じものかを判別する。
#判断したい式.almost_equals(基準にする式)という式で、bool型の値を返す
#例えば、論文に登場した式が、あるデータベースに含まれるか判断したいときに使える
print(mat_4_1.almost_equals(mat_4_1))

#出力結果：Ti0.58431373 O1
#Ti0.6 O1
#True

かなりモッさい書き方になった。また後で同じようなことをするが、その時にはもうちょっとスマートに書く予定。