# PMF from SMD simulations using constant force option **URL:** \ **Category:** User discussions\ **Tags:** mdp-parameters, analysis-tools, simulation-setup\ **Created:** [July 2, 2024, 6:02am UTC](https://gromacs.bioexcel.eu/t/pmf-from-smd-simulations-using-constant-force-option/9552 "2024-07-02T06:02:31Z")\ **Posts on this page:** 7\ **Page:** 1
**Author:** ![Mohan](https://avatars.discourse-cdn.com/v4/letter/m/a88e57/32.png) [@Mohan](https://gromacs.bioexcel.eu/u/Mohan)\ **Post date:** [July 2, 2024, 6:02am UTC](https://gromacs.bioexcel.eu/t/pmf-from-smd-simulations-using-constant-force-option/9552/1 "2024-07-02T06:02:31Z")
GROMACS version: GROMACS modification: Yes/No Here post your question Hello Gromacs users, I would like to perform steered MD simulations using gromacs. I have followed umbrella sampling tutorial in which first the protein is pulled to generate conformations which then ran and analysed using wham. I want to pull small molecule through a channel, for which I want to use option constant force for my pulling, which gives two output xvg files force.xvg, position.xvg. My question is, how to get PMF profile from these output files? ---
**Author:** ![milosz.wieczor](https://avatars.discourse-cdn.com/v4/letter/m/b19c9b/32.png) [@milosz.wieczor](https://gromacs.bioexcel.eu/u/milosz.wieczor)\ **Post date:** [July 3, 2024, 11:40am UTC](https://gromacs.bioexcel.eu/t/pmf-from-smd-simulations-using-constant-force-option/9552/2 "2024-07-03T11:40:08Z")
Are you doing umbrella sampling or steered MD? Steered MD is often used as a preparatory step for US, but you cannot really get any reasonable free energy estimate from a single steered MD run (even with Jarzynski you’d need at least dozens, if not thousands). ---
**Author:** ![Mohan](https://avatars.discourse-cdn.com/v4/letter/m/a88e57/32.png) [@Mohan](https://gromacs.bioexcel.eu/u/Mohan)\ **Post date:** [July 4, 2024, 3:58am UTC](https://gromacs.bioexcel.eu/t/pmf-from-smd-simulations-using-constant-force-option/9552/3 "2024-07-04T03:58:16Z")
Thank you so much for reply. I am doing steered MD simulations, pulling a small molecule through a channel. I have the following question, 1. My reaction coordinate is com of small molecule and com of dummy kept at just below the center of the channel. The initial distance between them is 38 Angstorms, i have used following pull.mdp file integrator = md dt = 0.002 nsteps = 1500000 nstxout-compressed = 500 ;nstvout = 50000 ;nstfout = 50000 nstcalcenergy = 5000 nstenergy = 5000 nstlog = 5000 ; cutoff-scheme = Verlet nstlist = 20 rlist = 1.2 vdwtype = Cut-off vdw-modifier = Force-switch rvdw\_switch = 1.0 rvdw = 1.2 coulombtype = PME rcoulomb = 1.2 ; tcoupl = nose-hoover tc\_grps = Protein wat\_cap tau\_t = 1.0 1.0 ref\_t = 310 310 constraints = h-bonds constraint\_algorithm = LINCS continuation = yes ;; Pull code pull = yes pull-coord1-dim = N Y N pull\_ncoords = 1 ; only one reaction coordinate pull\_ngroups = 2 ; one c6 and three dummies pull\_coord1\_groups = 1 2 pull\_group1\_name = HEP pull\_group2\_name = DUM pull\_coord1\_type = umbrella ; harmonic potential pull\_coord1\_geometry = distance ; simple angle increase pull-nstfout = 10 pull-nstxout = 10 pull-print-components = no pull-print-ref-value = yes pull\_coord1\_start = yes ; define initial COM distance \> 0 pull\_coord1\_rate = -0.001 ; 0.01 nm per ps = 10 nm per ns pull\_coord1\_k = 1000 ; kJ mol^-1 nm^-2 pull-coord1-vec = 0 -1 0 pull-pbc-ref-prev-step-com = yes according to the above mdp file, the reaction coordinate should cover 3 nm but when I use wham to get pmf. I get pmf from 3.8 nm to 4.2 nm, rest of the reaction coordinate is not covered? Could you please let me know why is it not covering entire reaction coordinate? ---
**Author:** ![milosz.wieczor](https://avatars.discourse-cdn.com/v4/letter/m/b19c9b/32.png) [@milosz.wieczor](https://gromacs.bioexcel.eu/u/milosz.wieczor)\ **Post date:** [July 5, 2024, 10:25am UTC](https://gromacs.bioexcel.eu/t/pmf-from-smd-simulations-using-constant-force-option/9552/4 "2024-07-05T10:25:46Z")
Your vector is towards negative values of Y, and your pulling rate is negative, so you end up pulling the molecule towards higher values of Y. Remove one of the minus signs and the direction should be fine. But again, steered MD by itself does not give you the PMF, make sure you run a free energy method such as US after your SMD. ---
**Author:** ![hess](https://dub1.discourse-cdn.com/flex017/user_avatar/gromacs.bioexcel.eu/hess/32/416_2.png) [@hess](https://gromacs.bioexcel.eu/u/hess)\ **Post date:** [July 5, 2024, 10:33am UTC](https://gromacs.bioexcel.eu/t/pmf-from-smd-simulations-using-constant-force-option/9552/5 "2024-07-05T10:33:37Z")
With steered MD the only option to get a PMF is using Jarzynski’s equality. You need (very) many simulations though. ---
**Author:** ![Mohan](https://avatars.discourse-cdn.com/v4/letter/m/a88e57/32.png) [@Mohan](https://gromacs.bioexcel.eu/u/Mohan)\ **Post date:** [July 12, 2024, 7:06am UTC](https://gromacs.bioexcel.eu/t/pmf-from-smd-simulations-using-constant-force-option/9552/6 "2024-07-12T07:06:47Z")
Thank you so much. I am running simulation and will let you know ---
**Author:** ![Mohan](https://avatars.discourse-cdn.com/v4/letter/m/a88e57/32.png) [@Mohan](https://gromacs.bioexcel.eu/u/Mohan)\ **Post date:** [July 12, 2024, 7:09am UTC](https://gromacs.bioexcel.eu/t/pmf-from-smd-simulations-using-constant-force-option/9552/7 "2024-07-12T07:09:37Z")
Thank you so much for a reply. If I Run 20 simulations using the above pull code for 1 ns each. after obtaining force.xvg, pos.xvg files as output. How can we convert using an equation? could you please provide stepwise details or direct me to any of the papers which contain step-by-step details?