# How to apply constant pull force at two ends of a molecule?

**URL:** https://gromacs.bioexcel.eu/t/how-to-apply-constant-pull-force-at-two-ends-of-a-molecule/10056
**Category:** User discussions
**Tags:** mdp-parameters, pull-code, simulation-setup
**Created:** [August 31, 2024, 11:19am UTC](https://gromacs.bioexcel.eu/t/how-to-apply-constant-pull-force-at-two-ends-of-a-molecule/10056 "2024-08-31T11:19:50Z")
**Posts on this page:** 3
**Page:** 1

<div class="post-metadata">

### Author: ![liuwq](https://avatars.discourse-cdn.com/v4/letter/l/c37758/32.png) [@liuwq](https://gromacs.bioexcel.eu/u/liuwq)
#### Post date: [August 31, 2024, 11:19am UTC](https://gromacs.bioexcel.eu/t/how-to-apply-constant-pull-force-at-two-ends-of-a-molecule/10056/1 "2024-08-31T11:19:50Z")

</div>

GROMACS version: 2018.8  
GROMACS modification: No  
I have a molecule, and I need to apply pulling forces with opposite directions to a head atom and a tail atom to stretch the molecule to an almost straight configuration, similar to that used in Figure 1 of the paper [J. Am. Chem. Soc. 2019, 141, 29, 11603–11613](https://pubs.acs.org/doi/abs/10.1021/jacs.9b04383).  
I tried with these pull settings, but it does not work as I expected.

```auto
pull = yes 
pull-ngroups = 2 
pull-group1-name = Left
pull-group2-name = Right
pull-ncoords = 2 
pull-coord1-type = constant-force
pull-coord1-geometry = direction-periodic
pull_coord1_vec = -1 0 0 
pull-coord1-dim = Y Y Y
pull_coord1_k = 10000
pull-coord1-groups = 1 2 
pull-coord2-type = constant-force
pull-coord2-geometry = direction-periodic
pull_coord2_vec = 1 0 0 
pull-coord2-dim = Y Y Y 
pull_coord2_k = 10000
pull-coord2-groups = 1 2 

```

The `Left` group is the head atom, the `Right` group is the tail atom.  
Could you please suggest to me the correct pull settings? So I can apply pulling forces with opposite directions at the two ends of a molecule.

---

<div class="post-metadata">

### Author: ![Karis](https://avatars.discourse-cdn.com/v4/letter/k/ed8c4c/32.png) [@Karis](https://gromacs.bioexcel.eu/u/Karis)
#### Post date: [September 2, 2024, 3:42am UTC](https://gromacs.bioexcel.eu/t/how-to-apply-constant-pull-force-at-two-ends-of-a-molecule/10056/2 "2024-09-02T03:42:22Z")

</div>

If the goal is to apply constant force to pull the 2 groups away from each other along a single direction, you should only need 1 pull coordinate (coord1 refers to the coordinate along which the groups are pulled, not the groups themselves). Also, if you’re using direction-periodic, you should only include pull\_coord1\_vec (setting pull\_coord1\_vec to -1 0 0 means group 1 is pulled in the negative x-direction and group 2 is pulled in the positive x-direction), no need to include pull\_coord1-dim (this is used for pull\_coord1\_geometry being set to distance, setting it to Y Y Y means it’s pulling in all 3 dimensions). For an example, see the mdp setup in this [tutorial](http://www.mdtutorials.com/gmx/umbrella/05_pull.html), and for a list of all possible mdp parameters, you can check the mdp options listed [here](https://manual.gromacs.org/documentation/2018/user-guide/mdp-options.html#com-pulling).

---

<div class="post-metadata">

### Author: ![liuwq](https://avatars.discourse-cdn.com/v4/letter/l/c37758/32.png) [@liuwq](https://gromacs.bioexcel.eu/u/liuwq)
#### Post date: [September 10, 2024, 6:58am UTC](https://gromacs.bioexcel.eu/t/how-to-apply-constant-pull-force-at-two-ends-of-a-molecule/10056/3 "2024-09-10T06:58:25Z")

</div>

Hello,  
Thank you very much for the detailed reply! I tried with your suggestions and made a new .mdp file, and it works.  
I have a molecule in a simulation box, and under NVT simulation, it will shrink or collapse. My goal is to pull the two ends of the molecule and to make the molecule straight again. In addition, I would like to know how much work is done by the pull force or the (potential) energy change of the molecule between the two states. However, I found it is hard to understand the output for me. Could you please help me further with this issue?  
The `COM Pull En.` from `gmx energy` is about `-335.808960 kJ/mol` at the last timestep. And from `pull_pullf.xvg` and `pull_pullx.xvg`, I found the force is recorded as `- 100 kJ/mol/nm` and the moving distance of the pulling force is `-3.35809 nm`. I think the work done by the pulling force is calculated by `F*r`, and it seems it is reflected by `COM Pull En.`. However, I do not understand why `COM Pull En.` is a negative value. The potential energy of the initial state is about 1050 kJ/mol, and the final value after pulling is about 700 kJ/mol, so roughly the new potential energy equals initial potential energy plus `COM Pull En.`, i.e., 700 = 1050 + (-335).  
Could you please tell me how can I calculate the work done by the pulling force or the potential energy difference between the state before pulling and after pulling? I expected the potential energy to increase after pulling, but it decreases according to the output from `gmx energy`.
