diff --git a/package/MDAnalysis/core/__init__.py b/package/MDAnalysis/core/__init__.py index d5e239a4b30..1e089e232d1 100644 --- a/package/MDAnalysis/core/__init__.py +++ b/package/MDAnalysis/core/__init__.py @@ -444,5 +444,5 @@ class flagsDocs(object): import groups -import Selection +import selection import Timeseries diff --git a/package/MDAnalysis/topology/TPRParser.py b/package/MDAnalysis/topology/TPRParser.py index 016b635565c..2428895bb6a 100644 --- a/package/MDAnalysis/topology/TPRParser.py +++ b/package/MDAnalysis/topology/TPRParser.py @@ -183,20 +183,13 @@ def parse(self): U.fver_err(V) if th.bTop: - tpr_top = U.do_mtop(data, V, self._u) - structure = { - 'atoms': tpr_top.atoms, - 'bonds': tpr_top.bonds, - 'angles': tpr_top.angles, - 'dihedrals': tpr_top.dihe, - 'impropers': tpr_top.impr - } + tpr_top = U.do_mtop(data, V) else: msg = "{0}: No topology found in tpr file".format(self.filename) logger.critical(msg) raise IOError(msg) - return structure + return tpr_top # THE FOLLOWING CODE IS WORKING FOR TPX VERSION 58, BUT SINCE THESE INFO IS # NOT INTERESTED, SO IT'S NOT COVERED IN ALL VERSIONS. PARSING STOPS HERE. diff --git a/package/MDAnalysis/topology/tpr/obj.py b/package/MDAnalysis/topology/tpr/obj.py index a403641ed88..84fa54d7a45 100644 --- a/package/MDAnalysis/topology/tpr/obj.py +++ b/package/MDAnalysis/topology/tpr/obj.py @@ -34,7 +34,6 @@ "bIr", "bTop", "bX", "bV", "bF", "bBox"]) Box = namedtuple("Box", "size rel v") Mtop = namedtuple("Mtop", "nmoltype moltypes nmolblock") -TPRTopology = namedtuple("TPRTopology", "atoms, bonds, angles, dihe, impr") Params = namedtuple("Params", "atnr ntypes functype reppow fudgeQQ iparams") Atom = namedtuple("Atom", ["m", "q", "mB", "qB", "tp", "typeB", "ptype", "resind", "atomnumber"]) Atoms = namedtuple("Atoms", "atoms nr nres type typeB atomnames resnames") diff --git a/package/MDAnalysis/topology/tpr/utils.py b/package/MDAnalysis/topology/tpr/utils.py index 2ce350ebb8e..799dc8d1ae8 100644 --- a/package/MDAnalysis/topology/tpr/utils.py +++ b/package/MDAnalysis/topology/tpr/utils.py @@ -45,6 +45,24 @@ from . import obj from . import setting as S +import numpy as np +from ..base import TopologyReader, squash_by +from ...core.topologyattrs import ( + Atomids, + Atomnames, + Atomtypes, + Masses, + Charges, + Resids, + Resnames, + Segids, + Bonds, + Angles, + Dihedrals, + Impropers +) +from ...core.topology import Topology + def ndo_int(data, n): """mimic of gmx_fio_ndo_real in gromacs""" return [data.unpack_int() for i in xrange(n)] @@ -146,7 +164,7 @@ def extract_box_info(data, fver): return obj.Box(box, box_rel, box_v) -def do_mtop(data, fver, u): +def do_mtop(data, fver): # mtop: the topology of the whole system symtab = do_symtab(data) do_symstr(data, symtab) # system_name @@ -162,7 +180,19 @@ def do_mtop(data, fver, u): mtop = obj.Mtop(nmoltype, moltypes, nmolblock) - ttop = obj.TPRTopology(*[[] for i in xrange(5)]) + bonds = [] + angles = [] + dihedrals = [] + impropers = [] + + atomids = [] + segids = [] + resids = [] + resnames = [] + atomnames = [] + atomtypes = [] + charges = [] + masses = [] atom_start_ndx = 0 res_start_ndx = 0 @@ -175,20 +205,20 @@ def do_mtop(data, fver, u): for j in xrange(mb.molb_nmol): mt = mtop.moltypes[mb.molb_type] # mt: molecule type for atomkind in mt.atomkinds: - ttop.atoms.append(Atom(atomkind.id + atom_start_ndx, - atomkind.name, - atomkind.type, - atomkind.resname, - atomkind.resid + res_start_ndx, - segid, - atomkind.mass, - atomkind.charge, - universe=u)) + atomids.append(atomkind.id + atom_start_ndx) + segids.append(segid) + resids.append(atomkind.resid + res_start_ndx) + resnames.append(atomkind.resname) + atomnames.append(atomkind.name) + atomtypes.append(atomkind.type) + charges.append(atomkind.charge) + masses.append(atomkind.mass) + # remap_ method returns [blah, blah, ..] or [] - ttop.bonds.extend(mt.remap_bonds(atom_start_ndx)) - ttop.angles.extend(mt.remap_angles(atom_start_ndx)) - ttop.dihe.extend(mt.remap_dihe(atom_start_ndx)) - ttop.impr.extend(mt.remap_impr(atom_start_ndx)) + bonds.extend(mt.remap_bonds(atom_start_ndx)) + angles.extend(mt.remap_angles(atom_start_ndx)) + dihedrals.extend(mt.remap_dihe(atom_start_ndx)) + impropers.extend(mt.remap_impr(atom_start_ndx)) atom_start_ndx += mt.number_of_atoms() res_start_ndx += mt.number_of_residues() @@ -201,7 +231,39 @@ def do_mtop(data, fver, u): # mtop_ffparams_cmap_grid_grid_spacing = 0.1 # mtop_ffparams_cmap_grid_cmapdata = 'NULL' # do_groups(data, symtab) - return ttop + + atomids = Atomids(np.array(atomids, dtype=np.int32)) + atomnames = Atomnames(np.array(atomnames, dtype=object)) + atomtypes = Atomtypes(np.array(atomtypes, dtype=object)) + charges = Charges(np.array(charges, dtype=np.float32)) + masses = Masses(np.array(masses, dtype=np.float32)) + + segids = np.array(segids, dtype=object) + resids = np.array(resids, dtype=np.int32) + resnames = np.array(resnames, dtype=object) + (residx, new_resids, + (new_resnames, perres_segids)) = squash_by(resids, resnames, segids) + residueids = Resids(new_resids) + residuenames = Resnames(new_resnames) + + segidx, perseg_segids = squash_by(perres_segids)[:2] + segids = Segids(perseg_segids) + + top = Topology(len(atomids), len(new_resids), len(perseg_segids), + attrs=[atomids, atomnames, atomtypes, + charges, masses, + residueids, residuenames, + segids], + atom_resindex=residx, + residue_segindex=segidx) + top.add_TopologyAttr(Bonds([bond for bond in bonds if bond])) + top.add_TopologyAttr(Angles([angle for angle in angles if angle])) + top.add_TopologyAttr(Dihedrals([dihedral for dihedral in dihedrals + if dihedral])) + top.add_TopologyAttr(Impropers([improper for improper in impropers + if improper])) + + return top def do_symstr(data, symtab): diff --git a/testsuite/MDAnalysisTests/test_tprparser.py b/testsuite/MDAnalysisTests/test_tprparser.py index 10d6c3dfb03..393e0dee297 100644 --- a/testsuite/MDAnalysisTests/test_tprparser.py +++ b/testsuite/MDAnalysisTests/test_tprparser.py @@ -19,239 +19,125 @@ TPR460, TPR461, TPR502, TPR504, TPR505, TPR510, TPR510_bonded from numpy.testing import TestCase, dec -from test_topology import _TestTopology +from MDAnalysisTests.topology.base import ParserBase import MDAnalysis.topology.TPRParser import functools -@dec.slow -class RefTPR(object): +class TestTPR(ParserBase): """ this test the data/adk_oplsaa.tpr which is of tpx version 58 """ - topology = TPR + filename = TPR parser = MDAnalysis.topology.TPRParser.TPRParser - ref_n_atoms = 47681 - ref_numresidues = 11302 - ref_proteinatoms = 3341 - - -class TestTPR(_TestTopology, RefTPR): - """Testing TPR version 58""" + expected_attrs = ['ids', 'names', 'resids', 'resnames'] + expected_n_atoms = 47681 + expected_n_residues = 11302 + expected_n_segments = 3 # The follow test the same system grompped by different version of gromacs # FORMAT: TPRABC, where numbers ABC indicates the version of gromacs that # generates the corresponding tpr file -class TPRBase(object): +class TPRBase(ParserBase): parser = MDAnalysis.topology.TPRParser.TPRParser - ref_n_atoms = 2263 - ref_numresidues = 230 - ref_proteinatoms = 1962 - - -@dec.slow -class TPR400(TPRBase): - topology = TPR400 - - -class TestTPR400(_TestTopology, TPR400): - """Testing TPR version 58""" - - -@dec.slow -class TPR402(TPRBase): - topology = TPR402 - - -class TestTPR402(_TestTopology, TPR402): - """Testing TPR version 58""" - - -@dec.slow -class TPR403(TPRBase): - topology = TPR403 - - -class TestTPR403(_TestTopology, TPR403): - """Testing TPR version 58""" - - -@dec.slow -class TPR404(TPRBase): - topology = TPR404 - - -class TestTPR404(_TestTopology, TPR404): - """Testing TPR version 58""" - - -@dec.slow -class TPR405(TPRBase): - topology = TPR405 - - -class TestTPR405(_TestTopology, TPR405): - """Testing TPR version 58""" - - -@dec.slow -class TPR406(TPRBase): - topology = TPR406 - - -class TestTPR406(_TestTopology, TPR406): - """Testing TPR version 58""" - - -@dec.slow -class TPR407(TPRBase): - topology = TPR407 - + expected_attrs = ['ids', 'names', 'resids', 'resnames'] + expected_n_atoms = 2263 + expected_n_residues = 230 + expected_n_segments = 2 -class TestTPR407(_TestTopology, TPR407): - """Testing TPR version 58""" +# All these classes should be generated in a loop. Yet, nose test generation +# seems to work only with functions, and not with classes. +class TestTPR400(TPRBase): + filename = TPR400 -@dec.slow -class TPR450(TPRBase): - topology = TPR450 +class TestTPR402(TPRBase): + filename = TPR402 +class TestTPR403(TPRBase): + filename = TPR403 -class TestTPR450(_TestTopology, TPR450): - """Testing TPR version 73""" +class TestTPR404(TPRBase): + filename = TPR404 +class TestTPR405(TPRBase): + filename = TPR405 -@dec.slow -class TPR451(TPRBase): - topology = TPR451 +class TestTPR406(TPRBase): + filename = TPR406 +class TestTPR450(TPRBase): + filename = TPR450 -class TestTPR451(_TestTopology, TPR451): - """Testing TPR version 73""" +class TestTPR451(TPRBase): + filename = TPR451 +class TestTPR452(TPRBase): + filename = TPR452 -@dec.slow -class TPR452(TPRBase): - topology = TPR452 +class TestTPR453(TPRBase): + filename = TPR453 +class TestTPR454(TPRBase): + filename = TPR454 -class TestTPR452(_TestTopology, TPR452): - """Testing TPR version 73""" +class TestTPR455(TPRBase): + filename = TPR455 - -@dec.slow -class TPR453(TPRBase): - topology = TPR453 - - -class TestTPR453(_TestTopology, TPR453): - """Testing TPR version 73""" - - -@dec.slow -class TPR454(TPRBase): - topology = TPR454 - - -class TestTPR454(_TestTopology, TPR454): - """Testing TPR version 73""" - - -@dec.slow -class TPR455(TPRBase): - topology = TPR455 - - -class TestTPR455(_TestTopology, TPR455): - """Testing TPR version 73""" - - -@dec.slow -class TPR455Double(object): +class TPRDouble(ParserBase): parser = MDAnalysis.topology.TPRParser.TPRParser - ref_n_atoms = 21692 - ref_numresidues = 4352 - ref_proteinatoms = 0 # no protein, but DOPC, DPPC, CHOL, SOL - topology = TPR455Double - + expected_attrs = ['ids', 'names', 'resids', 'resnames'] + expected_n_atoms = 21692 + expected_n_residues = 4352 + expected_n_segments = 7 -class TestTPR455Double(_TestTopology, TPR455Double): - """Testing TPR version 73, double precision""" +class TestTPR455Double(TPRDouble): + filename = TPR455Double -class TPR46xBase(object): +class TPR46xBase(ParserBase): parser = MDAnalysis.topology.TPRParser.TPRParser - ref_n_atoms = 44052 - ref_numresidues = 10712 - ref_proteinatoms = 1885 - - -@dec.slow -class TPR460(TPR46xBase): - topology = TPR460 - - -class TestTPR460(_TestTopology, TPR460): - """Testing TPR version 83""" - - -@dec.slow -class TPR461(TPR46xBase): - topology = TPR461 - - -class TestTPR461(_TestTopology, TPR461): - """Testing TPR version 83""" - -@dec.slow -class TPR502(TPRBase): - topology = TPR502 - - -class TestTPR502(_TestTopology, TPR502): - """Testing TPR version 100""" - - -@dec.slow -class TPR504(TPRBase): - topology = TPR504 - + expected_attrs = ['ids', 'names', 'resids', 'resnames'] + expected_n_atoms = 44052 + expected_n_residues = 10712 + expected_n_segments = 8 -class TestTPR504(_TestTopology, TPR504): - """Testing TPR version 100""" -@dec.slow -class TPR505(TPRBase): - topology = TPR505 +class TestTPR460(TPR46xBase): + filename = TPR460 +class TestTPR461(TPR46xBase): + filename = TPR461 -class TestTPR505(_TestTopology, TPR505): - """Testing TPR version 100""" -@dec.slow -class TPR510(TPRBase): - topology = TPR510 +class TestTPR502(TPRBase): + filename = TPR502 +class TestTPR504(TPRBase): + filename = TPR504 -class TestTPR510(_TestTopology, TPR510): - """Testing TPR version 103""" +class TestTPR505(TPRBase): + filename = TPR505 +class TestTPR510(TPRBase): + filename = TPR510 def _test_is_in_topology(name, elements, topology_section, topology_path): """ Test if an interaction appears as expected in the topology """ universe = MDAnalysis.Universe(topology_path) + parser = MDAnalysis.topology.TPRParser.TPRParser(topology_path) + top = parser.parse() for element in elements: - assert element in universe._topology[topology_section], \ + assert element in getattr(top, topology_section).values, \ 'Interaction type "{}" not found'.format(name) def test_all_bonds(): - """Test that all bond types are parsed as expected""" topology = TPR510_bonded bonds = {'BONDS':[(0, 1)], 'G96BONDS':[(1, 2)], 'MORSE':[(2, 3)], 'CUBICBONDS':[(3, 4)], 'CONNBONDS':[(4, 5)], 'HARMONIC':[(5, 6)], diff --git a/testsuite/MDAnalysisTests/topology/base.py b/testsuite/MDAnalysisTests/topology/base.py index ccbfce38663..108c05203a0 100644 --- a/testsuite/MDAnalysisTests/topology/base.py +++ b/testsuite/MDAnalysisTests/topology/base.py @@ -45,9 +45,15 @@ def test_expected_attributes(self): def test_size(self): """Check that the Topology is correctly sized""" - assert_(self.top.n_atoms == self.expected_n_atoms) - assert_(self.top.n_residues == self.expected_n_residues) - assert_(self.top.n_segments == self.expected_n_segments) + assert_(self.top.n_atoms == self.expected_n_atoms, + '{} atoms read, {} expected' + .format(self.top.n_atoms, self.expected_n_atoms)) + assert_(self.top.n_residues == self.expected_n_residues, + '{} residues read, {} expected' + .format(self.top.n_residues, self.expected_n_residues)) + assert_(self.top.n_segments == self.expected_n_segments, + '{} segment read, {} expected' + .format(self.top.n_segments, self.expected_n_segments)) def test_tt_size(self): """Check that the transtable is appropriately sized"""