1.关于NPT系综
NPT系综中需要设置压力PSTRESS，在OUTCAR中，搜索关键词pressure会出现三个值
其中，in kB行是压力矩阵的元素，前三个是主对角元元素，后面三个是非对角元元素
Pullay stress 是设置的PSTRESS的值，也就是外部压力（pstress)
external pressure是in kB行三个主对角元素平均值和Pullay stress的差值，对应于压力的波动量，这是结构的压力。
由于有动能压力，导致内部会有一个内应力，动能压力总是正的，内应力也就是负数，对应于使得体积膨胀
kinetic pressure 动能也会有一个压力，total pressure是最终的总压力，检查NPT是否平衡，计算这个值随步数的波动情况就可以
当有温度时，kinetic pressure不为0 ，且温度越高，kinetic pressure越大，external pressrure的值倾向于保持为负的kinetic pressure，这样，保持total pressure = pullay stress，也就会导致热膨胀的出现
不同的计算精度（ENCUT和KPOINTS对同一个结构计算的external pressure不同，尤其是对于大体系来说ENCUT, 如果要计算NPT来得到热膨胀效应，就必须要用相同的参数，不然可能会出现没有热膨胀，反而在设置的温度下体积相比结构优化的0 K下降的现象
提取pressure的值
grep 'total pressure' OUTCAR | awk '{print NR, $4}' > pressure.dat
import numpy as np
import matplotlib.pyplot as plt
import matplotlib
matplotlib.use('Agg')
import argparse
parser = argparse.ArgumentParser(description='deal two column data')
parser.add_argument('filename')
parser.add_argument('xlabel')
parser.add_argument('ylabel')
args = parser.parse_args()
data_txt = args.filename
x_label = args.xlabel
y_label = args.ylabel
# print('='*40)
# filename = input('请输入文件名：').strip()
# x_label = input("请输入x轴标签（直接回车可跳过）: ").strip()
# y_label = input("请输入y轴标签（直接回车可跳过）: ").strip()
# print("="*40)
# data_txt = filename
# x_label = x_label
# y_label = y_label
data = np.loadtxt(data_txt)
plt.rcParams['font.size'] = 16
plt.rcParams['axes.linewidth'] = 2
plt.rcParams['axes.titlesize'] = 16
plt.rcParams['axes.labelsize'] = 20
plt.rcParams['axes.labelweight'] = 'bold'
plt.rcParams['xtick.labelsize'] = 16
plt.rcParams['ytick.labelsize'] = 16
plt.rcParams['legend.fontsize'] = 16
figure = plt.figure(figsize=(8,6),dpi = 300)
grid = figure.add_gridspec(1,1)
ax1 = figure.add_subplot(grid[0,0])
x = data[:,0]
y = data[:,1]
ax1.plot(x,y,linewidth=1.8)
plt.xlabel(x_label)
plt.ylabel(y_label)
plt.xticks(fontweight = 'bold')
plt.yticks(fontweight = 'bold')
plt.legend(prop = {'weight' :'bold'})
plt.tight_layout()
plt.savefig('plot.png',dpi = 300)
转载请注明来源 有问题可通过github提交issue