Planned maintenance on Wednesday, 2021-01-20, 17:00-18:00. Expect some interruptions during that time

Commit ae6b4ca4 authored by Carl Poelking's avatar Carl Poelking

Trajectory, configuration, system referencing.

parent 9b09c724
......@@ -2,7 +2,7 @@ import os
import sys
import re
import json
import logging
#import logging
import setup_paths
import numpy as np
......@@ -71,25 +71,32 @@ def parse(output_file_name):
terminal_ctrls = DlPolyControls(osio)
terminal_ctrls.ParseControls(ctrl_file_name)
# PARSE OUTPUT / TOPOLOGY ...
output_file_name = os.path.join(base_dir, 'OUTPUT')
terminal = DlPolyParser(osio)
terminal.ParseOutput(output_file_name)
# PARSE TRAJECTORY ...
# PARSE CONFIG ...
cfg_file_name = os.path.join(base_dir, 'CONFIG')
terminal_trj = DlPolyConfig(osio)
terminal_trj.ParseConfig(cfg_file_name)
terminal_cfg = DlPolyConfig(osio)
terminal_cfg.ParseConfig(cfg_file_name)
# PARSE TRAJECTORY
trj_file_name = os.path.join(base_dir, 'HISTORY')
terminal_trj = DlPolyTrajectory(osio)
terminal_trj.ParseTrajectory(trj_file_name)
# SUMMARIZE KEY-TABLE DEFAULTS ...
terminal.SummarizeKeyDefaults()
terminal.topology.SummarizeKeyDefaults()
terminal_ctrls.SummarizeKeyDefaults()
terminal_cfg.SummarizeKeyDefaults()
terminal_trj.SummarizeKeyDefaults()
# ABBREVIATE ...
ctr = terminal_ctrls
out = terminal
top = terminal.topology
trj = terminal_trj
cfg = terminal_cfg
trj = terminal_trj
ofs = open('parser.keys.log', 'w')
terminals = [ctr, out, top, trj]
terminals = [ctr, out, top, cfg, trj]
for t in terminals:
keys = sorted(t.data.keys())
for key in keys:
......@@ -168,7 +175,9 @@ def parse(output_file_name):
push_array_values(jbe, np.asarray(atoms_to_molidx_atomidx), 'atom_to_molecule')
# SAMPLING-METHOD SECTION
with open_section(jbe, 'section_sampling_method'):
sec_sampling_method_ref = None
with open_section(jbe, 'section_sampling_method') as gid_sec_sampling_method:
sec_sampling_method_ref = gid_sec_sampling_method
# Ensemble
ensemble = push(jbe, out, 'ensemble_type', lambda s: s.As().split()[0].upper())
# Method
......@@ -184,9 +193,55 @@ def parse(output_file_name):
push(jbe, out, 'barostat_target_pressure', lambda s: s.As(float))
push(jbe, out, 'barostat_tau', lambda s: s.As(float))
pass
# TODO Store state variables in frames/system description (temperature, pressure)
# TODO Store interactions
# TODO Add section_method (!= section_sampling_method)
# SYSTEM DESCRIPTION
refs_system_description = []
all_frames = [cfg] + trj.frames # <- Initial config + trajectory
for frame in all_frames:
with open_section(jbe, 'section_system_description') as gid:
refs_system_description.append(gid)
# Configuration core
atom_labels = np.array([ atom['atom_name'].As() for atom in frame.atoms ])
push_array_values(jbe, atom_labels, 'atom_label')
push_array_values(jbe, frame.position_matrix, 'atom_position')
push_array_values(jbe, frame.box_matrix, 'simulation_cell')
push_array_values(jbe, frame.pbc_booleans, 'configuration_periodic_dimensions')
if frame.has_velocities:
push_array_values(jbe, frame.velocity_matrix, 'atom_velocities')
if frame.has_forces:
# TODO Wouldn't it be nicer if forces were added here?
pass
pass
# SINGLE CONFIGURATIONS
refs_single_configuration = []
i_frame = -1
for frame in all_frames:
i_frame += 1
with open_section(jbe, 'section_single_configuration_calculation') as gid:
refs_single_configuration.append(gid)
# Reference system description section
ref_system = refs_system_description[i_frame]
push_value(jbe, ref_system, 'single_configuration_calculation_to_system_description_ref')
# Forces
if frame.has_forces:
push_array_values(jbe, frame.force_matrix, 'atom_forces')
# TODO Set section_method reference
pass
# FRAME-SEQUENCE SECTION
with open_section(jbe, 'section_frame_sequence'):
push_value(jbe, len(all_frames), 'number_of_frames_in_sequence')
# Reference configurations and sampling method
push_value(jbe, sec_sampling_method_ref, 'frame_sequence_to_sampling_ref')
refs_config = np.array(refs_single_configuration)
push_array_values(jbe, refs_config, 'frame_sequence_local_frames')
# TODO Push this to frame_sequence_time
time_values = [ frame['time_value'].As(float) for frame in trj.frames ]
pass
jbe.finishedParsingSession("ParseSuccess", None)
......
......@@ -420,13 +420,135 @@ class DlPolyControls(DlPolyParser):
pass
ifs.close()
return
class DlPolyFrame(DlPolyParser):
def __init__(self, log=None):
super(DlPolyFrame, self).__init__(log)
self.atoms = []
self.has_velocities = False
self.has_forces = False
self.position_matrix = None
self.pbc_booleans = None
self.box_matrix = np.zeros((3,3))
def ParseAtoms(self, ifs, n_atoms, log_level):
# Log level
if log_level > 0:
self.has_velocities = True
if log_level > 1:
self.has_forces = True
# Create atoms
for i in range(n_atoms):
atom_name, atom_id = tuple(ifs.ln().split()[0:2]) # This leaves out weight, charge, rsd
xyz = map(float, ifs.ln().split())
records = [atom_name, atom_id, xyz]
record_labels = ['atom_name', 'atom_id', 'xyz']
if log_level > 0:
vel = map(float, ifs.ln().split())
records.append(vel)
record_labels.append('vel')
if log_level > 1:
force = map(float, ifs.ln().split())
records.append(force)
record_labels.append('force')
new_atom = DlPolyAtom(records, record_labels, self)
self.atoms.append(new_atom)
assert len(self.atoms) == n_atoms
# Position matrix
self.position_matrix = np.zeros((n_atoms, 3))
for i in range(n_atoms):
self.position_matrix[i] = np.array(self.atoms[i]['xyz'].As())
# Velocity matrix
if self.has_velocities:
self.velocity_matrix = np.zeros((n_atoms, 3))
for i in range(n_atoms):
self.velocity_matrix[i] = np.array(self.atoms[i]['vel'].As())
if self.has_forces:
self.force_matrix = np.zeros((n_atoms, 3))
for i in range(n_atoms):
self.force_matrix[i] = np.array(self.atoms[i]['force'].As())
return
def ParseBox(self, ifs, pbc_type):
# PBC booleans
if pbc_type == 6:
self.pbc_booleans = np.array([ True, True, False ])
elif pbc_type > 0:
self.pbc_booleans = np.array([ True, True, True ])
else:
self.pbc_booleans = np.array([ False, False, False ])
# Box
if pbc_type > 0:
a = map(float, ifs.ln().split())
b = map(float, ifs.ln().split())
c = map(float, ifs.ln().split())
self.Set('box_a', a)
self.Set('box_b', b)
self.Set('box_c', c)
# Box matrix
# ATTENTION First index: x,y,z; second index: a,b,c
for i in range(3):
self.box_matrix[i][0] = a[i]
self.box_matrix[i][1] = b[i]
self.box_matrix[i][2] = c[i]
return
def ParseFrame(self, ifs):
ln = ifs.ln()
directives = ln.split()
print directives
assert 'timestep' == directives[0]
# Frame meta-info
self.Set('timestep', directives[1])
self.Set('n_atoms', directives[2])
self.Set('log_level', directives[3])
self.Set('pbc_type', directives[4])
self.Set('dt', directives[5])
self.Set('time_value', self['timestep'].As(int)*self['dt'].As(float))
n_atoms = self['n_atoms'].As(int)
pbc_type = self['pbc_type'].As(int)
# Logging level
log_level = self['log_level'].As(int)
# Box
self.ParseBox(ifs, pbc_type)
# Atoms
self.ParseAtoms(ifs, n_atoms, log_level)
return
class DlPolyTrajectory(DlPolyParser):
def __init__(self, log=None):
super(DlPolyTrajectory, self).__init__(log)
self.logtag = 'trj'
self.frames = [] # List of DlPolyFrame's
def ParseTrajectory(self, trj_file):
if self.log:
self.log << self.log.mg << "Start trajectory ..." << self.log.endl
if not os.path.isfile(trj_file):
trj_file = '__nofile__'
pass
else:
ifs = FileStream(trj_file)
title = ifs.ln().replace('\n','').strip()
# Title
self.Set('title', title)
# Directives: logging, pbc
directives = ifs.ln().split()
self.Set('log_level', directives[0])
self.Set('pbc_type', directives[1])
# Frames
while not ifs.all_read():
new_frame = DlPolyFrame(self.log)
new_frame.ParseFrame(ifs)
self.frames.append(new_frame)
if self.log:
self.log << "Read %d frames from %s" % (len(self.frames), trj_file) << self.log.endl
return
class DlPolyConfig(DlPolyParser):
class DlPolyConfig(DlPolyFrame):
def __init__(self, log=None):
super(DlPolyConfig, self).__init__(log)
DlPolyFrame.__init__(self, log)
#super(DlPolyConfig, self).__init__(log)
self.logtag = 'cfg'
self.atoms = []
return
def ParseConfig(self, trj_file):
if self.log:
......@@ -435,38 +557,20 @@ class DlPolyConfig(DlPolyParser):
# Title
title = ifs.ln().replace('\n','').strip()
self.Set('title', title)
self.Set('time_value', 0.)
# Directives: logging, pbc
directives = ifs.ln().split()
self.Set('log_level', directives[0]) # 0 -> 1 -> 2: coords -> + vel. -> + forces
self.Set('pbc_type', directives[1]) # 0 / ... / 6: no / cubic / orthorhom. / par.-epiped / xy
self.Set('n_atoms', directives[2])
# Box
if self['pbc_type'].As(int) > 0:
a = map(float, ifs.ln().split())
b = map(float, ifs.ln().split())
c = map(float, ifs.ln().split())
self.Set('box_a', a)
self.Set('box_b', b)
self.Set('box_c', c)
# Atom records
n_atoms = self['n_atoms'].As(int)
pbc_type = self['pbc_type'].As(int)
# Logging level
log_level = self['log_level'].As(int)
for i in range(n_atoms):
atom_name, atom_id = tuple(ifs.ln().split())
xyz = map(float, ifs.ln().split())
records = [atom_name, atom_id, xyz]
record_labels = ['atom_name', 'atom_id', 'xyz']
if log_level > 0:
vel = map(float, ifs.ln().split())
records.append(vel)
record_labels.append('vel')
if log_level > 1:
force = map(float, ifs.ln().split())
records.append(force)
record_labels.append('force')
new_atom = DlPolyAtom(records, record_labels, self)
self.atoms.append(new_atom)
assert len(self.atoms) == n_atoms
# Box
self.ParseBox(ifs, pbc_type)
# Atom records
self.ParseAtoms(ifs, n_atoms, log_level)
return
# =======================
......
This diff is collapsed.
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment