Commit f8317c70 authored by Berk Onat's avatar Berk Onat
Browse files

Flow of the section open/closes corrected with topology/system addition

parent 8cc4a62c
......@@ -583,6 +583,9 @@ def get_updateDictionary(self, defname):
'activeSections' : ['x_amber_mdin_method']
}
# ---------------------------------------------------------------
# Definitions of meta data values for section_sampling_method
# ---------------------------------------------------------------
sampling = {
'ensemble_type' : MetaInfoMap(startpage,
depends=[
......@@ -770,6 +773,9 @@ def get_updateDictionary(self, defname):
)
}
# ------------------------------------------------------------
# Definitions for section_single_configuration_calculation
# ------------------------------------------------------------
singleconfcalc = {
#'atom_forces_type' : MetaInfoMap(startpage,
# depends=[{'assign' : 'Amber Force Field'}],
......@@ -813,6 +819,7 @@ def get_updateDictionary(self, defname):
'stress_tensor_value' : MetaInfoMap(startpage)
}
# section_single_energy_van_der_Waals of section_single_configuration_calculation
singlevdw = {
'energy_van_der_Waals_value' : MetaInfoMap(startpage,
depends=[{'value' : 'VDWAALS'}],
......@@ -822,6 +829,9 @@ def get_updateDictionary(self, defname):
),
}
# ------------------------------------------
# Definitions for section_frame_sequence
# ------------------------------------------
frameseq = {
'frame_sequence_conserved_quantity_frames' : MetaInfoMap(startpage,
depends=[{'store' : 'NSTEP'}],
......@@ -940,49 +950,123 @@ def get_updateDictionary(self, defname):
'previous_sequence_ref' : MetaInfoMap(startpage)
}
topology = {
'molecule_to_molecule_type_map' : MetaInfoMap(startpage),
'number_of_topology_atoms' : MetaInfoMap(startpage,
# ----------------------------------
# Definitions for section_system
# ----------------------------------
# section_system
sec_system = {
#'topology_ref' : MetaInfoMap(startpage),
'atom_velocities' : MetaInfoMap(startpage),
'configuration_raw_gid' : MetaInfoMap(startpage),
'local_rotations' : MetaInfoMap(startpage),
'number_of_atoms' : MetaInfoMap(startpage,
depends=[{'value' : 'NATOM'}],
lookupdict=self.parmDict
),
'number_of_topology_molecules' : MetaInfoMap(startpage),
'topology_force_field_name' : MetaInfoMap(startpage,
value='Amber Force Field',
)
'number_of_sites' : MetaInfoMap(startpage),
'number_of_symmetry_operations' : MetaInfoMap(startpage),
'reduced_symmetry_matrices' : MetaInfoMap(startpage),
'reduced_symmetry_translations' : MetaInfoMap(startpage),
'SC_matrix' : MetaInfoMap(startpage),
'spacegroup_3D_choice' : MetaInfoMap(startpage),
'spacegroup_3D_hall' : MetaInfoMap(startpage),
'spacegroup_3D_international' : MetaInfoMap(startpage),
'spacegroup_3D_number' : MetaInfoMap(startpage),
'spacegroup_3D_origin_shift' : MetaInfoMap(startpage),
'spacegroup_3D_pointgroup' : MetaInfoMap(startpage),
'spacegroup_3D_std_lattice' : MetaInfoMap(startpage),
'spacegroup_3D_std_positions' : MetaInfoMap(startpage),
'spacegroup_3D_std_types' : MetaInfoMap(startpage),
'spacegroup_3D_trasformation_matrix' : MetaInfoMap(startpage),
'spacegroup_3D_wyckoff' : MetaInfoMap(startpage),
'symmorphic' : MetaInfoMap(startpage),
'system_name' : MetaInfoMap(startpage,
subfunction=self.topology_system_name,
#functionbase=self
),
'time_reversal_symmetry' : MetaInfoMap(startpage)
}
atom_list = {
'atom_to_molecule' : MetaInfoMap(startpage,
#subfunction=func_atom_to_mol()
),
# section_configuration_core of section_system
configuration_core = {
'number_of_electrons' : MetaInfoMap(startpage,
value=0,
),
'atom_labels' : MetaInfoMap(startpage,
#subfunction=func_atom_labels()
#subfunction=self.system_atom_labels()
),
'atom_positions' : MetaInfoMap(startpage,
#subfunction=func_atom_positions()
#subfunction=self.system_atom_positions()
),
'configuration_periodic_dimensions' : MetaInfoMap(startpage,
#subfunction=func_pbc()
depends=[
{'test' : [['ntb', '== 0']],
'assign' : np.asarray([False, False, False])},
{'test' : [['ntb', '> 0']],
'assign' : np.asarray([True, True, True])}
],
lookupdict=self.cntrlDict,
),
'embedded_system' : MetaInfoMap(startpage),
'lattice_vectors' : MetaInfoMap(startpage,
#subfunction=func_lat_vec()
#subfunction=self.system_lattice_vectors()
),
'simulation_cell' : MetaInfoMap(startpage,
#subfunction=func_unitcell()
#subfunction=self.system_simulation_cell()
)
}
# section_spacegroup_3D_operation of section_system
spacegroup_op = {
'spacegroup_3D_rotation' : MetaInfoMap(startpage),
'spacegroup_3D_translation' : MetaInfoMap(startpage)
}
# section_system_to_system_refs of section_system
sys_to_sys = {
'system_to_system_kind' : MetaInfoMap(startpage),
'system_to_system_ref' : MetaInfoMap(startpage)
}
# --------------------------------------------------------
# Definitions of meta data values for section_topology
# --------------------------------------------------------
# section_topology of section_run
topology = {
'atom_to_molecule' : MetaInfoMap(startpage,
subfunction=self.topology_atom_to_mol,
#functionbase=self
),
'molecule_to_molecule_type_map' : MetaInfoMap(startpage),
'number_of_topology_atoms' : MetaInfoMap(startpage,
depends=[{'value' : 'NATOM'}],
lookupdict=self.parmDict
),
'number_of_topology_molecules' : MetaInfoMap(startpage,
subfunction=self.topology_num_topo_mol
),
'topology_force_field_name' : MetaInfoMap(startpage,
value='Amber Force Field',
)
}
# section_atom_type of section_topology
atom_type = {
'atom_type_charge' : MetaInfoMap(startpage),
'atom_type_mass' : MetaInfoMap(startpage),
'atom_type_name' : MetaInfoMap(startpage)
}
# section_constraint of section_topology
constraint = {
'constraint_atoms' : MetaInfoMap(startpage),
'constraint_kind' : MetaInfoMap(startpage),
'constraint_parameters' : MetaInfoMap(startpage),
'number_of_atoms_per_constraint' : MetaInfoMap(startpage),
'number_of_constraints' : MetaInfoMap(startpage)
}
# section_interaction of section_topology
interaction = {
'interaction_atoms' : MetaInfoMap(startpage),
'interaction_kind' : MetaInfoMap(startpage),
......@@ -991,13 +1075,27 @@ def get_updateDictionary(self, defname):
'number_of_interactions' : MetaInfoMap(startpage)
}
# -------------------------------------------------------------
# Definitions of meta data values for section_molecule_type
# -------------------------------------------------------------
# section_molecule_type of section_topology
mol_type = {
'molecule_type_name' : MetaInfoMap(startpage),
'number_of_atoms_in_molecule' : MetaInfoMap(startpage),
'settings_atom_in_molecule' : MetaInfoMap(startpage)
}
molecule_interaction = {
# section_molecule_constraint of section_molecule_type
mol_constraint = {
'molecule_constraint_atoms' : MetaInfoMap(startpage),
'molecule_constraint_kind' : MetaInfoMap(startpage),
'molecule_constraint_parameters' : MetaInfoMap(startpage),
'number_of_atoms_per_molecule_constraint' : MetaInfoMap(startpage),
'number_of_molecule_constraints' : MetaInfoMap(startpage)
}
# section_molecule_interaction of section_molecule_type
mol_interaction = {
'molecule_interaction_atoms' : MetaInfoMap(startpage),
'molecule_interaction_kind' : MetaInfoMap(startpage),
'molecule_interaction_parameters' : MetaInfoMap(startpage),
......@@ -1005,6 +1103,7 @@ def get_updateDictionary(self, defname):
'number_of_molecule_interactions' : MetaInfoMap(startpage)
}
# section_atom_in_molecule of section_molecule_type
atom_in_mol = {
'atom_in_molecule_charge' : MetaInfoMap(startpage),
'atom_in_molecule_name' : MetaInfoMap(startpage),
......
......@@ -5,9 +5,10 @@ import setup_paths
import numpy as np
import nomadcore.ActivateLogging
from nomadcore.caching_backend import CachingLevel
from nomadcore.simple_parser import AncillaryParser, mainFunction, ParsingContext
from nomadcore.simple_parser import mainFunction, ParsingContext
from nomadcore.simple_parser import SimpleMatcher as SM
from AMBERDictionary import set_excludeList, set_includeList, get_updateDictionary, getList_MetaStrInDict, getDict_MetaStrInDict
from contextlib import contextmanager
import AMBERCommon as AmberC
import trajectory_reader as TrajRead
import logging
......@@ -24,6 +25,12 @@ LOGGER = logging.getLogger("nomad.AMBERParser")
#PRINTABLE = re.compile(r"\W+")
@contextmanager
def open_section(parser, name):
gid = parser.openSection(name)
yield gid
parser.closeSection(name, gid)
class AMBERParser(AmberC.AMBERParserBase):
"""Context for parsing AMBER main file.
......@@ -84,13 +91,18 @@ class AMBERParser(AmberC.AMBERParserBase):
"""
self.secMethodGIndex = None
self.secSystemGIndex = None
self.secTopologyGIndex = None
self.secSamplingGIndex = None
self.secSingleGIndex = None
self.secVDWGIndex = None
self.secAtomType = None
self.inputMethodIndex = None
self.mainMethodIndex = None
self.mainCalcIndex = None
self.MD = True
self.topologyDict = None
self.topologyTable = None
self.topologyBonds = None
self.topology = None
self.topologyFormat = None
self.topologyFile = None
......@@ -111,6 +123,8 @@ class AMBERParser(AmberC.AMBERParserBase):
self.logFileName = None
self.lastfInLine = None
self.lastfInMatcher = None
self.secOpen = open_section
self.superP = self.parser.backend.superBackend
def startedParsing(self, fInName, parser):
"""Function is called when the parsing starts.
......@@ -172,6 +186,7 @@ class AMBERParser(AmberC.AMBERParserBase):
traj_loaded = self.trajectory.load()
self.atompositions = self.trajectory.iread()
if self.atompositions is not None:
self.topology = self.trajectory.get_topology()
return fileFormat
def initializeFileHandlers(self):
......@@ -186,7 +201,10 @@ class AMBERParser(AmberC.AMBERParserBase):
topoformat = self.topologyFileHandler(fileItem)
# Second check trajectory file
trajformat = self.trajectoryFileHandler(fileItem, topoformat)
print(self.atompositions)
if self.topology:
self.topologyTable, self.topologyBonds = self.topology.to_dataframe()
self.topologyDict = self.topologyTable.to_dict(orient='list')
#print(self.atompositions)
self.topologyFormat = topoformat
self.trajectoryFormat = trajformat
#if trajformat or topoformat:
......@@ -194,6 +212,164 @@ class AMBERParser(AmberC.AMBERParserBase):
#else:
# return False
def topology_to_dictionary(self):
""" This function generates self.topologyDict dictionary
and converts/stores all the topology information
according to the meta info definitions
"""
newDict = {}
return newDict
def topology_num_topo_mol(self, itemdict):
""" Function to generate data for number_of_topology_molecules
"""
if self.topology:
residueList = self.topologyDict["resSeq"]
return False, residueList[len(residueList)-1], itemdict
else:
return False, None, itemdict
def topology_system_name(self, itemdict):
""" Function to generate data for atom_to_molecule
"""
system_name = self.fName.split('.')[-1]
if self.topology:
residueList = ','.join(list(set(self.topologyDict["resName"])))
return False, system_name, itemdict
def topology_atom_to_mol(self, itemdict):
""" Function to generate data for atom_to_molecule
"""
if self.topology:
residueList = self.topologyDict["resSeq"]
atomIndex = np.arange(len(residueList))
atom_to_residue = np.zeros((len(residueList), 2), dtype=int)
atom_to_residue[:,0] = atomIndex
atom_to_residue[:,1] = residueList
return False, atom_to_residue, itemdict
else:
return False, None, itemdict
def topology_atom_type_and_interactions(self):
""" Function to generate data for atom_to_molecule
"""
if self.topology:
atom_type_list=list(set(self.topologyDict["name"]))
atom_type_dict = {}
massesDict = {}
elementDict = {}
radiusDict = {}
for ielm in range(len(atom_type_list)):
elm = atom_type_list[ielm]
atom_type_dict.update({elm : ielm})
for atom in self.topology.atoms:
if elm == atom.name:
massesDict.update({atom.name : atom.element.mass})
elementDict.update({atom.name : atom.element.symbol})
radiusDict.update({atom.name : atom.element.radius})
massesList = list(massesDict.values())
elementList = list(elementDict.values())
radiusList = list(radiusDict.values())
interDict = {}
interTypeDict = {}
for ielm in range(len(atom_type_list)):
elm = atom_type_list[ielm]
bondList = []
typeList = []
bondid = 0
for bond in self.topology.bonds:
molname1, molatom1 = str(bond[0]).split('-')
molname2, molatom2 = str(bond[1]).split('-')
if (elm == str(molatom1) or elm == str(molatom2)):
bondList.append(list(self.topologyBonds[bondid]))
typeList.append(list([
atom_type_dict[str(molatom1)],
atom_type_dict[str(molatom2)]
]))
#str(molatom1)+':'+str(ielm+1),
#str(molatom2)+':'+str(ielm+1)
#]))
interDict.update({elm : bondList})
interTypeDict.update({elm : typeList})
bondid += 1
for elm in range(len(atom_type_list)):
with self.secOpen(self.superP, 'section_atom_type'):
### !!! ----------------------------------------------- !!!
### !!! Need to check the definition of atom type name.
### !!! Which one is better? C1, CA, C:1 or C !!!
### !!! ----------------------------------------------- !!!
# Atom name? Atom element?
self.superP.addValue('atom_type_name', str(atom_type_list[elm]) + ':' + str(elm+1))
# Atomic mass
self.superP.addValue('atom_type_mass', massesList[elm])
# Atomic element/symbol
self.superP.addValue('x_amber_atom_type_element', elementList[elm])
# Atomic van der Waals radius
self.superP.addValue('x_amber_atom_type_radius', radiusList[elm])
# Atomic charge
#self.superP.addValue('atom_type_charge', atom_charge)
pass
for elm in atom_type_list:
with self.secOpen(self.superP, 'section_interaction'):
# atom indexes of bound pairs for a specific atom type
self.superP.addArrayValues('interaction_atoms', np.asarray(interDict[elm]))
# number of bonds for this type
self.superP.addValue('number_of_interactions', len(interDict[elm]))
# number of atoms involved (2 for covalent bonds)
self.superP.addValue('number_of_atoms_per_interaction', len(interDict[elm][0]))
#if bondFunctional:
# self.superP.addValue('interaction_kind', bondFunctional) # functional form of the interaction
# this points to the relative section_atom_type
self.superP.addArrayValues('x_lammps_interaction_atom_to_atom_type_ref',
np.asarray(interTypeDict[elm]))
# interaction parameters for the functional
#self.superP.addValue('interaction_parameters', bondParameters)
# for i in range(len(moleculeTypeInfo)):
#
# with o(p, 'section_molecule_type'):
# # gindex = 0
# p.addValue('molecule_type_name', 'molecule'+'_'+str(moleculeTypeInfo[i][0]))
# p.addValue('number_of_atoms_in_molecule', len(moleculeTypeInfo[i][1]))
#
# p.addArrayValues('atom_in_molecule_to_atom_type_ref', np.asarray([x-1 for x in moleculeTypeInfo[i][1]]))
#
#
# residueList = self.topologyDict["resSeq"]
# atomIndex = np.arange(len(residueList))
# atom_to_residue = np.zeros((len(residueList), 2), dtype=int)
# atom_to_residue[:,0] = atomIndex
# atom_to_residue[:,1] = residueList
# return atom_to_residue
# else:
# return None
# mol_to_mol_bond = []
# mol_to_mol_bond_num = []
# atom_in_mol_bond = []
# atom_in_mol_bond_num = []
# index=0
# print(len([res for res in mytopology.bonds]))
# for res in mytopology.bonds:
# molname1, molatom1 = str(res[0]).split('-')
# molname2, molatom2 = str(res[1]).split('-')
# if molname1 in molname2:
# atom_in_mol_bond.append(res)
# atom_in_mol_bond_num.append(bonds[index])
# else:
# mol_to_mol_bond.append(res)
# mol_to_mol_bond_num.append(bonds[index])
# index += 1
# print(mol_to_mol_bond)
# print(np.array(mol_to_mol_bond_num))
def onClose_section_run(self, backend, gIndex, section):
"""Trigger called when section_run is closed.
......@@ -219,15 +395,6 @@ class AMBERParser(AmberC.AMBERParserBase):
backend.superBackend.addValue(k, v[-1])
# reset all variables
self.initialize_values()
# for k,v in section.simpleValues.items():
# if k == 'x_amber_mdin_file_mdin':
# print('RUN YES',v)
if section['x_amber_mdin_imin'] is not None:
# print('YEEEEEEEEEESSSSSS',section['x_amber_mdin_imin'])
if section['x_amber_mdin_imin'][-1] == '1':
self.MD = True
# else:
# self.MD = False
# if self.MD:
# sampling_method = "molecular_dynamics"
......@@ -249,13 +416,6 @@ class AMBERParser(AmberC.AMBERParserBase):
backend.addArrayValues("frame_sequence_local_frames_ref", np.asarray(self.singleConfCalcs))
backend.closeSection("section_frame_sequence", frameSequenceGIndex)
# def onClose_x_amber_section_MD_detect(self, backend, gIndex, section):
# """Trigger called when fhi_aims_section_MD_detect is closed.
#
# Detect if a MD run was performed.
# """
# self.MD = True
def onClose_x_amber_section_input_output_files(self, backend, gIndex, section):
"""Trigger called when x_amber_section_input_output_files is closed.
......@@ -302,47 +462,8 @@ class AMBERParser(AmberC.AMBERParserBase):
"""Trigger called when section_method is closed.
"""
# input method
if gIndex == self.inputMethodIndex:
self.closingInputMethodSection(backend, gIndex, section)
# for k,v in section.simpleValues.items():
# if (k.startswith('x_amber_mdin_') or
# k.startswith('x_amber_mdout_') or
# k.startswith('x_amber_parm_')):
# backend.superBackend.addValue(k, v[-1])
def closingInputMethodSection(self, backend, gIndex, section):
"""Called when section_method that should contain the main input is closed.
Write the keywords from AMBER output from the parsed mdin, which belong to section_method.
ATTENTION
backend.superBackend is used here instead of only the backend to write the JSON values,
since this allows to bybass the caching setting which was used to collect the values for processing.
However, this also bypasses the checking of validity of the metadata name by the backend.
The scala part will check the validity nevertheless.
"""
# check if control keywords were found or verbatim_writeout is false
#verbatim_writeout = True
counter = 0
exclude_list = set_excludeList(self)
include_list = set_includeList()
#for name in self.metaInfoEnv.infoKinds:
# if name.startswith('x_fhi_aims_mdin_'):
# exclude_list.append(name)
# write settings of aims output from the parsed mdin
for k,v in section.simpleValues.items():
if (k.startswith('x_amber_mdin_') or
# k.startswith('x_amber_mdout_') or
k.startswith('x_amber_parm_')):
# if k in exclude_list and k not in include_list:
# continue
# # default writeout
# else:
backend.superBackend.addValue(k, v[-1])
if k.startswith('x_amber_mdin_imin'):
# print("-----IMIN----",k,v)
if int(v[-1]) == 1:
self.MD = False
#if gIndex == self.inputMethodIndex:
pass
def onOpen_section_sampling_method(self, backend, gIndex, section):
# keep track of the latest sampling description section
......@@ -366,75 +487,78 @@ class AMBERParser(AmberC.AMBERParserBase):
startsection=['section_sampling_method'],
autoopenclose=False)
#backend.closeSection("section_sampling_method", self.secSamplingGIndex)
def onOpen_section_topology(self, backend, gIndex, section):
# keep track of the latest topology description section
if (gIndex is None or gIndex == -1 or gIndex == "-1"):
self.secTopologyGIndex = backend.superBackend.openSection("section_topology")
else:
self.secTopologyGIndex = gIndex
def onClose_section_topology(self, backend, gIndex, section):
"""Trigger called when section_topology is closed.
"""
section_topology_Dict = get_updateDictionary(self, 'topology')
updateDict = {
'startSection' : [
['section_topology']],
#'muteSections' : [['section_method']],
'dictionary' : section_topology_Dict
}
self.metaStorage.update(updateDict)
self.metaStorage.updateBackend(backend.superBackend,
startsection=['section_topology'],
autoopenclose=False)
self.topology_atom_type_and_interactions()
if (gIndex is None or gIndex == -1 or gIndex == "-1"):
backend.superBackend.closeSection("section_topology", self.secTopologyGIndex)
def onOpen_section_system(self, backend, gIndex, section):
# keep track of the latest system description section
self.secSystemGIndex = gIndex
if (gIndex is None or gIndex == -1 or gIndex == "-1"):
self.secSystemGIndex = backend.superBackend.openSection("section_system")
else:
self.secSystemGIndex = gIndex
def onClose_section_system(self, backend, gIndex, section):
"""Trigger called when section_system is closed.
Writes atomic positions, atom labels and lattice vectors.
"""
# check if control keywords were found
#section_sampling_Dict = get_updateDictionary(self, 'sampling')
#updateDict = {
# 'startSection' : [['section_sampling_method']],
# #'muteSections' : [['section_system']],
# 'dictionary' : section_sampling_Dict
# }
#self.secSamplingGIndex = backend.openSection("section_sampling_method")
#self.metaStorage.update(updateDict)
#self.metaStorage.updateBackend(backend,
# startsection=['section_sampling_method'],
# autoopenclose=False)
#backend.closeSection("section_sampling_method", self.secSamplingGIndex)
#self.secMethodGIndex = backend.openSection("section_method")
#backend.closeSection("section_method", self.secMethodGIndex)
# Our special recipe for the confused backend and more
# See our other recipes at the back of this parser: MetaInfoStorage! :)
if (gIndex is None or gIndex == -1 or gIndex == "-1"):
SloppyBackend = backend.superBackend
else:
SloppyBackend = backend
if self.topology:
if (self.secTopologyGIndex is None or
(self.secTopologyGIndex == -1 or
self.secTopologyGIndex == "-1")):
self.onOpen_section_topology(backend, None, None)
self.onClose_section_topology(backend, None, None)
SloppyBackend.addValue("topology_ref", self.secTopologyGIndex)
counter = 0
exclude_list = set_excludeList(self)
include_list = set_includeList()
#for name in self.metaInfoEnv.infoKinds:
# if name.startswith('x_fhi_aims_mdin_'):
# exclude_list.append(name)
# write settings of aims output from the parsed mdin
for k,v in section.simpleValues.items():
if (k.startswith('x_amber_mdin_') or
# k.startswith('x_amber_mdout_') or
k.startswith('x_amber_parm_')):