Packing fixes #385

mattwthompson · 2017-09-07T19:07:41Z

Boxes can be packed with an aspect ratio, i.e. generating a rectangular box with L x L x 3L. This is for the case in which a target density and number of compounds is specified. Previously, a cubic box was generated.
Boxes can now be packed with a mixture of compounds, given the target composition. This is for the case in which a target box size and density are specified. Previously, this case was only supported with one compound. Other cases already supported mixtures.
There is now an argument to create a buffer at the edge of the box (see Add fill_region function #311 (comment)). The default value is 0.2 nm.

Travis will probably fail due to an unrelated bug (fixed in #384) found after this branch was created.

ctk3b · 2017-09-08T17:43:42Z

mbuild/tests/test_packing.py

    def test_fill_box_density_n_compounds(self, h2o):
        filled = mb.fill_box(h2o, density=1000,
                             box=mb.Box([3.1042931, 3.1042931, 3.1042931]))
        assert filled.n_particles == 3000

+    def test_fill_box_compound_ratio(self, h2o, ethane):
+        filled = mb.fill_box(compound=[h2o, ethane], density=800,
+                compound_ratio=[2, 1], box=[2, 2, 2, 4, 4, 4])


Can you also assert here that you indeed have placed h2o/ethane in the correct ratios?

Well that test looks silly. But it would have failed anway if it that assertion existed. Now should be passing!

ctk3b · 2017-09-08T17:44:22Z

mbuild/packing.py

@@ -233,6 +258,9 @@ def solvate(solute, solvent, n_solvent, box, overlap=0.2, seed=12345):
    box_maxs = box.maxs * 10
    overlap *= 10
    center_solute = (box_maxs + box_mins) / 2
+
+    # Apply edge buffer 
+    box_maxs = box_maxs - edge * 10


Do we always want to add this buffer? This is changing prior behavior so need to be a bit careful

I think it should exist as default. Currently it is easy to create a box that has overlaps across the PBC because packmol doesn't. This should really be a default behavior in packmol itself, imo, but they instead give an optional add_box_sides argument. I prefer a buffer approach because if we have a volume inside we we want molecules to be placed, it makes sense to add the buffer inside of the box than increasing the box size. If we want to hash out a list of uses cases of this function I could dig into that, but this handles all the cases that I am working with at the moment and otherwise come to mind.

ctk3b · 2017-09-08T17:44:42Z

mbuild/packing.py


    # In angstroms for packmol.
    box_mins = box.mins * 10
    box_maxs = box.maxs * 10
    overlap *= 10
+
+    # Apply edge buffer 
+    box_maxs = box_maxs - edge * 10


Same question here as below about always applying edge buffer.

ctk3b · 2017-09-08T17:45:11Z

mbuild/packing.py

+                    prototype_mass += r * np.sum([a.mass for a in c.to_parmed().atoms])
+                # Conversion from kg/m^3 / amu * nm^3 to dimensionless units
+                n_prototypes = int(density/prototype_mass*np.prod(box.lengths)*.60224)
+                print(n_prototypes)


Stray debug output (add these checks to the tests)

ctk3b · 2017-09-08T17:46:07Z

mbuild/packing.py

+            else:
+                if compound_ratio is None:
+                    msg = ("Determing `n_compounds` from `density` and `box` "
+                    "for systems with more than one compound type requires"


Could you add indentation here? Makes it visually clear that it's a continuation of the previous line even though it's not a syntax error in this case.

mattwthompson · 2017-09-09T14:57:58Z

Looks like Travis is being stubborn - don't know if a restart would fix it. But the tests are passing for me locally. As we discussed, mbuild/tests/test_coordinate_transform.py::TestCoordinateTransform::test_spin_relative_compound_coordinates is failing but should be fixed with #384 even though it fails here.

ctk3b · 2017-09-10T03:51:52Z

Ok let's let all the tests run to completion though and make sure the #384 error is the only one.

One last request: could you add a few sanity checks on all the various errors that can be raised from incorrect combinations of input values? Can probably be achieved somewhat elegantly with a pytest parameterize that accepts combos and the expected errors (or lack thereof).

summeraz · 2017-09-19T19:02:00Z

mbuild/packing.py

    overlap : float, units nm
+    edge : float, units nm
    density : float, units kg/m^3


Could you add descriptions for each of these arguments in the docstring (I know several of these existed already without a description)?

summeraz · 2017-09-19T19:03:40Z

Looks good. Just a little bit more description in the docstrings I think would be helpful and then this should be ready to merge.

…into pack-fixes

mattwthompson · 2017-09-21T21:16:12Z

Will get to Chrisoph's comments when I get around to figuring out how pytest works

mattwthompson · 2017-10-14T18:26:19Z

Okay, finally got around to this again. As usual, I did weird and bad stuff with git but I think it turned out ok.

Added some tests for bad combinations of n_compounds, density, and box, and also added an error for the case of len(n_compounds) != len(compound), which I did catch previously. I don't know if these should be ValueError or MBuildError

If it looks good, might want to hold off on merging now to see if we can fix #389 #388 in the near week or so.

…into pack-fixes-2

ctk3b · 2017-10-20T02:03:56Z

mbuild/packing.py

@@ -37,7 +39,8 @@
 """


-def fill_box(compound, n_compounds=None, box=None, density=None, overlap=0.2, seed=12345):
+def fill_box(compound, n_compounds=None, box=None, aspect_ratio=None,
+        density=None, overlap=0.2, seed=12345, write_tempfile=False):


Drive by comment: let's append the new args to the end just to avoid breaking any code that called this with positional args.

ctk3b · 2017-10-20T02:05:07Z

mbuild/packing.py

@@ -137,15 +193,20 @@ def fill_box(compound, n_compounds=None, box=None, density=None, overlap=0.2, se
    return filled


-def fill_region(compound, n_compounds, region, overlap=0.2, seed=12345):
+def fill_region(compound, n_compounds, region, overlap=0.2, edge=0.2,
+                seed=12345, write_tempfile=False):


Ditto about arg order

ctk3b · 2017-10-20T02:05:47Z

mbuild/packing.py

@@ -203,16 +277,22 @@ def fill_region(compound, n_compounds, region, overlap=0.2, seed=12345):
    return filled


-def solvate(solute, solvent, n_solvent, box, overlap=0.2, seed=12345):
+def solvate(solute, solvent, n_solvent, box, overlap=0.2, edge=0.2,
+            seed=12345, write_tempfile=False):


One more args order!

mattwthompson · 2017-10-20T02:13:52Z

Okay ... there's been a bit of a bloodbath in the git history but I think it's at least all in one place. (The actual history is messier than what git shows, somehow.)

I have added a flag write_tempfile that optionally copies over the PDB file that packmol writes. This resolves #346

We have also come across a probable solution to #388 #389. Here I will try and explain some of what @summeraz and @tcmoore3 and I worked on today. When packmol fails, there are two results (maybe more).

It can exit with an error, i.e. if the box dimensions are invalid. However, the error messages are passed to STDOUT, not STDERR. So the existing code was not properly called because if err was never satisfied. This is now properly raises the error using the existing _packmol_error helper function. If this happens, packmol writes an empty PDB file. When mdtraj reads an empty PBD file, it doesn't notice that it has no atoms, box, etc. until here, where it's trying to assign disulfide bonds. So most of the cases we encountered this error it was probably the fault of reading in a empty file, not the non-standard residue names. (We should probably do a sanity check when mb.load() is called on an empty file or add this to mdtraj.)
Packmol often exits "without perfect packing," which can be reproduced in at least two different ways
2.1 For systems with a minimum tolerance that is too large (say in excess of 0.5 nm for a simple liquid), it will run the algorithm for a while without really be able to satisfy the constraints, and just end when it hits the iteration limit (default 200).
2.2 This result can also happen when a system too large for the default # iterations to satisfy the constraints. I figure the algorithm scales decently but large systems with a million atoms probably need more than 200 iterations.
In both cases, it writes two files: foo.pdb and foo.pdb_FORCED. Often the first is very far away from the origin, i.e. has values less than -1000. This results is weird stuff in the PDB file, such as a bunch of asterisks where coordinates should be. The latter file represents the system that has been packed but not packed to totally satisfy the constraints. Either way it happens, a warning is read and coordinates from the foo.pdb file is used.

To summarize ... there were problems with the mbuild code but also some pretty silly problems with packmol and mdtraj.

I will look through this again tomorrow to see if the tests are adequate and any cleanup is needed, but this should comprise everything I wanted to get done for a pre-AIChE release.

ctk3b · 2017-10-20T14:18:04Z

Nice work! This looks great.

One thought for a future PR: there's a fair bit of duplicated code flying around this module now which we could probably factor out into a base _run_packmol function of some sort.

mattwthompson · 2017-10-20T14:19:41Z

👍 Andrew and I were just talking about that

ctk3b · 2017-10-24T19:47:01Z

Can you pin the networkX version to the latest 1.X release? We can move to 2.0 after building new conda packages.

summeraz · 2017-10-24T19:54:51Z

The networkX changes have already been made to Foyer, so we can just build that package first and rerun the failed builds for this PR (and the other failing ones) before merging. I'll try to get the new Foyer package built tomorrow.

ctk3b · 2017-11-09T19:17:58Z

Is this good to go since we're punting on the py35 travis issues?

justinGilmer · 2017-11-09T19:21:21Z

@ctk3b Were wanting Matt to merge the current fixes and repush so we can rerun the builds to ensure that we only have the 3.5 issues.

mattwthompson · 2017-11-09T19:33:29Z

In my view, this is either good to go or I can take a stab or your _run_packmol idea now.

ctk3b · 2017-11-09T19:37:44Z

mbuild/packing.py

@@ -124,9 +170,20 @@ def fill_box(compound, n_compounds=None, box=None, density=None, overlap=0.2, se

    proc = Popen(PACKMOL, stdin=PIPE, stdout=PIPE, stderr=PIPE, universal_newlines=True)
    out, err = proc.communicate(input=input_text)
-    if err:
+    #import pdb; pdb.set_trace()


Stray trace

ctk3b · 2017-11-09T19:38:15Z

mbuild/packing.py

    overlap : float
+        Minimum separation between atoms of different molecules.


Could you update the parameters here with the new options?

ctk3b · 2017-11-09T19:38:29Z

mbuild/packing.py

    overlap : float
+        Minimum separation between atoms of different molecules.


Could you update the parameters here with the new options?

ctk3b · 2017-11-09T19:38:52Z

mbuild/packing.py

    density : float, units kg/m^3
+        Target density for the system.


Could you update the parameters here with the new options?

ctk3b · 2017-11-09T19:40:07Z

Ok just did a quick final pass. If you want to take a stab at refactoring it, now's probably the time but it's obviously not crucia. Could you leave a file an issue if you don't?

summeraz · 2017-11-09T20:09:25Z

Also if you could merge in the changes from the latest master that would be great. Just want to make sure all of our tests are passing (except for Python 3.5 on OSX - we still can't seem to figure out that error, but that's a Travis issue).

mattwthompson · 2017-11-11T19:53:33Z

I moved some repeated lines of code into new helper functions: _run_packmol and _check_packmol.

I also changed the writing of temporary files a little bit. Instead of a boolean flag that can write out a pack_temp.pdb file, that argument is now a filename to which it will be saved. This is more natural, adds a tiny amount of functionality, and actually allowed us to remove a couple lines of code. Default is still None. Thanks to @summeraz for the suggestion.

Upon review of these new commits, it should be good to go. It looks like travis is only failing on the known issue with OS X and python3.5. @summeraz and/or @justinGilmer might also want to review this more thoroughly before merging.

summeraz

Looks good to me!

justinGilmer · 2017-11-13T22:49:43Z

Same here!

* Add edge buffer in packing functions * Support aspect ratio for generated boxes * Start working with compound_ratio * Support compound_ratio and add test * Cleanup * Actually build n_compounds and test it * Ensure n_compounds contains ints if list * Add docstrings for packing functions * Test for bad combinations of input parameters * Raise an error if compound and n_compounds are of different len * Add optional flag for writing packmol output * Copy packmol PDB directly; don't write from mbuild compound * First attempt at handling packmol errors * Reorder arguments to match previous versions * Force ruby 2.3 for homebrew (#392) * Test packmol error and warning * Small requested fixes * Add _run_packmol helper function and rework tempfile * Add new arguments to some docstrings * Add _check_packmol helper function * Test temp_file behavior in each packing function * Cleanup some docstrings * Cleanup in fill_region: applying edge buffer & docstrings

mattwthompson added 5 commits September 7, 2017 10:35

Add edge buffer in packing functions

7056e92

Support aspect ratio for generated boxes

3bef5a8

Start working with compound_ratio

a2d0b7d

Support compound_ratio and add test

7b5472f

Cleanup

1f8504e

ctk3b requested changes Sep 8, 2017

View reviewed changes

Actually build n_compounds and test it

ea5ef80

Ensure n_compounds contains ints if list

dff4b1c

summeraz requested changes Sep 19, 2017

View reviewed changes

mattwthompson added 2 commits September 21, 2017 15:16

Add docstrings for packing functions

85af741

Merge branch 'pack-fixes' of https://github.com/mattwthompson/mbuild …

20ae518

…into pack-fixes

summeraz approved these changes Sep 21, 2017

View reviewed changes

mattwthompson added 2 commits October 14, 2017 13:12

Test for bad combinations of input parameters

ac28876

Raise an error if compound and n_compounds are of different len

3af9378

mattwthompson force-pushed the pack-fixes branch from 20ae518 to 3af9378 Compare October 14, 2017 18:25

ctk3b approved these changes Oct 17, 2017

View reviewed changes

mattwthompson added 5 commits October 19, 2017 11:44

Add optional flag for writing packmol output

b806d3c

Copy packmol PDB directly; don't write from mbuild compound

5d7ed28

First attempt at handling packmol errors

87c4739

Merge branch 'pack-fixes-2' of https://github.com/mattwthompson/mbuild …

91a203e

…into pack-fixes-2

Merge branch 'pack-tempfile' into pack-fixes-2

ae715c7

ctk3b reviewed Oct 20, 2017

View reviewed changes

mattwthompson added 2 commits October 19, 2017 21:22

Reorder arguments to match previous versions

19e8940

Force ruby 2.3 for homebrew (mosdef-hub#392)

ac7938a

Test packmol error and warning

08480eb

ctk3b requested changes Nov 9, 2017

View reviewed changes

mattwthompson added 8 commits November 9, 2017 14:12

Merge branch 'master' into pack-fixes

a87c47a

Small requested fixes

c38dc36

Add _run_packmol helper function and rework tempfile

f2747db

Add new arguments to some docstrings

018aafb

Add _check_packmol helper function

94aec64

Test temp_file behavior in each packing function

3c312d5

Cleanup some docstrings

b583953

Cleanup in fill_region: applying edge buffer & docstrings

cc2a436

ctk3b approved these changes Nov 13, 2017

View reviewed changes

summeraz approved these changes Nov 13, 2017

View reviewed changes

summeraz merged commit 7e3b096 into mosdef-hub:master Nov 13, 2017

This was referenced Nov 15, 2017

PDB reader creates bonds based on standard residue templates #389

Closed

Optional flag to write intermediate packmol files #346

Closed

Provide feedback

Saved searches

Use saved searches to filter your results more quickly

Packing fixes #385

Packing fixes #385

mattwthompson commented Sep 7, 2017

ctk3b Sep 8, 2017

mattwthompson Sep 9, 2017

ctk3b Sep 8, 2017

mattwthompson Sep 9, 2017

ctk3b Sep 8, 2017

ctk3b Sep 8, 2017

ctk3b Sep 8, 2017

mattwthompson commented Sep 9, 2017

ctk3b commented Sep 10, 2017

summeraz Sep 19, 2017

summeraz commented Sep 19, 2017

mattwthompson commented Sep 21, 2017

mattwthompson commented Oct 14, 2017

ctk3b Oct 20, 2017

ctk3b Oct 20, 2017

ctk3b Oct 20, 2017

mattwthompson commented Oct 20, 2017

ctk3b commented Oct 20, 2017

mattwthompson commented Oct 20, 2017

ctk3b commented Oct 24, 2017

summeraz commented Oct 24, 2017

ctk3b commented Nov 9, 2017

justinGilmer commented Nov 9, 2017

mattwthompson commented Nov 9, 2017

ctk3b Nov 9, 2017

ctk3b Nov 9, 2017

ctk3b Nov 9, 2017

ctk3b Nov 9, 2017

ctk3b commented Nov 9, 2017

summeraz commented Nov 9, 2017

mattwthompson commented Nov 11, 2017 •

edited

Loading

summeraz left a comment

justinGilmer commented Nov 13, 2017

		overlap : float
		Minimum separation between atoms of different molecules.

Packing fixes #385

Packing fixes #385

Conversation

mattwthompson commented Sep 7, 2017

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

mattwthompson commented Sep 9, 2017

ctk3b commented Sep 10, 2017

Choose a reason for hiding this comment

summeraz commented Sep 19, 2017

mattwthompson commented Sep 21, 2017

mattwthompson commented Oct 14, 2017

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

mattwthompson commented Oct 20, 2017

ctk3b commented Oct 20, 2017

mattwthompson commented Oct 20, 2017

ctk3b commented Oct 24, 2017

summeraz commented Oct 24, 2017

ctk3b commented Nov 9, 2017

justinGilmer commented Nov 9, 2017

mattwthompson commented Nov 9, 2017

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

ctk3b commented Nov 9, 2017

summeraz commented Nov 9, 2017

mattwthompson commented Nov 11, 2017 • edited Loading

summeraz left a comment

Choose a reason for hiding this comment

justinGilmer commented Nov 13, 2017

mattwthompson commented Nov 11, 2017 •

edited

Loading