Reality check: time needed for FES convergence with WT Metadynamics

[David M. Sherman]

So, I'm looking at metal-chloride complexation in NPT_I aqueous solutions (55 H2O in box). I'm using the PBE functional with Grimme dispersion. Because I'm doing NPT, I use a plane-wave cutoff to 400/40; the barostat performs very nicely with average pressures fluctuating within a few % of my target (simulations run from 1 to 3000 bars with temperatures ranging from 50 to 350C).  My colvars are M-Cl CN and M-O CN with SCALE small enough to nicely resolve the FES mimina.  I'm using well-tempered metadynamics with a bias factor (T+DeltaT/T) of 10.  I'm adding hills every 150 steps (0.5 fs stepsize).  I'm using 10 walkers.  From the differences in FES minima, I calculate logK values for the various M+nCl = MCln reactions. Qualitatively, the results are pretty good: even after only 50ps per walker, most of the the logK values are within 1-2 units of experiment at 50C, 1 bar and show a nice systematic (and physically correct) change with pressure and temperature.  However, even after each walker has run for  150 ps (i.e., 1.5 ns for each P,T run), I still have the logK values drifting over  +/-  1-2 log units every 20ps or so; this corresponds to fluctuations in the differences between  FES minima of 10-20 kJ/mole.  After 150ps per walker, the added hill height is still about 0.00002 AU (0.05 k/mole).  Given that I'm using WT Metadynamics,  I know that, eventually,   the added hills will be 0, and the FES will become constant, but why is this taking so long?    I know I'm slowing things down by also including the oxygen (water) CN as a collective variable.   Indeed, it seems like the significant hills  that are being added after 150 ps are just for non-important hydration complexes; still, the FES mimima for the MCl, MCl2, MCl3 complexes  haven't really converged either (the free surfaces are clearly on their way to convergence; however the energy differences between the free energy minima are still too noisy).  I know that if I have poorly defined colvars, I will never get FES convergence, but these colvars are well-defined relative to the system and, indeed, give very clear minima in the FES with a barrier height between the minima of around 10-30 kJ/mole. So, my question is: should we expect a highly converged (+/- 5kJ/mole="chemical accuracy") FES in a system like this to require more than 1-2 nanoseconds (100-200 ps/walker with 10 walkers)?  Is everything here normal? Alternatively, Is my bias factor too high?  Should I increase NT-HILLS (am I "hill surfing"?). Thanks!

[David M. Sherman]

Well, I answered my own question: for a system like this, it took about 220 ps for each walker (2.2 ns total) to get log K values converged to +/-0.5 units. The results are pretty interesting and show where both DFT and dielectric continuum models fail (and succeed).  Paper in geochemical literature is forthcoming....