In [9]:
import pandas as pd
import exa, exatomic
from exatomic import qe
from exatomic.base import resource, list_resources # Helper to find example files available in exatomic
In [10]:
# list_resources shows all available example files; here are some QE related ones
[i for i in list_resources() if "qe" in i]
Out[10]:
In [33]:
inppath = resource("qe.inp")
pospath = resource("qe.pos")
velpath = resource("qe.vel")
forpath = resource("qe.for")
evppath = resource("qe.evp")
In [12]:
exa.Editor(velpath).head() # We can use the Editor to view the actual file contents
In [22]:
exa.Editor(inppath)[34:38] # Likewise for the input file...
# Note that the atom order of the pos, vel, etc files is the
# order of the ATOMIC_SPECIES block
Out[22]:
In [ ]:
# This is one way to access the docs of a specific function
# in the Jupyter notebook or IPython terminal
qe.parse_symbols_from_input?
In [24]:
# Here we get and order the symbol list
xyz = qe.parse_symbols_from_input(inppath)
xyz # pandas DataFrame
Out[24]:
In [27]:
# Turns out the symbols are in the right order but if we needed to reorder them...
def mapper(sym):
return {'H': 0, 'O': 1}[sym]
xyz['order'] = xyz['symbol'].map(mapper)
symbols = xyz.sort_values('order')['symbol'].tolist() # symbols = xyz['symbol'].tolist()
symbols
Out[27]:
In [28]:
pos = qe.parse_xyz(pospath, symbols)
vel = qe.parse_xyz(velpath, symbols)
vel.columns = ("vx", "vy", "vz", "frame", "symbol")
force = qe.parse_xyz(forpath, symbols)
force.columns = ("fx", "fy", "fz", "frame", "symbol")
In [29]:
pos.head()
Out[29]:
In [30]:
vel.head()
Out[30]:
In [31]:
force.head()
Out[31]:
In [32]:
atom = pd.concat((pos, vel[["vx", "vy", "vz"]], force[["fx", "fy", "fz"]]), axis=1)
In [36]:
# Some extra steps to set the column names
evp = qe.parse_evp(evppath, symbols)
columns = evp.iloc[0, :].tolist()[1:] + ["None"] # Add a dummy column
evp = evp.drop("#", axis=0) # Delete first row
evp.columns = columns
evp['atom_count'] = len(symbols)
evp = evp.dropna(how='all', axis=1) # Remove dummy column
evp.head()
Out[36]:
In [37]:
# Make a universe
uni = exatomic.Universe(atom=atom, frame=evp)
In [38]:
len(uni)
Out[38]:
In [39]:
uni.atom.dtypes
Out[39]:
In [40]:
%time uni.compute_atom_two() # Compute interatomic distances
In [41]:
uni.atom.shape
Out[41]:
In [42]:
uni.atom_two.shape
Out[42]:
In [43]:
uni.molecule.shape # Automatically computes molecules from atom-atom distances
Out[43]:
In [ ]:
exatomic.UniverseWidget(uni) # Visualize the universe!
In [ ]:
In [ ]:
In [ ]:
In [ ]: