Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Parmed conversions #334

Merged
merged 26 commits into from
Mar 16, 2020
Merged

Parmed conversions #334

merged 26 commits into from
Mar 16, 2020

Conversation

daico007
Copy link
Member

@daico007 daico007 commented Mar 4, 2020

Something happens with my original Parmed conversions PR so I open a new one here. This PR adds to_parmed method, allowing the conversion from gmso.Topology to parmed.Structure. There is an argument to also convert types and parameters information (default to True).
Basically, I am mapping:

@daico007 daico007 added enhancement New feature or request ready for review labels Mar 4, 2020
@codecov
Copy link

codecov bot commented Mar 4, 2020

Codecov Report

Merging #334 into master will decrease coverage by 0.43%.
The diff coverage is 89.53%.

Impacted file tree graph

@@            Coverage Diff             @@
##           master     #334      +/-   ##
==========================================
- Coverage   95.06%   94.63%   -0.44%     
==========================================
  Files          66       66              
  Lines        3873     4193     +320     
==========================================
+ Hits         3682     3968     +286     
- Misses        191      225      +34
Impacted Files Coverage Δ
gmso/external/__init__.py 100% <100%> (ø) ⬆️
gmso/tests/base_test.py 100% <100%> (ø) ⬆️
gmso/tests/test_convert_parmed.py 90.9% <88.13%> (-9.1%) ⬇️
gmso/external/convert_parmed.py 90.26% <89.32%> (-7.18%) ⬇️
gmso/core/topology.py 92.68% <0%> (+0.48%) ⬆️

Continue to review full report at Codecov.

Legend - Click here to learn more
Δ = absolute <relative> (impact), ø = not affected, ? = missing data
Powered by Codecov. Last update feee88f...8346081. Read the comment docs.

Copy link
Member

@mattwthompson mattwthompson left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Thanks for this! This is no small feat. I will let the linters and CI bots suggest their own changes, but I think the missing pieces here are just more tests, maybe a few along the lines of

  • A system with other elements than C & H

  • A system with multiple residues (i.e. two waters in one structure)

  • Round-trip coversions

  • Checking
    I wonder if each of the _consolidate_core_type functions in from_parmed would be better placed in discrete internal functions, since each block of code produces a single output (the dict). I noticed that you have done this in the opposite direction to_parmed

Comment on lines 36 to 38
# simple check if our pmd.Structure is fully parametrized
is_parametrized = True if isinstance(structure.atoms[0].atom_type,
pmd.AtomType) else False
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Can this this be done more comprehensively? Looking at more atom types and some connection types. I believe pmd can operate with partially-parametrized systems. Unless the point here is to check if anything is parametrized at all

Comment on lines 49 to 50
'sigma': (atom_type.sigma * u.angstrom).in_units(u.nm),
'epsilon': atom_type.epsilon * u.Unit('kcal / mol')},
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Stuff like this is probably worth discussing in the context of #5 #315 - I don't think there's necessarily a downside to converting it into our internal units but it does seem like an unnecessary function call

Comment on lines +67 to +70
assert struc.bond_types == struc_from_top.bond_types
assert struc.angle_types == struc_from_top.angle_types
assert struc.dihedral_types == struc_from_top.dihedral_types
assert struc.rb_torsion_types == struc_from_top.rb_torsion_types
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

How are pmd's equality overrides done? Does this go through and check each attribute of the object against the other? My first guess in a round-trip conversion would be to compare the GMSO objects but this could be easier (and I don't see any reason to believe it tests different behavior, unless there was a reason to do two full loops)

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

For Parmed, two bond types are equal if their parameters are equal. Parmed Structure does not store atom_types, so I cannot compare that one

@daico007
Copy link
Member Author

daico007 commented Mar 4, 2020

I will try to work on breaking _consolidate_core_type into smaller chunks

Break up _consolidate_core_type into bite size chunks,
Add more unit tests and pytest.fixture,
Update documentation.
Copy link
Collaborator

@ahy3nz ahy3nz left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Check your parmed_loop test since that gets skipped.

Watch the the use of default mutable arguments for Angle and Dihedral constructors

import mbuild as mb
import mbuild.recipes
import foyer
from gmso.external.convert_parmed import from_parmed, to_parmed

pentane = mb.recipes.Alkane(5)
pentane.name = "PENT"
box = mb.fill_box(pentane, n_compounds=3, box=[5,5,5])
struc = foyer.forcefields.load_OPLSAA().apply(box, residues="PENT")
print(struc.residues)
top = from_parmed(struc)
print(top.subtops)
new_struc = to_parmed(top)
print(new_struc.residues)

Subtops don't quite match the residues. FF params and bondtypes etc look pristine though

for dihedraltype in structure.dihedral_types:
dihedral_params = {
'k': (dihedraltype.phi_k * u.Unit('kcal / mol')),
'phi_eq': (dihedraltype.phase * u.degree),
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Is there anything about the different interpretations of phi_eq that apply here?

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Do you mean the different names that variable can have?

Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

#31

phi_eq as compared to cis or trans I think. If we take phi_eq as trans, and parmed takes phi_eq as cis, then we're pi radians off

pmd_top_dihedraltypes = {}
for dihedraltype in structure.dihedral_types:
dihedral_params = {
'k': (dihedraltype.phi_k * u.Unit('kcal / mol')),
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Is there a factor of 2 we need to worry about?

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Do you know where I can look up parmed default equation for these Types? I try to look for these in their documentation but could not find anything helpful. Right now, i am basing this conversion (as in if we need a factor of 2) on your earlier version of from_parmed (from your other 2 PRs).

Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

You kind of have to work backwards from how these parameters are stored in OpenMM or gromacs etc, and then figure out parmed digests those parameters to get an idea of what parmed's form is

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I checked the Parmed gro writer and there is no factor of 2 for dihedral's k

Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Okay then I think the gmso equation for periodic torsions matches openmm/gromacs so you're good that no factor of 2 conversion is necessary

@ahy3nz
Copy link
Collaborator

ahy3nz commented Mar 10, 2020

Residue translation from the earlier example looks good now

Copy link
Member

@mattwthompson mattwthompson left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I think these changes are big enough we should try to merge soon, knowing that we'll be tinkering around the edges of ParmEd conversions in perpituity 🙂 (one that comes to mind is impropers, which we just added, but that is appropriate for another PR)

@daico007 daico007 merged commit a452689 into mosdef-hub:master Mar 16, 2020
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
enhancement New feature or request ready for review
Projects
None yet
Development

Successfully merging this pull request may close these issues.

None yet

3 participants