Phono3py hdf5文件数据读取与处理

Phono3py是一个主要用python写的声子-声子相互作用相关性质的模拟包，可以基于有限位移算法实现三阶力常数和晶格热导率的计算过程，同时输出包括声速，格林奈森常数，声子寿命和累积晶格热导率等参量。

相关介绍和安装请参考往期推荐。

phono3py快速安装教程

QE+phono3py 热导率计算

VASP+phono3py:快速计算晶格热导率

理论到实践：VASP+Phono3py计算Phonon Lifetime

phonopy 和 phono3py 安装教程Phonopy-Spectroscopy计算材料红外和Raman光谱

在phonopy和phono3py的计算过程中，会将过程文件和结果文件写入到hdf5的格式文件中，以此来方便程序读入和写入，phonopy有关的诸多输出量可以通过额外设置获得文本文件，并用于数据处理，phono3py则不能自主输出文本数据。但是在phono3py介绍网址上有如何读取hdf5文件的说明。

http://phonopy.github.io/phono3py/hdf5_howto.html

本文则以此来介绍读取phono3py计算完晶格热导率后的数据读出。

生成有限位移文件后，并通过VASP计算获得力常数fc3.hdf5和fc2.hdf5后，进行热导率计算，采用11×11×11的网格，

phono3py-load --mesh 11 11 11 --br

有关计算结果会被写入到 kappa-m111111.hdf5中，首先通过pip安装 ipython

pip install ipython 

随后打开ipython编辑器

ipython

随后输入（In [1]:是ipython的输入提示，只需要输入冒号后面的部分就行）

In [1]: import h5pyIn [2]: f = h5py.File("kappa-m111111.hdf5")In [3]: list(f)

则输出kappa-m111111.hdf5中包含的数据的title

Out[3]:['frequency', 'gamma', 'group_velocity', 'gv_by_gv', 'heat_capacity', 'kappa', 'kappa_unit_conversion', 'mesh', 'mode_kappa', 'qpoint', 'temperature', 'weight']

获得热导率张量

In [4]: f['kappa'].shapeOut[4]: (101, 6)

In [5]: f['kappa'][:]Out[5]:array([[  0.00000000e+00,   0.00000000e+00,   0.00000000e+00,          0.00000000e+00,   0.00000000e+00,   0.00000000e+00],       [  2.11702476e+05,   2.11702476e+05,   2.11702476e+05,          6.64531043e-13,   6.92618921e-13,  -1.34727352e-12],       [  3.85304024e+04,   3.85304024e+04,   3.85304024e+04,          3.52531412e-13,   3.72706406e-13,  -7.07290889e-13],       ...,       [  2.95769356e+01,   2.95769356e+01,   2.95769356e+01,          3.01803322e-16,   3.21661793e-16,  -6.05271364e-16],       [  2.92709650e+01,   2.92709650e+01,   2.92709650e+01,          2.98674274e-16,   3.18330655e-16,  -5.98999091e-16],       [  2.89713297e+01,   2.89713297e+01,   2.89713297e+01,          2.95610215e-16,   3.15068595e-16,  -5.92857003e-16]])

获得温度的张量

In [6]: f['temperature'][:]Out[6]:array([    0.,    10.,    20.,    30.,    40.,    50.,    60.,    70.,          80.,    90.,   100.,   110.,   120.,   130.,   140.,   150.,         160.,   170.,   180.,   190.,   200.,   210.,   220.,   230.,         240.,   250.,   260.,   270.,   280.,   290.,   300.,   310.,         320.,   330.,   340.,   350.,   360.,   370.,   380.,   390.,         400.,   410.,   420.,   430.,   440.,   450.,   460.,   470.,         480.,   490.,   500.,   510.,   520.,   530.,   540.,   550.,         560.,   570.,   580.,   590.,   600.,   610.,   620.,   630.,         640.,   650.,   660.,   670.,   680.,   690.,   700.,   710.,         720.,   730.,   740.,   750.,   760.,   770.,   780.,   790.,         800.,   810.,   820.,   830.,   840.,   850.,   860.,   870.,         880.,   890.,   900.,   910.,   920.,   930.,   940.,   950.,         960.,   970.,   980.,   990.,  1000.])

获得特定温度热导率和模式热导率（需提前定义）

In [7]: f['kappa'][30]Out[7]:array([  1.09089896e+02,   1.09089896e+02,   1.09089896e+02,         1.12480528e-15,   1.19318349e-15,  -2.25126057e-15])In [8]: f['mode_kappa'][30, :, :, :].sum(axis=0).sum(axis=0) / weight.sum()Out[8]:array([  1.09089896e+02,   1.09089896e+02,   1.09089896e+02,         1.12480528e-15,   1.19318349e-15,  -2.25126057e-15])

计算声子寿命

In [9]: g = f['gamma'][30]In [10]: import numpy as npIn [11]: g = np.where(g > 0, g, -1)In [12]: lifetime = np.where(g > 0, 1.0 / (2 * 2 * np.pi * g), 0)

声子群速和热容等数据处理方式同理。

具体修改请查看h5py使用方法和hdf5格式文件读取以及ipython命令使用方法。