Skip to content

Commit

Permalink
backported fix for DCD and XTC/TRR pickling
Browse files Browse the repository at this point in the history
- fix issue #2878
- add tests to the base tests (but xfail for any format that does not
  implement pickling yet)
  • Loading branch information
yuxuanzhuang authored and orbeckst committed Aug 13, 2020
1 parent 7313cf1 commit 9e89878
Show file tree
Hide file tree
Showing 4 changed files with 45 additions and 4 deletions.
3 changes: 2 additions & 1 deletion package/CHANGELOG
Original file line number Diff line number Diff line change
Expand Up @@ -13,7 +13,7 @@ The rules for this file:
* release numbers follow "Semantic Versioning" http://semver.org

------------------------------------------------------------------------------
??/??/20 orbeckst, VOD555, lilyminium
??/??/20 orbeckst, VOD555, lilyminium, yuxuanzhuang

* 1.0.1

Expand All @@ -24,6 +24,7 @@ Fixes
density=True; the keyword was available since 0.19.0 but with incorrect
semantics and not documented and did not produce correct results (Issue
#2811, PR #2812)
* pickle correct frames of DCD, TRR, and XTC (#2878)


06/09/20 richardjgowers, kain88-de, lilyminium, p-j-smith, bdice, joaomcteixeira,
Expand Down
3 changes: 2 additions & 1 deletion package/MDAnalysis/lib/formats/libdcd.pyx
Original file line number Diff line number Diff line change
Expand Up @@ -265,7 +265,8 @@ cdef class DCDFile:
return

current_frame = state[1]
self.seek(current_frame)
self.seek(current_frame - 1)
self.current_frame = current_frame


def tell(self):
Expand Down
3 changes: 2 additions & 1 deletion package/MDAnalysis/lib/formats/libmdaxdr.pyx
Original file line number Diff line number Diff line change
Expand Up @@ -306,7 +306,8 @@ cdef class _XDRFile:

# where was I
current_frame = state[1]
self.seek(current_frame)
self.seek(current_frame - 1)
self.current_frame = current_frame

def seek(self, frame):
"""Seek to Frame.
Expand Down
40 changes: 39 additions & 1 deletion testsuite/MDAnalysisTests/coordinates/base.py
Original file line number Diff line number Diff line change
Expand Up @@ -25,6 +25,7 @@
import numpy as np
import pytest
from six.moves import zip, range
from six.moves import cPickle as pickle
from unittest import TestCase
from numpy.testing import (assert_equal, assert_almost_equal,
assert_array_almost_equal, assert_allclose)
Expand Down Expand Up @@ -419,12 +420,22 @@ def test_transformations_copy(self,ref,transformed):
ideal_coords = ref.iter_ts(i).positions + v1 + v2
assert_array_almost_equal(ts.positions, ideal_coords, decimal = ref.prec)


def test_add_another_transformations_raises_ValueError(self, transformed):
# After defining the transformations, the workflow cannot be changed
with pytest.raises(ValueError):
transformed.add_transformations(translate([2,2,2]))

def test_pickle_reader(self, reader):
try:
reader_p = pickle.loads(pickle.dumps(reader))
except (TypeError, ValueError):
pytest.xfail("Pickling not yet implemented for format {}".format(
reader.format))
assert_equal(len(reader), len(reader_p))
assert_equal(reader.ts, reader_p.ts,
"Timestep is changed after pickling")


class MultiframeReaderTest(BaseReaderTest):
def test_last_frame(self, ref, reader):
ts = reader[-1]
Expand Down Expand Up @@ -490,6 +501,33 @@ def test_iter_as_aux_lowf(self, ref, reader):
ref.iter_ts(ref.aux_lowf_frames_with_steps[i]),
decimal=ref.prec)

# To make sure we not only save the current timestep information,
# but also maintain its relative position.
def test_pickle_next_ts_reader(self, reader):
try:
reader_p = pickle.loads(pickle.dumps(reader))
except (TypeError, ValueError):
pytest.xfail("Pickling not yet implemented for format {}".format(
reader.format))
assert_equal(next(reader), next(reader_p),
"Next timestep is changed after pickling")

# To make sure pickle works for last frame.
def test_pickle_last_ts_reader(self, reader):
# move current ts to last frame.
reader[-1]
try:
reader_p = pickle.loads(pickle.dumps(reader))
except (TypeError, ValueError):
pytest.xfail("Pickling not yet implemented for format {}".format(
reader.format))
assert_equal(len(reader), len(reader_p),
"Last timestep is changed after pickling")
assert_equal(reader.ts, reader_p.ts,
"Last timestep is changed after pickling")




class BaseWriterTest(object):
@staticmethod
Expand Down

0 comments on commit 9e89878

Please sign in to comment.