FORCE_CONSTANTS中3阶力常数大小与原子间距的分析脚本

shengBTE所需的3阶力常数文件格式为FORCE_CONSTANTS_3RD, 想要得知所有力常数数值与原子间距的变化关系，或者输出更小截断半径下新的FORCE_CONSTANTS_3RD以验证热导率随截断半径的变化，需要分析FORCE_CONSTANTS_3RD内数据，这里提供一个PYTHON脚本，能够读取FORCE_CONSTANTS_3RD文件，并绘制力常数强度和原子间距的散点图，以及输出新的截断半径下的力常数文件，供大家参考，注意python的缩进。

#####################para_constants
bohr2ang=0.529177210
rhy2ev=13.60569193 #1 Ry = 1/2 Ha = 13.60569193 eV
#####################################
import os
import numpy as np
import matplotlib.pyplot as plt
######################
#输入1:力常数文件所在目录
filepath=R"填入总路径"
ifc3_std_filename=filepath+'/'+'FORCE_CONSTANTS_3RD'
#输入2：为了计算间距，需要输入uc内原子坐标索引(为XYZ下单位angstrom),同vesta用.xyz格式
uc_xyz=[[5.654830126,1.395610194,1.885231820], \
[7.273416028,4.199308702,1.885231820],\
[4.036048945,4.199421446,1.885231820],\
[5.654700060,5.133874941 , 3.114744400],\
[4.035983972,2.330176433, 3.114744400],\
[7.273416149 , 2.330176433 , 3.114744400],\
[1.885025911, 1.088320217 , 2.486000000],\
[9.424374089 , 5.441164917 , 2.514000000] \
]
#输入3：设置新的cutoff,输出相应的3阶力常数
ifc3_new_cut=6.0 #新的截断半径
ifc3_cut_filename=filepath+'/'+'FORCE_CONSTANTS_3RD_cut' #新的力常数文件
#############################
ifc3_std=[]
#3阶力常数的读取
ifc3_file=open(ifc3_std_filename,'r')
ifc3_line=ifc3_file.readlines()
ifc3_file.close()
#
#格式1行为总块数，空行，块(1序号2行3列超胞索引1行原胞原子索引27行4值)，每个块前有空行
#块数目
temp_line=ifc3_line[0].split()
ifc3_b_num=int(temp_line[0])
#读取块，转为同输入1的ifc_res列表
ifc3_sc=[] #字符串，1重为所有块，2重为2个相对晶格矢的行字符串，两个字符串
ifc3_name=[] #为原胞3原子索引，字符串
ifc3_val=[] #27数值
for m1 in range(len(ifc3_line)-1): #排除第一行的总数目
temp_line=ifc3_line[m1+1].split()
if temp_line and len(temp_line)==1:
#块的开始
temp_line1=ifc3_line[m1+1+3].strip()
ifc3_name.append(temp_line1)
temp_line2=ifc3_line[m1+1+1].strip()
temp_line3=ifc3_line[m1+1+2].strip()
ifc3_sc.append([temp_line2,temp_line3])
val_temp=[]
for m2 in range(27):
temp_line4=ifc3_line[m1+1+4+m2].split()
val_temp.append(float(temp_line4[3]))
ifc3_val.append(val_temp)
###
ifc3_std=[ifc3_b_num,ifc3_sc,ifc3_name,ifc3_val]
######################分析标准文件间距与
#ifc3_std=[ifc3_b_num,ifc3_sc,ifc3_name,ifc3_val]
ifc3_len=[] #1重为块索引，2重为2和3原子索引以及最大值，共3值
ifc3_max_val=[]
for k1 in range(ifc3_std[0]):
temp1=ifc3_std[1][k1][0].split()
temp2=ifc3_std[1][k1][1].split()
a2_sc=[float(temp1[0]),float(temp1[1]),float(temp1[2])]
a3_sc=[float(temp2[0]),float(temp2[1]),float(temp2[2])]
temp3=ifc3_std[2][k1].split()
a1_ind=int(temp3[0])-1
a2_ind=int(temp3[1])-1
a3_ind=int(temp3[2])-1
a21_xyz=[uc_xyz[a2_ind][0]-uc_xyz[a1_ind][0],uc_xyz[a2_ind][1]-uc_xyz[a1_ind][1],uc_xyz[a2_ind][2]-uc_xyz[a1_ind][2]]
a31_xyz=[uc_xyz[a3_ind][0]-uc_xyz[a1_ind][0],uc_xyz[a3_ind][1]-uc_xyz[a1_ind][1],uc_xyz[a3_ind][2]-uc_xyz[a1_ind][2]]
d21=[a2_sc[0]+a21_xyz[0],a2_sc[1]+a21_xyz[1],a2_sc[2]+a21_xyz[2]]
d31=[a3_sc[0]+a31_xyz[0],a3_sc[1]+a31_xyz[1],a3_sc[2]+a31_xyz[2]]
d21_len=(d21[0]**2+d21[1]**2+d21[2]**2)**0.5
d31_len=(d31[0]**2+d31[1]**2+d31[2]**2)**0.5
d_max=max([abs(d21_len),abs(d31_len)])
d_now=[d21_len,d31_len,d_max]
ifc3_len.append(d_now)
###找到每个块27值的最大值
val=[abs(x) for x in ifc3_std[3][k1]]
val_max=max(val)
ifc3_max_val.append(val_max)
###
######绘图
x_d=[x[2] for x in ifc3_len]
y_d=ifc3_max_val
cut_res=[x_d,y_d]
plt.scatter(x_d,y_d,s=5,c="r")
plt.xlabel('std_res:x_max(d12 vs d13) VS y_max(27-val)')
plt.show()
############筛选ifc3_std的三阶力常数，用ifc3_new_cut
#ifc3_std=[ifc3_b_num,ifc3_sc,ifc3_name,ifc3_val]
ifc3_cut_num=0
ifc3_sc_cut=[]
ifc3_name_cut=[]
ifc3_val_cut=[]
#ifc3_len--------ifc3_max_val
for k1 in range(ifc3_std[0]):
if ifc3_len[k1][2]<=ifc3_new_cut:
ifc3_sc_cut.append(ifc3_std[1][k1])
ifc3_name_cut.append(ifc3_std[2][k1])
ifc3_val_cut.append(ifc3_std[3][k1])
ifc3_cut_num=ifc3_cut_num+1
ifc3_cut=[ifc3_cut_num,ifc3_sc_cut,ifc3_name_cut,ifc3_val_cut]
#输出3阶力常数
ifc3_new_file=open(ifc3_cut_filename,'w')
#ifc3_tot=[ifc3_b_num,ifc3_sc,ifc3_name,ifc3_val]
print(" %d "%(ifc3_cut[0]),file=ifc3_new_file)
#print(" ",file=ifc3_new_file)
for kk1 in range(ifc3_cut[0]):
print(" ",file=ifc3_new_file)
print(" %d"%(kk1+1),file=ifc3_new_file)
print("%s"%(ifc3_cut[1][kk1][0]),file=ifc3_new_file)
print("%s"%(ifc3_cut[1][kk1][1]),file=ifc3_new_file)
print(" %s"%(ifc3_cut[2][kk1]),file=ifc3_new_file)
#27值
pp0=0
for pp1 in range(3):
for pp2 in range(3):
for pp3 in range(3):
be_str=' '+repr(pp1+1)+' '+repr(pp2+1)+' '+repr(pp3+1)+' '
be_val=ifc3_cut[3][kk1][pp0]
pp0=pp0+1
print("%s %.14f"%(be_str,be_val),file=ifc3_new_file)
ifc3_new_file.close()

FORCE_CONSTANTS中3阶力常数大小与原子间距的分析脚本相关推荐

R语言ggplot2可视化：jupyter中设置全局图像大小、jupyter中自定义单个ggplot2图像结果的大小
R语言ggplot2可视化:jupyter中设置全局图像大小.jupyter中自定义单个ggplot2图像结果的大小目录
html里字号rem,html中的em和rem到底该如何使用，自适应效果中如何确定文字大小/字号？...
如今手机屏幕繁多,自适应效果中如何确定文字大小/字号? em rem vm vw vh你都了解吗? 先说说em和rem em:继承父级的,假设html的font-size默认为16px,body字体大 ...
android auto 字体大小,Android中给TextView字体大小用dp还是sp？dp和sp有什么区别？
dp,sp都会根据屏幕ppi显示大小不同. ppi的运算方式是: PPI = (开根号(长度像素数² + 宽度像素数²)) / 屏幕对角线英寸数 dp:Density-independent pixe ...
oracle设置缓存大小设置多少,【数据库类※从V$DB_CACHE_ADVICE中设置数据缓冲大小※】...
[数据库类※从V$DB_CACHE_ADVICE中设置数据缓冲大小※] Oracle 9i引入了一个新的途径来预测数据缓冲cache中附加数据缓存的所带来的好处的多少.V$DB_CACHE_A ...
深度机器学习中的batch的大小对学习效果有何影响？
深度机器学习中的batch的大小对学习效果有何影响? 如题,在深度学习中,刚入门的小弟一直听闻一个batch中同时训练多个数据可以得到较好的效果,于是小弟在caffe上跑deepID的网络时对如何选取 ...
5 控件固定大小_【聊技术】在Android中实现自适应文本大小显示
本周的聊技术话题和大家说说如何在Android中实现自适应文本大小显示. 想象一下,在布局中,通常显示文本的区域大小是固定的,但是文本长度并不总是固定的.比如列表中的文章标题.界面下方的按钮文本等等. ...
在 sql server 中，查询数据库的大小和数据库中各表的大小
2019独角兽企业重金招聘Python工程师标准>>> 在 sql server 中,查询数据库的大小和数据库中各表的大小其实本来只想找一个方法能查询一下数据库的大小,没 ...
ITK:计算矢量图像中每个像素的大小以生成大小图像
ITK:计算矢量图像中每个像素的大小以生成大小图像内容提要输出结果 C++实现代码内容提要计算矢量图像中每个像素的大小以生成大小图像输出结果 C++实现代码 #include "i ...
发现dba_segments和dba_extents中统计段空间大小居然不一样
最近在测试系统上发现dba_segments和dba_extents中统计段空间大小居然不一样 SQL> select bytes,blocks,extents from dba_segment ...
用递归求最大值PHP,php如何使用递归来计算一个目录中所有文件的大小（代码）...
本篇文章给大家带来的内容是关于php如何使用递归来计算一个目录中所有文件的大小 (代码),有一定的参考价值,有需要的朋友可以参考一下,希望对你有所帮助.sudo find /private/etc - ...

FORCE_CONSTANTS中3阶力常数大小与原子间距的分析脚本

FORCE_CONSTANTS中3阶力常数大小与原子间距的分析脚本相关推荐

最新文章

热门文章