Project

General

Profile

Feature #2451

Linear virtual sites with fixed distance

Added by David van der Spoel about 1 year ago. Updated 6 months ago.

Status:
Closed
Priority:
Normal
Category:
mdrun
Target version:
-
Difficulty:
simple
Close

Description

In the Charmm force field linear virtual site are use for halogens to model halogen bonding. By design the distance between vsite particle and halogen atom is fixed at
CL-v: 0.016 nm,
BR-v: 0.019 nm, and
I-v: 0.022 nm

The present virtual_site2 does not support this since the position of the vsite is relative to two atoms and hence the distance can fluctuate. Here the proposal is make an additional virtual_site (type 2) with fixed distance to the halogen.

History

#1 Updated by Berk Hess about 1 year ago

Is the vsite in between the two atoms or outside? Should we support both (using the sign of the distance)?

#2 Updated by David van der Spoel about 1 year ago

Outside the halide atom I think, maybe Justin can comment, but this is a matter of definition for the implementation as well. I envision one could write
[ virtual_sites2 ]
; i j k type distance
Vsite Cl C 2 0.016

#3 Updated by Justin Lemkul about 1 year ago

David van der Spoel wrote:

Outside the halide atom I think, maybe Justin can comment, but this is a matter of definition for the implementation as well. I envision one could write
[ virtual_sites2 ]
; i j k type distance
Vsite Cl C 2 0.016

Yes, this is correct. In our convention, the construction is C-Cl-X. There could be other use cases, of course, and I agree that this could just be done with a sign convention (>0 for "outside" the bond, <0 for "inside").

#4 Updated by Berk Hess about 1 year ago

I would think the consistent sign definition with the normal 2-atom vsite would be negative outside.

#5 Updated by Justin Lemkul about 1 year ago

However we want to define it is fine, though it is a little confusing to have a negative number extend beyond the end of a bond. But I guess that's what documentation is for :)

#6 Updated by Berk Hess about 1 year ago

A normal 2-atom vsite the parameter a, as one would expect, interpolates linear between atom i and j. So to get an atom sticking out from i you need a negative a.

#7 Updated by David van der Spoel about 1 year ago

Just wondering whether fixing the 3FAD such that it can deal with 180 degree angles wouldn't be sufficient? It would do the trick for everything except diatomics like HF and CO.

#8 Updated by Justin Lemkul about 1 year ago

David van der Spoel wrote:

Just wondering whether fixing the 3FAD such that it can deal with 180 degree angles wouldn't be sufficient? It would do the trick for everything except diatomics like HF and CO.

From my perspective, 3fad would cover things. I'm not sure how many existing force fields will require this kind of construction for diatomics. I know CHARMM/CGenFF do not.

#9 Updated by Berk Hess about 1 year ago

Looking at the manual I realized that you can already achieve what you want by using 3fd with j=k (and arbitrary a, e.g. a=0).
Implementing 2fd is then as simple as cut-and-pasting the 3fd code and removing k and a.
This also means that for consistency we should not have a minus sign for the "out" location.

#10 Updated by Mark Abraham 7 months ago

  • Status changed from New to Closed
  • Target version deleted (2019)

Sounds like there is a solution that doesn't need new GROMACS code.

#11 Updated by Justin Lemkul 6 months ago

To follow up on this issue, I have tried to do what Berk suggested with 3fd construction, and it indeed works just fine. The only issue is that grompp throws a warning about a duplicate atom number in the virtual site construction. That's fair, but should we add an exception, e.g. in the case of a=0, to allow for such odd behavior? I can get around this with -maxwarn 1, but I generally think doing so is poor practice and encourages users to ignore more significant issues.

#12 Updated by Justin Lemkul 6 months ago

Even worse, with version 2018, the virtual site construction throws an error that cannot be circumvented. I haven't tried the 2019 beta, but should I expect different behavior? This lone pair construction is essential for CGenFF going forward and our users are seeking solutions.

#13 Updated by David van der Spoel 6 months ago

Then let's get back to the original plan an implement a new type, 2fd. So we need a new type and a description here:

http://manual.gromacs.org/documentation/2019-beta2/reference-manual/functions/interaction-methods.html#virtual-interaction-sites

Also available in: Atom PDF