-
Notifications
You must be signed in to change notification settings - Fork 0
/
plot_xyz.py
executable file
·81 lines (68 loc) · 1.94 KB
/
plot_xyz.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
#!/bin/env python3
import numpy as np
from matplotlib import rc
import matplotlib.pyplot as plt
import matplotlib.ticker as mtick
import sys
from plot import *
show = True # show plots?
dfpath = "vi.xyz"
tn = 10000 # number of time steps
dt = 1 # read frequency
if len(sys.argv) > 1:
dfpath = sys.argv[1]
if len(sys.argv) > 2:
tn = int(sys.argv[2])
if len(sys.argv) > 3:
dt = int(sys.argv[3])
if len(sys.argv) > 2:
show = bool(int(sys.argv[4]))
def parsexyz(path, tn, dt):
df = open(path, "r")
start = False
atoms = 0 # number of atoms
time = [] # time
pos = [] # xyz position
i = 0
j = 0
tj = 0
t = 0
posj = []
for l in df: # parse dump file
if "Time=" in l:
if start:
j += 1
tj = j * dt
time.append(t)
pos.append(np.array(posj))
posj = []
start = False
continue
if i == tj:
ts = l.split(",")[0].split("=")[1]
t = float(ts)
start = True
i += 1
continue
i += 1
if len(l.split()) == 1:
atoms = int(l)
elif start: # append data to array
words = l.split()
try:
w0 = int(words[0])
w1 = float(words[1])
w2 = float(words[2])
w3 = float(words[3])
w4 = float(words[4])
except Exception as e:
print(e)
posj.append([w2, w3, w4])
return np.array(time), pos
times, r = parsexyz(dfpath, tn, dt)
msd = []
for i in range(len(times)):
ri = r[i]
msd.append(np.mean(np.square(ri[:, 0])+np.square(ri[:, 1])+np.square(ri[:, 2])))
msd = np.array(msd)
linescatter([np.array([times, msd]).T], titles=["Mean squared displacement", "$t$ (fs)", "MSD (Å$^2$)"], fpath="msd.pdf")