Dependece of diffusion constant on force constant in AWH simulation

GROMACS version: 2024.4
GROMACS modification: No
Hi,

I am simulating a particle permeating a membrane using the martini2 force field and wanted to compute the permeability of the membrane. To do this, I am using AWH to obtain the PMF and diffusion constant along Z.

While i was doing a test run on diffusion coefficient estimation in an isotrophic water box containing the solute particle, I found that the diffusion constant D_z from the simulation is highly dependent on the chosen force constant of the potential. The D_z decreases with increasing force constant. In my expectation, the estimation of D_z should be relatively stable regardless of the choice of force constant.

Any comments or suggestions are appreciated.Thanks in advance for any help.

The attached is the results and the settings for the AWH simulation is as follows:
pull = yes
pull_ncoords = 1 ; only one reaction coordinate
pull_ngroups = 2 ; two groups defining one reaction coordinate
pull_group1_name = W
;pull-group1-pbcatom = 1198
;pull-group2-pbcatom = 12211
;pull-pbc-ref-prev-step-com = yes
pull_group2_name = EOL
pull_coord1_type = external-potential
pull-coord1-potential-provider = AWH
pull_coord1_geometry = direction ; simple distance increase
;pull-cylinder-r = 2
pull_coord1_dim = N N Y
pull_coord1_vec = 0 0 1
pull_coord1_groups = 1 2
pull_coord1_start = yes ; define initial COM distance > 0
pull-nstxout = 1000
pull-nstfout = 0

awh = yes
awh-nstout = 50000
awh-nbias = 1
awh-nstsample = 100
awh-nsamples-update = 10
awh-potential = convolved

awh1-ndim = 1
awh1-dim1-coord-provider = pull
awh1-dim1-coord-index = 1
awh1-dim1-start = -2
awh1-dim1-end = 0
awh1-dim1-force-constant = 10000
awh1-dim1-diffusion = 4e-3

Firstly, I recommend that you use geometry=cylinder. Now the whole membrane might deform, which can give artifacts in the PME and diffusion constant.

The diffusion constant will depend (strongly) on the force constant when the force constant is comparable to or smaller than the curvature of the free-energy landscape. So the question is what curvature you have. I would guess though that a force constant of 10000 should be sufficient, but there is still a strong change when going to 20000.

Thanks Hess for the detailed exaplanation. I will use the geometry=cylinder as you suggest for membrane systems.

The test systems and the above results are obtained in a simple anisotropic water box containing a solute without the membrane. The a Z span of ~ 1 nm is used for running AWH for the solute. The simulation quickly exit the initial stage and I did not run too long. Might it help if I run longer?

Thanks again.

But then you are measuring something very strange. I don’t understand what you are pulling against then.

Ah, I was pulling using an absolute coordinat, which may not be optimal. I also observed this dependence in a simulation with membrane. The solute in the bulk water (not interacting with membrane) has quite different D_z when using different force constants. I’ll double check and thanks!