Skip to content

Fix parsing of GLYCAM topology files#407

Open
lohedges wants to merge 8 commits intodevelfrom
fix_glycam
Open

Fix parsing of GLYCAM topology files#407
lohedges wants to merge 8 commits intodevelfrom
fix_glycam

Conversation

@lohedges
Copy link
Contributor

This PR fixes several bugs affecting systems parameterised with the GLYCAM force field (Glycam-06), which uses full 1-4 non-bonded interactions (SCEE=1.0, SCNB=1.0) rather than the standard AMBER scaling (SCEE=1.2, SCNB=2.0). GLYCAM topologies from CHARMM-GUI represent these interactions using [pairs] funct=2 with explicit parameters, rather than funct=1 which applies the global fudgeLJ/fudgeQQ scaling from [defaults].

Bugs fixed:

GROMACS read (grotop.cpp): [pairs] entries were listed in missed_tags and completely ignored on read. All 1-4 pairs were assigned the global fudgeQQ/fudgeLJ scaling factors, so GLYCAM's full 1-4 interactions were silently replaced with standard AMBER scaling. Fixed by parsing funct=2 entries and storing them as explicit pair overrides applied after the gen-pairs CLJNBPairs is constructed.

GROMACS write (grotop.cpp): CLJScaleFactor(1.0, 1.0) pairs fell through the write_pairs logic without being written. Since nrexcl=3 excludes all 1-4 pairs by default, omitting them from [pairs] gives zero 1-4 interaction. Fixed by adding a funct=2 write branch that computes combined LJ parameters from atom types using the topology's combining rule.

AMBER write (amberparams.cpp): Two conditions incorrectly skipped CLJScaleFactor(1.0, 1.0) pairs, causing them to be silently dropped and written as SCEE=0/SCNB=0 in the prm7 file. Fixed by correcting both conditions to exclude only fully-excluded (0.0, 0.0) pairs.

GROMACS read (gro87.cpp, grotop.cpp): Residue numbering that restarts at 1 (common in protein+glycan complexes) caused duplicate ResNum exceptions or multi-minute hangs due to O(n²) conflict resolution. Fixed by tracking residue boundaries via name+number rather than number alone, with O(1) unique numbering.

BioSimSpace merge (_merge.py, biosimspace.{h,cpp}): The merged molecule intrascale properties were constructed from global forcefield scale factors rather than the actual per-pair values stored in the molecule's intrascale property. For GLYCAM this caused an energy error in the merged end states. Fixed via a new C++ function SireIO::mergeIntrascale that directly copies per-pair scale factors into the merged molecule's atom index space using a two-pass approach, without relying on CLJNBPairs::merge which reconstructs pairs from global ff parameters.

  • I confirm that I have merged the latest version of devel into this branch before issuing this pull request (e.g. by running git pull origin devel): [y]
  • I confirm that I have added a changelog entry to the changelog (we will add a link to this PR as part of the review): [y]
  • I confirm that I have permission to release this code under the GPL3 license: [y]

@lohedges
Copy link
Contributor Author

It looks like the reindexing fix has broken one of my tests. It also appears that I need to make changes in AmberParams to account for the GLYCAM updates.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

bug Something isn't working

Projects

None yet

Development

Successfully merging this pull request may close these issues.

1 participant