test for NCRST

Author: jeremyleung521

testsuite/MDAnalysisTests/coordinates/test_netcdf.py

@@ -41,6 +41,8 @@
     DLP_CONFIG,
     CPPTRAJ_TRAJ_TOP,
     CPPTRAJ_TRAJ,
+    PRMNCRST,
+    TRJNCRST,
 )
 from MDAnalysisTests.coordinates.test_trj import _TRJReaderTest
 from MDAnalysisTests.coordinates.reference import RefVGV, RefTZ2
@@ -389,6 +391,72 @@ def test_warn_user_no_time_information(self, u):
             u2 = mda.Universe(CPPTRAJ_TRAJ_TOP, CPPTRAJ_TRAJ)
 
 
+class TestNCDFReader5(object):
+    """NCRST Restart File with positions and forces, exported by CPPTRAJ.
+
+    Contributed by Jeremy M. G. Leung
+    """
+
+    prec = 6
+
+    @pytest.fixture(scope="class")
+    def u(self):
+        return mda.Universe(PRMNCRST, TRJNCRST)
+
+    def test_positions(self, u):
+        """Check positions on first frame"""
+        u.trajectory[0]
+        ref_1 = np.array(
+            [
+                [-1.1455358, -2.0177484, -0.55771565],
+                [-0.19042611, -2.2511053, -1.0282656],
+                [ 0.53238064, -1.5778863, -0.56737846],
+            ],
+            dtype=np.float64,
+        )
+        assert_almost_equal(ref_1, u.atoms.positions[:3], self.prec)
+
+    def test_velocities(self, u):
+        """Check forces on first frame"""
+        u.trajectory[0]
+        ref_1 = np.array(
+            [
+                [11.86471367, 31.22108269, -4.03538418],
+                [7.36676359, -4.68035316, 1.78124952],
+                [12.86675262, 1.39324546, -14.97190762],
+            ],
+            dtype=np.float64,
+        )
+        assert_almost_equal(ref_1, u.atoms.velocities[:3], self.prec)
+
+    def test_forces(self, u):
+        """Check forces on first frame"""
+        u.trajectory[0]
+        ref_1 = np.array(
+            [
+                [-9.7262249, -0.37627634, -24.79876137],
+                [-32.43384552, 37.47783279, -12.45422745],
+                [24.10939026, -19.80618095, -17.09523582],
+            ],
+            dtype=np.float64,
+        )
+        assert_almost_equal(ref_1, u.atoms.forces[:3], self.prec)
+
+    def test_time(self, u):
+        """Check time on first frame"""
+        ref = 5.0
+        assert_almost_equal(ref, u.trajectory[0].time, self.prec)
+
+    def test_dt(self, u):
+        """Default 1.0 fs"""
+        ref = 1.0
+        assert_almost_equal(ref, u.trajectory.dt, self.prec)
+        assert_almost_equal(ref, u.trajectory.ts.dt, self.prec)
+
+    def test_box(self, u):
+        assert u.trajectory[0].dimensions is None
+
+
 class _NCDFGenerator(object):
     """A class for generating abitrary ncdf files and exhaustively test
     edge cases which might not be found in the wild"""

testsuite/MDAnalysisTests/data/Amber/ace_mbondi3.ncrst

@@ -191,7 +191,8 @@
     "CPPTRAJ_TRAJ_TOP",
     "CPPTRAJ_TRAJ",  # Amber ncdf extracted from CPPTRAJ without time variable
     "PRMcs",  # Amber (format, Issue 1331)
-    "PRMNCRST",  # Amber ncrst with positions/forces/velocities
+    "PRMNCRST",  # Amber topology for ncrst with positions/forces/velocities
+    "TRJNCRST",  # Amber ncrst with positions/forces/velocties
     "PRM_NCBOX",
     "TRJ_NCBOX",  # Amber parm7 + nc w/ pos/forces/vels/box
     "PRMNEGATIVE",  # Amber negative ATOMIC_NUMBER (Issue 2306)
@@ -662,6 +663,7 @@
 PRMcs = (_data_ref / "Amber/chitosan.prmtop").as_posix()
 
 PRMNCRST = (_data_ref / "Amber/ace_mbondi3.parm7").as_posix()
+TRJNCRST = (_data_ref / "Amber/ace_mbondi3.ncrst").as_posix()
 
 PRM_NCBOX = (_data_ref / "Amber/ace_tip3p.parm7").as_posix()
 TRJ_NCBOX = (_data_ref / "Amber/ace_tip3p.nc").as_posix()