Skip to content

Commit

Permalink
fep.Gsolv automatically converts EDR to XVG
Browse files Browse the repository at this point in the history
- completely resolves #82
- test added
- CHANGES updated
  • Loading branch information
orbeckst committed May 2, 2017
1 parent 6b90ac9 commit c3b7fbd
Show file tree
Hide file tree
Showing 3 changed files with 44 additions and 16 deletions.
3 changes: 3 additions & 0 deletions CHANGES
Original file line number Diff line number Diff line change
Expand Up @@ -20,6 +20,9 @@ orbeckst, iorga, ianmkenney, rhheilma
* FEP failed with IndexError when lambdas where presented as a list
(#78)
* EDR files are now output by default instead of XVG (#75)
* fep.Gsolv.analyze() automatically converts EDR to XVG.bz2 if
necessary and can process the XVG/MBAR format of Gromacs 5.x
for TI (issue #82)
* deprecated mdpow-ghyd script (issue #14)


Expand Down
10 changes: 9 additions & 1 deletion mdpow/fep.py
Original file line number Diff line number Diff line change
Expand Up @@ -896,7 +896,15 @@ def analyze(self, force=False, stride=None, autosave=True, ncorrel=25000):
stride = stride or self.stride

if force or not self.has_dVdl():
self.collect(stride=stride, autosave=False)
try:
self.collect(stride=stride, autosave=False)
except IOError as err:
if err.errno == errno.ENOENT:
self.convert_edr()
self.collect(stride=stride, autosave=False)
else:
logger.exception()
raise
else:
logger.info("Analyzing stored data.")

Expand Down
47 changes: 32 additions & 15 deletions mdpow/tests/test_analysis.py
Original file line number Diff line number Diff line change
Expand Up @@ -52,32 +52,24 @@ def fix_manifest(topdir):

# session scope if read-only use

@pytest.fixture(scope="session")
@pytest.fixture(scope="function")
def fep_benzene_directory(tmpdir_factory):
topdir = tmpdir_factory.mktemp('analysis')
m = pybol.Manifest(fix_manifest(topdir).strpath)
m.assemble('FEP', topdir.strpath)
return topdir.join("benzene")

class TestAnalyze(object):
def test_TI(self, fep_benzene_directory):
gsolv = fep_benzene_directory.join("FEP", "water", "Gsolv.fep")
def get_Gsolv(self, pth):
gsolv = pth.join("FEP", "water", "Gsolv.fep")
G = pickle.load(gsolv.open())
# patch paths
G.basedir = fep_benzene_directory.strpath
G.basedir = pth.strpath
G.filename = gsolv.strpath
return G

# convert EDR to XVG.bz2; if the fixture is session scoped then other
# workers will pick up these files. Make sure that only one runs convert because
# there is no file locking, if in doubt, make fep_benzene_directory locally scoped
G.convert_edr()

try:
G.analyze(force=True, autosave=False)
except IOError as err:
raise AssertionError("Failed to convert edr to xvg: {0}: {1}".format(
err.strerror, err.filename))

@staticmethod
def assert_DeltaA(G):
DeltaA = G.results.DeltaA
assert_array_almost_equal(DeltaA.Gibbs.astuple(),
(-3.7217472974883794, 2.3144288928034911),
Expand All @@ -89,3 +81,28 @@ def test_TI(self, fep_benzene_directory):
(-4.6128782195215781, 2.1942144688960972),
decimal=6)


def test_convert_edr(self, fep_benzene_directory):
G = self.get_Gsolv(fep_benzene_directory)
try:
G.analyze(force=True, autosave=False)
except IOError as err:
raise AssertionError("Failed to auto-convert edr to xvg: {0}: {1}".format(
err.strerror, err.filename))
self.assert_DeltaA(G)


def test_TI(self, fep_benzene_directory):
G = self.get_Gsolv(fep_benzene_directory)
# ensure conversion EDR to XVG.bz2; if the fixture is session scoped
# then other workers will pick up these files. Make sure that only one
# runs convert because there is no file locking, if in doubt, make
# fep_benzene_directory locally scoped
G.convert_edr()
try:
G.analyze(force=True, autosave=False)
except IOError as err:
raise AssertionError("Failed to convert edr to xvg: {0}: {1}".format(
err.strerror, err.filename))
self.assert_DeltaA(G)

0 comments on commit c3b7fbd

Please sign in to comment.