Coulomb PME convergence criteria no specified?

[mrshirts]

I don't think PME convergence criteria are specified (gaussaina smoothing, number of grid points, etc). Most codes different in what convergence criteria they use, so one could get relatively different answers without specifying this.

[jchodera]

I believe our general philosophy is to define the force field for the convergence limit for (Particle Mesh) Ewald-based methods, and allow some level of tolerance in implementations as long as it does not deviate too much from the intended energies.

However, we don't specify how much deviation is "too much" or what metrics of deviation are to be measured. I think this is an area where we could develop some useful tools in the future.

For now, we could place a limit on the energy deviation stddev for a benchmark set of systems, for example. Or, as you suggest, specify a PME cutoff tolerance.

[mattwthompson]

SMIRNOFF currently offers no guidance beyond

a method like particle mesh Ewald should be used

In the current implementation .setEwaldErrorTolerance(1.0e-4) is used in the OpenMM export. The reason for this choice is historical and I can't argue for it in favor of other values or possible combinations of settings. These settings are not part of the Amber or GROMACS exports as the Ewald settings are not specified in files provided. The case of LAMMPS is a little different but not the most important one to consider here.

[mattwthompson]

I run into this once every few months and it causes a minor headache each time, purely due to 0.005, 0.0001, and even 0.00001 being used almost interchangeably by different tools.

With all due respect to the complexity here, having this specified in the force field would make it less awkward for me to change the defaults, since that qualifies as a breaking change that concretely impacts energy evaluations (and there potentially any of the physics).

[mikemhenry]

We should use whatever was used to parameterize the force field (or used to validate it).

[mrshirts]

I think the thing is that all of these choices should result in very similar energies, so people don't police what this number is very much.

[jchodera]

I think the thing is that all of these choices should result in very similar energies, so people don't police what this number is very much.

That's not quite correct. In particular, 0.005 = 5e-3 is far too lax, and would likely result in significant issues in free energy calculations. The potential energy stddev for the OpenMM default, 5e-4, was >> 1 kT for solvated protein:ligand systems.

The best approach here is to have the force field be built for Ewald (the PME accuracy limit) with an "acceptable error" limit on how much one can relax the tolerance dependent on application and system size. For example, for reliable free energy calculations, you'd want to tune the PME tolerance so that your potential energy stddev is << 1 kT to ensure that the sampled probability density is not significantly perturbed from the target density. This will depend on system size, composition, etc.

I realize this is different from current practice, but current practice has no actual rationale underlying it that I can discern.