Fixing Issue#550

package/CHANGELOG

@@ -53,6 +53,7 @@ Changes
 
 Fixes
 
+  * GROWriter now truncates atom numbers over 99999 (Issue #550)
   * AMBER netcdf writer now correctly uses float32 precision (Issue #518)
   * Fixed a numpy incompatibility in `analysis.leaflet.LeafletFinder`.
     (Issue #533)

package/MDAnalysis/analysis/psa.py

@@ -1075,7 +1075,7 @@ class PSAnalysis(object):
     """
     def __init__(self, universes, reference=None, ref_select='name CA',
                  ref_frame=0, path_select=None, labels=None,
-                 targetdir=None):
+                 targetdir=os.path.curdir):
         """Setting up Path Similarity Analysis.
 
         The mutual similarity between all unique pairs of trajectories

package/MDAnalysis/coordinates/GRO.py

@@ -286,11 +286,12 @@ def write(self, selection, frame=None):
             output_gro.write(self.fmt['n_atoms'].format(len(atoms)))
             # Atom descriptions and coords
             for atom_index, atom in enumerate(atoms):
+                truncated_atom_index = int(str(atom_index + 1)[-5:])
                 if has_velocities:
                     output_gro.write(self.fmt['xyz_v'].format(
                         resid=atom.resid,
                         resname=atom.resname,
-                        index=atom_index+1,
+                        index=truncated_atom_index,
                         name=atom.name,
                         pos=coordinates[atom_index],
                         vel=velocities[atom_index],
@@ -299,7 +300,7 @@ def write(self, selection, frame=None):
                     output_gro.write(self.fmt['xyz'].format(
                         resid=atom.resid,
                         resname=atom.resname,
-                        index=atom_index+1,
+                        index=truncated_atom_index,
                         name=atom.name,
                         pos=coordinates[atom_index]
                     ))

testsuite/MDAnalysisTests/analysis/test_psa.py

@@ -27,16 +27,17 @@
 from MDAnalysisTests.datafiles import PSF, DCD, DCD2
 
 
-
 class TestPSAnalysis(TestCase):
     def setUp(self):
+        self.tmpdir = tempfile.mkdtemp()
         self.iu1 = np.triu_indices(3, k=1)
         self.universe1 = MDAnalysis.Universe(PSF, DCD)
         self.universe2 = MDAnalysis.Universe(PSF, DCD2)
         self.universe_rev = MDAnalysis.Universe(PSF, DCD)
         self.universes = [self.universe1, self.universe2, self.universe_rev]
         self.psa = MDAnalysis.analysis.psa.PSAnalysis(self.universes,           \
-                                               path_select='name CA')
+                                               path_select='name CA',           \
+                                               targetdir=self.tmpdir)
         self.psa.generate_paths(align=True)
         self.psa.paths[-1] = self.psa.paths[-1][::-1,:,:] # reverse third path
         self._run()
@@ -54,7 +55,8 @@ def tearDown(self):
         del self.universe2
         del self.universe_rev
         del self.psa
-        shutil.rmtree('psadata') # remove psa data directory
+        if self.tmpdir:
+            shutil.rmtree(self.tmpdir)
 
     def test_hausdorff_bound(self):
         err_msg = "Some Frechet distances are smaller than corresponding "      \

testsuite/MDAnalysisTests/coordinates/test_gro.py

@@ -1,14 +1,15 @@
 import MDAnalysis as mda
 import numpy as np
 import os
+import bz2
 
 from nose.plugins.attrib import attr
 from numpy.testing import (assert_equal, assert_almost_equal, dec,
                            assert_array_almost_equal, assert_raises)
 from unittest import TestCase
 import tempdir
 
-from MDAnalysisTests.datafiles import (GRO, GRO_velocity)
+from MDAnalysisTests.datafiles import (GRO, GRO_velocity, GRO_large)
 from MDAnalysisTests.coordinates.reference import RefAdK
 from MDAnalysisTests.coordinates.base import BaseTimestepTest
 
@@ -235,6 +236,31 @@ def test_check_coordinate_limits_max_noconversion(self):
         del u
 
 
+class TestGROWriterLarge(TestCase, tempdir.TempDir):
+    def setUp(self):
+        self.tmpdir = tempdir.TempDir()
+        self.large_universe = mda.Universe(GRO_large)
+
+    def tearDown(self):
+        del self.tmpdir
+        del self.large_universe
+
+    @dec.slow
+    @attr('issue')
+    def test_writer_large(self):
+        """Test that atom numbers are truncated for large
+        GRO files (Issue 550)."""
+        outfile = self.tmpdir.name + '/outfile1.gro'
+        self.large_universe.atoms.write(outfile)
+        with open(outfile, 'r') as mda_output:
+            with bz2.BZ2File(GRO_large, 'r') as expected_output:
+                produced_lines = mda_output.readlines()[1:]
+                expected_lines = expected_output.readlines()[1:]
+                assert_equal(produced_lines,
+                             expected_lines,
+                             err_msg="Writing GRO file with > 100 000 "
+                                 "coords does not truncate properly.")
+
 class TestGROWriterVels(object):
     def setUp(self):
         self.tmpdir = tempdir.TempDir()

testsuite/MDAnalysisTests/datafiles.py

@@ -47,6 +47,7 @@
     "NUCL",  # nucleic acid (PDB)
     "INC_PDB",  # incomplete PDB file (Issue #396)
     "PDB", "GRO", "XTC", "TRR", "TPR", "GRO_velocity",  # Gromacs (AdK)
+    "GRO_large", #atom number truncation at > 100,000 particles, Issue 550
     "PDB_xvf", "TPR_xvf", "TRR_xvf",  # Gromacs coords/veloc/forces (cobrotoxin, OPLS-AA, Gromacs 4.5.5 tpr)
     "PDB_xlserial",
     "TPR400", "TPR402", "TPR403", "TPR404", "TPR405", "TPR406", "TPR407",
@@ -130,6 +131,7 @@
 
 GRO = resource_filename(__name__, 'data/adk_oplsaa.gro')
 GRO_velocity = resource_filename(__name__, 'data/sample_velocity_file.gro')
+GRO_large = resource_filename(__name__, 'data/bigbox.gro.bz2')
 PDB = resource_filename(__name__, 'data/adk_oplsaa.pdb')
 XTC = resource_filename(__name__, 'data/adk_oplsaa.xtc')
 TRR = resource_filename(__name__, 'data/adk_oplsaa.trr')