Support 4-site water models in Amber export

[mattwthompson]

Description
This should be simpler than implementing all virtual site types i.e. for ligands, I think.

[mattwthompson]

I'm trying to figure out how this is implemented in Amber files and I'm getting pretty confused pretty quickly.

• There's a NUMEXTRA flag, but it's not clear what it's used for. It doesn't seem to be used anywhere but the top of the file. Maybe the integrator just looks for massless particles and assumes they're all virtual sites.
• Amber seems to want to add valence terms between the virtual site and its parent atoms and orientation atoms, with no obvious indication that the participating atom is a virtual site.
  • Confusingly, this seems to be how fixed water models are implemented - by writing a bond between the hydrogens
• Exclusions for virtual sites seem to also treat them just the same as any other atom.

I seem to be in good company with my confusion. Section 17.5 of the Amber22 manual begins

It is not completely feasible to perform molecular dynamics with massless particles. However, for many useful
cases in which the locations of massless particles are determined by the locations of two or more atoms with mass,
it is possible to perform dynamics by using the chain rule to transfer forces from the “virtual sites” to the massive
particles. These constructions, enumerated below, provide a means for breaking out of the “one atom, one site”
paradigm that has dominated classical molecular dynamics. The prmtop format utilized by the sander and pmemd
programs does not always provide a straightforward means of expressing the relationships between virtual sites
and their parent (or “frame”) atoms. In Amber20, the sander and pmemd programs only support the most widely
used cases of virtual sites (e.g. TIP4P and TIP5P water), but efforts are underway to support a broader variety of
these virtual sites

There's a lot of useful information in that chapter on how mdgx can be used to fit force fields with virtual sites, but that's not what I'm trying to do here. I'm just trying to find a way to get virtual sites into the relevant input files and get sander/pmemd to understand that they're there.

[mattwthompson]

I'm still super perplexed at how the geometry of the virtual site is kept in place. There's a soft convention of adding a bond between the virtual site and its parent atom, which might constrain that particular distance, but what about the rest of the geometry? In other words, what's keeping the virtual site in 4-site water in the plane?

[mrshirts]

I wonder if we can ask Jason Swails this?

[mattwthompson]

@swails could you suggest a starting point for us to understand how various engines within Amber expect to see virtual sites structured in files? I can't find a lot on the Amber website.

I know mdgx can do plenty of fancy stuff in the context of parameter development, but for now the goal is just supporting 4-site water for which we already know the parameters. Writing out files without virtual sites and a "rule" file and then calling a separate program to smoosh them together is an option I'd like to avoid if possible.