思路:
1.首先在一个文件夹下不同frame下的、该体系内的不同种类的原子的成键信息
2.利用python安装numpy、numba包
3.将第0 frame下的原子成键信息导入到数组中(本例为二维数组)
4. 不同frame中的成键信息与第0 frame进行对比,判断那些键依然存在,那些键随着时间推移已经断裂、不在成键。
5. 将不同frame与第0frame的对比情况(第4步),写入到文件之中
import numpy as np
from numba import jit
import os
# 文件保存路径
path = "E:/MSD/Hbonds/CAOtest/cao"
path1 = "E:/MSD/Hbonds/CAOtest/repeat"
path2 = "E:/MSD/Hbonds/CAOtest/norepeat"
path_list = os.listdir(path)
path_list.sort(key=lambda x: int(x[0:10])) # 对文件进行排序,确保按照顺序进行读写
# path_list.sort(key=lambda x: int(x.split("_")[1][:-5]))
print(path_list)
# path_list1 = [path_list[i] for i in range(10, len(path_list), 10)]
# 读入第0 frame下的成键情况
datao = np.genfromtxt(os.path.join(path, "0"), dtype=None) # 将文件中数据加载到data数组里
data0 = np.array(datao.tolist()).astype(int) # data0为二维数组
Ca0 = data0[:, 0]
O0 = data0[:, 1]
# print(Ca0)
# print(O0)
for filename in path_list:
data = np.genfromtxt(os.path.join(path, filename), dtype=None) # 将文件中数据加载到data数组里
data1 = np.array(data.tolist()).astype(int)
Ca = data1[:, 0]
O = data1[:, 1]
@jit(nopython=True)
def cao():
repeat = []
# norepeat =[]
for i in range(len(Ca0)):
for j in range(len(Ca)):
if Ca0[i] == Ca[j]:
if O0[i] == O[j]:
repeat.append([Ca0[i], O0[i]])
if Ca0[i] == O[j]:
if O0[i] == Ca[j]:
repeat.append([Ca0[i], O0[i]])
print((len(repeat)))
return repeat
t = cao()
np.savetxt(os.path.join(path1, filename), t) #将仍然存在的键写入到文件中
# print(t1)
t = np.array(t)
# np.savetxt(os.path.join(path2, filename), np.setdiff1d(t,data0),fmt='%s')
# print(np.setdiff1d(t,data0))
norepeat = []
for i in range(len(data0)):
if (data0[i]==t).all(1).any():
continue
else:
norepeat.append(data0[i])
print(len(norepeat))
np.savetxt(os.path.join(path2,filename),norepeat) # 将已经改变的、断裂的键写入到文件中