AWH simulation with cylinder pull geometry

GROMACS version: 2025.4
GROMACS modification: No

Hi all,

we are simulating the membrane insertion of a cyclic peptide (with a backbone diameter of 0.8 nm) into a lipid bilayer using AWH with cylinder pulling.

The initial geometry is:

  • Z COM Peptide: 9.78 nm
  • Z COM Membrane: 5.2 nm
  • Z COM Lipids below peptide: 5.3 nm

Therefore, the expected distance between the peptide and the local membrane is approximately 4.5 nm. Since it was recommended at several posts in this forum, we would like to use the cylinder pull coordinate. However, we see strange results depending on the radius of the cylinder. We would assume naively that an increasing radius converges to the distance calculated with the direction pull coordinate. But apparently, it is vice versa. Did we understand it correctly that the first group in pull-coord1-groups must be the membrane?

Geometry pull-cylinder-r Distance at start
cylinder 5.0 nm 1.178 nm
cylinder 1.5 nm 1.644 nm
cylinder 0.1 nm 4.130 nm
direction -- 4.482 nm

Below I add the code for cylinder and radius. Why would this happen in the cylinder simulations?

#----------------------------------CYLINDER---------------------------

Radius 5.0 nm:

; Pull code
pull = yes
pull-ncoords = 1 ; only one reaction coordinate
pull-ngroups = 2 ; two groups defining one reaction coordinate
pull-group1-name = onlyMEMB ; pull towards membrane
pull-group2-name = LUGDBB ; pull the peptide COM
pull-coord1-type = external-potential ; Apply the bias using an external module.
pull-coord1-potential-provider = AWH ; The external module is called AWH!
pull-coord1-geometry = cylinder
pull-cylinder-r = see above
pull-coord1-dim = N N Y
pull-coord1-groups = 1 2
pull-coord1-vec = 0 0 -1 ; pull in z direction
pull-coord1-start = yes
pull-group1-pbcatom = 5902 ; roughly in the middle- took into account only for first step so not super-important
pull-pbc-ref-prev-step-com = yes ; the ref is adapted every step- necessary for large groups like membranes

#----------------------------------DIRECTION---------------------------

pull = yes
pull-ncoords = 1 ; only one reaction coordinate
pull-ngroups = 2 ; two groups defining one reaction coordinate
pull-group1-name = onlyMEMB ; pull towards lipid A phosphates
pull-group2-name = LUGDBB ; pull the peptide COM
pull-coord1-type = external-potential ; Apply the bias using an external module.
pull-coord1-potential-provider = AWH ; The external module is called AWH!
pull-coord1-geometry = direction
pull-coord1-dim = N N Y
pull-coord1-groups = 1 2
pull-coord1-vec = 0 0 -1 ; pull in z direction
pull-coord1-start = yes
pull-group1-pbcatom = 5902 ; roughly in the middle- took into account only for first step so not super-important
pull-pbc-ref-prev-step-com = yes ; the ref is adapted every step- necessary for large groups like membranes

Thanks a lot!

Your results look very strange. Unless the membrane is not flat, you should get the COM distance between the two groups independent of the radius. You have pull-pbc-ref-prev-step-com set, which I think should not be necessary. Could you try without that?

That is what we also expected. We have a flat membrane with a size of 8 nm x 8 nm, so there should not be a lot of curvature.

We get the same results without using pull-pbc-ref-prev-step-com.

I have no clue about what could be going wrong. Could you provide all input files to run grompp?

Thank you for your help! Can we send it to you via email?

If everything fits in a mail, yes.

We think we found the problem. It seems like we did not have enough solvent above the peptide.

If we reduce the peptide Z COM position to 8 nm, we do not get these issues.

Is there some way to get around this without increasing the amount of solvent around the membrane? What is the reason for this happening?

What do you mean with “solvent”? I suppose you mean periodic distance. You can set the start distance manually, I think that should solve it. But you want you peptide far enough away from the membrane over the periodic boundary.

Yes, the problem probably originates from the box being too small in the z-dimension.
Could be it be related to this note that we got?

NOTE 1 [file walker1/awh1.mdp, line 104]:
With pull type ‘external-potential’ and geometry ‘cylinder’, the distance
component along the cylinder axis between atoms in the cylinder group and
the COM of the pull group should be smaller than half the box length

The COM of the membrane and the peptide is below half the box height (10.8 nm / 2 = 5.4 nm > 4.5 nm). However, some atoms in the opposite leaflet will most likely be further away than this. If we move the peptide to an initial COM-COM distance of 3.6 nm the output of GROMPP is fine:

From GROMPP:

Cylinder radius = 1.5 nm

Pull group natoms pbc atom distance at start reference at t=0
1 25644 5148
2 27 73832 -3.234 nm -3.234 nm

COM-COM distance calculated with ‘gmx traj’ using the whole membrane as reference: 3.56376 nm

The problem is that with a COM-COM distance of 3.6, the peptide isn’t far enough from the membrane to generate a good PMF curve.
Is this an intrinsic limitation of this method, or does it sound like we have a different issue in our setup?

Thank you for your input!

I think it might work if you turn of start and set the reference distance manually. Then the pull code will search for the peptide atoms near the reference, if I am not mistaken.

If this works, we should add this to the manual. This issue comes up more often.

I am not sure if I understood correctly, but I made the following changes:

In the .mdp file:
pull-coord1-start = no
pull-coord1-init = -4.5

Then, when calling grompp I get the output:

Pull group natoms pbc atom distance at start reference at t=0
1 25644 5902
2 27 73823 -1.644 nm -4.500 nm

So the start distance is still incorrect even though the reference is set.

I sent a zip file to your email address with a minimal reproduction, I hope this makes the issue clear.

I tried your system and looked in the code and I see what the issue is now. Indeed my suggestion works when only using the pull code. But not with AWH. A comment in the code says:

/* With an external pull potential, the reference value loses
* it’s meaning and should not be used. Setting it to zero
* makes any terms dependent on it disappear.
* The only issue this causes is that with cylinder pulling
* no atoms of the cylinder group within the cylinder radius
* should be more than half a box length away from the COM of
* the pull group along the axial direction.
*/

I must have written this code and comment. But that was a long time ago. I need to see if there is a solution for this.

Thank you very much for you time, I am glad you could find the issue.

Since our membrane is quite small we will just keep using direction pulling for now.

I have found a solution. If you want to use it you can copy the two changed files into a GROMACS 2026 distribution. Only applying the changes in the files should also work in 2025.

We incorporated the changes into the 2026.1 GROMACS distribution, and the cylinder pulling now works as expected with the following mdp options:

; Pull code
pull                           = yes
pull-ncoords                   = 1                 
pull-ngroups                   = 2                  
pull-group1-name               = onlyMEMB            
pull-group2-name               = LUGDBB              
pull-coord1-type               = external-potential 
pull-coord1-potential-provider = AWH                
pull-coord1-geometry           = cylinder
pull-cylinder-r                = 1.5
pull-coord1-dim                = N N Y
pull-coord1-groups             = 1 2
pull-coord1-vec                = 0 0 -1           
pull-coord1-start              = yes

pull-group1-pbcatom            =  5902 
pull-pbc-ref-prev-step-com     = no 

It was important that pull-pbc-ref-prev-step-com not be set to yes.

We initially tried to incorporate the changes in the 2025.4 version, but there was some issue, but since it is working in 2026.1, it is fine.

Thank you for the patch!

1 Like