Water box

[chrisiacovella]

PR Summary:

This adds a recipe to efficiently initialize a box of 3-site water molecules at ~1000 kg/m^3. This works by reading in a 4nm x 4nm x 4 nm relaxed water configuration (run in gromacs with tip3p and saved as a .gro file ) and performing any duplication/cropping of the configuration as necessary.

Since initializing water is a very common function, it seems to be a reasonable case to have a prebuilt function, as compared to calling packmol every time, since packmol: (1) will be much slower and (2) will not necessarily generate a low energy state for water.

A few features worth noting:

• Users can specify a value for the "edge" just like the packmol fill box, to give some space to prevent overlaps at period boundaries. Note, it is of course possible to remove overlapping molecules in periodic dimensions, but in the current state this would become quite inefficient as system size grows; this feature can be easily added in after I get the cell list PR submitted/merged.
• Users can also specify which 3-site model to use. This code uses the force_overlap routine to orient the model to each molecule in the configuration and updates coordinates and names. Note, because this uses force_overlap, the exact distance and angle specified in the model will be preserved (even if not exactly preserved in the simulation .gro file read in). If a model isn't specified, the tip3p model will be used. Note, even though the .gro file was generated using tip3p, the code still does the force_overlap step such the mbuild configuration that is generated is as ideal as possible. Note, 4 site models are not supported, but could be easily in a future PR after generating a reasonable configuration
• This also includes a simple masking function, where water molecules will only be added if they don't overlap with the mask Compound (or list of mask Compounds). The masking function will use the particle radii (if the element field is defined), to determine this. Users can also specify a dictionary ("radii_dict)" of radii where the Particle name is the key; this will allow masks to be non-atomistic species; this dictionary will take precedence over the element radii (thus allowing for fine tuning of masking). If the element field is not set and dictionary is not passed, a single value will be used (default 0.15 nm, can be set by a user via "radii_overlap" ). Additional scaling can be applied to the radii to provide some more control over the determination of overlaps. This will apply to all radii values, whether radii come from by the ele package, radii_dict, or radii_overlap.

usage is pretty simple:

`
from mbuild.lib.recipes.water_box import Water3SiteBox

new_box = mb.Box([7.0, 3.0, 6.0])
wb = Water3SiteBox(box=new_box)
`

This also addresses issue #1116 where unbounded sites in the compound got added twice after calling condense. The code will now check to ensure an unbonded particle is not part of a compound we have already added. This bug came up when I tried to initialize a tip4p configuration.

PR Checklist

---

• Includes appropriate unit test(s)
• Appropriate docstring(s) are added/updated
• Code is (approximately) PEP8 compliant
• Issue(s) raised/addressed?

[chrisiacovella]

Tests pass locally but fail on the cloud. I think it is related to using packmol; I'll revise the tests.

[daico007]

Said failed test:

self = <mbuild.tests.test_water_box.TestWaterBox object at 0x7f1600a97220>
ethane = <Ethane 8 particles, 7 bonds, non-periodic, id: 139731837261024>

    def test_mask(self, ethane):
        box_temp = mb.Box([2.0, 2.0, 1.0])
    
        ethane_system = mb.fill_box(
            ethane,
            n_compounds=50,
            overlap=0.22,
            edge=0.10,
            sidemax=box_temp.Lz + 1,
            box=box_temp,
            seed=1234,
        )
    
        wb = Water3SiteBox(box=mb.Box([2.0, 2.0, 2.0]), mask=ethane_system)
    
>       assert wb.n_particles == 285
E       assert 288 == 285
E        +  where 288 = <Water3SiteBox 288 particles, 192 bonds, System box: Box: Lx=2.000000, Ly=2.000000, Lz=2.000000, xy=0.000000, xz=0.000000, yz=0.000000, , id: 139731837262032>.n_particles

/home/runner/work/mbuild/mbuild/mbuild/tests/test_water_box.py:67: AssertionError

[daico007]

From the documentation about the mask option, this test will inherently flaky. I think we determined that PACKMOL on different machine (even if they are the same version), will not generate the exact final configurations.

[codecov]

Codecov Report

Patch coverage: 96.22% and project coverage change: +0.25 🎉

Comparison is base (b4d2082) 86.89% compared to head (5248f18) 87.15%.

Additional details and impacted files

@@            Coverage Diff             @@
##             main    #1115      +/-   ##
==========================================
+ Coverage   86.89%   87.15%   +0.25%     
==========================================
  Files          61       62       +1     
  Lines        6319     6422     +103     
==========================================
+ Hits         5491     5597     +106     
+ Misses        828      825       -3     

Impacted Files	Coverage Δ
mbuild/lib/recipes/water_box.py	95.78% <95.78%> (ø)
mbuild/compound.py	97.11% <100.00%> (+0.02%)	⬆️

... and 1 file with indirect coverage changes

☔ View full report in Codecov by Sentry.
📢 Do you have feedback about the report comment? Let us know in this issue.

[daico007]

LGTM! Only minor suggestions regarding docstrings consistency, code styles, etc., but main logic seems great.

[chrisiacovella]

Putting a note here in case we need to reference it later. Mamba started pulling protobuf 4.22.5 (rather than 4.21.12). 4.22.5 won't import properly unless you export PROTOCOL_BUFFERS_PYTHON_IMPLEMENTATION=python in the shell (telling it to use python, not c++). As such, it appears that the problem may actually lie in libprotobuf (the associated c++ library, based on reported issues from prior versions where this has happened). Co pinned 4.21 and all tests appear to pass without needed to do the export statement.

[chrisiacovella]

It looks like bleeding edge tests are failing because it can't find symengine. I'm guessing this package just hasn't been included in the testing yaml file.

[chrisiacovella]

I included symengine in the environment_dev.yml file (confirmed it got installed by mamba), but bleeding edge tests still seems to be failing because they can't find symengine. .

[daico007]

Symengine is technically a C++ package and will need its python API (python-symengine) to work. However, I don't thin we should add symengine to mbuild dev.yml since it technically a GMSO's dependency. Currently, the bleeding test is failing because symengine and python-symengine is not yet included as hard dependency of GMSO (pending the next release). I am more in favor of leaving it out and just have the bleeding test fail for now, since we are planning to release GMSO soon.

[chrisiacovella]

@daico007 I'll revert adding that to the dev.yml and then I think we are ready to merge.

[daico007]

Sounds cool!