Something to watch out for when we couple smirky moves with parameter moves

[davidlmobley]

One thing that just came up in discussions with Christopher Bayly is that we will have to be careful of something when we make moves with SMIRKY which relate to torsions and couple these to parameters. Particularly, consider a relatively generic torsion [*:1]~[#6:2]~[#6:3]~[*:4] and compare this to a (still fairly generic) torsion [#1:1]~[#6:2]~[#6:3]~[#1:4]. Imagine applying these to some ring system -- let's say cyclohexane for simplicity. Consider the torsions which apply to a C-C bond within the ring of cyclohexane. If we were considering the first pattern above, we would end up with many instances of it applying to that bond (H1-C-C-H1, H2-C-C-H1, H2-C-C-H2, H2-C-C-H1, C-C-C-H1, C-C-C-H2, H1-C-C-C, H2-C-C-C, ...), whereas the second pattern only applies to the H_-C-C-H_ instances. Notice that these two torsions apply to a different number paths involving the same two central atoms.

I haven't totally thought through the implications of this yet, but it strikes me that there is something to be careful of here. Particularly, I'm imagining we could relatively easily make moves from the generic to the more specialized torsion, but since we're changing how many paths the torsion applies to, I'm thinking that might mean we would need some corresponding change to the barrier heights? (i.e., this is what AMBER used the "idivf" factor to deal with, I think).

Anyway, just wanted to archive this year so I'll remember it if/when we get to this point.

[jchodera]

We discussed whether or not we should use for our new forcefield (1) explicitly specified idivf, (2) automatic division by the number of torsions along a bond, or (3) direct application of the torsion barrier as-is long ago. We selected option (3) for our new forcefield with support for (1) for backwards compatibility, but maybe we should rethink this. (1) would require moves that automatically scale torsions when idivf is changed, (3) would require some cleverness to not propose moves that change idivf (or else get high rejection rates), and (2) would not require any special considerations.

[davidlmobley]

I'm not totally sure I understand all the potential implications/issues here yet.

Chris had the idea that one solution to this would be to just not allow torsions like [*:1]~[#6:2]~[#6:3]~[*:4] but to instead require that each place have at least an element in it. However, it's not totally clear to me that this solves the problem either, as some torsions will certainly involve an "OR" of things so it seems to me that moves could still modify the number of terms applied.

Also, there's another place where I'm unclear. If I remember correctly, we require that EVERY set of atoms have a torsion, so in the example I gave above, we're not actually changing the number of torsions applied to that bond, we're just changing which torsions apply to that bond (i.e. we're changing a SMIRKS which applies to eight into one which applies to four, but at the same time there is ANOTHER torsion which I haven't listed which will have to pick up the other four). So, maybe there is no problem at all, in that we're just shifting the torsional coverage around.

On the other hand, it strikes me that this still may have implications for the types of moves we can propose/want to make in SMIRKS-space, since unless we do option (2) we might have to avoid having SMIRKS which cover patterns of different torsional multiplicity...

I still don't think option (1) is viable, since unless we make our SMIRKS particularly cumbersome, I don't think we'll necessarily KNOW the multiplicity of a given torsion (i.e. it might cover multiple multiplicities).

@jchodera -- thoughts?

I perhaps should discuss again with Christopher after he emerges a little more from parm@frosst torsions (900 to go). Though, perhaps this is urgent enough I should talk to him about it now. I'll try to do so.

[davidlmobley]

OK, I discussed this with Christopher and the take-away, basically, is that we are already doing things exactly the way we should be.

Some summary:

• Yes, my concern is absolutely right that if we change the connectivity of the central atoms in a torsion, we will change the number of terms we get applied and that would reduce the likelihood of being able to find good parameters for that. But, that's what we want.
• He has already baked this above concern into the move proposal probabilities he shared with us (low probability of proposing a change to the connectivity of central atoms in a torsion since this normally won't be a good thing)
• The other issues are all handled properly just by solution (3) noted above, which we have currently.

So, I think the main take-away is that he thinks that central atoms in a torsion should basically always have connectivity, and we probably don't want many torsions which don't contain at least some specification of element somewhere.

I'll mark this as resolved and just save for the record.

[jchodera]

Did you guys decide we don't want to change the current behavior? I think it may be worth adding support for the automated division by number of incident torsions along a bond.

[davidlmobley]

@jchodera - yes, the current behavior seems to be exactly what we want.

[davidlmobley]

I'll have to elaborate on this more later, but right before @cbayly13 left, we noticed that we really are going to strongly prefer option 2 here:

2. automatic division by the number of torsions along a bond

for reasons which only became obvious at the very end -- i.e. relating to changes in protonation of certain groups which could modify the number of torsions being applied in ways which they shouldn't. If we do this division (effectively, a type of averaging), it will avoid the problem in a sensible way.

I'm reopening this issue as we'll need to implement automatic division by the number of torsions.

There is some overlap between this and #51 because in #51 , I'm arguing for automatic division by the number of impropers applied (6) which would be consistent with automatic division by the number of propers applied, here.