It's not US if you have a non-zero pull-rate1 that gets applied during your simulation, although with your pull-geometry=distance I don't think that pull-rate is even used (therefore I think that part is technically correct, although your .mdp file is confusing on this issue).
Careful that X and Y are going to be free to diffuse introducing what will likely be a convergence nightmare for separation distances at which the drug and protein are just coming into contact. Specifically, while you may set up the system with the drug and protein having the same COM in X and Y, that will change over time (more or less, depending on what define = -DPOSRES_Protein does for your system). This likely isn't a problem for large separation distances (although there will be some convergence issues here due to the need to sample over all X/Y at some larger displacements for which coulomb interaction are still non-negligable). Neither is it likely to be a problem for umbrellas in which the protein constrains your drug within a favourable binding pocket. However, in umbrellas where the drug and protein just come into contact, you will probably start with contact, but then need to obtain the proper boltzmann populations for interactions that depend on X/Y. What's worse , if you end up predominantly sampling conformational space that overlaps along your order parameter (Z) but not in some other orthogonal degree of freedom (e.g., X/Y), then your PMF will be simply wrong (that will get fixed with extensive sampling, but the amount of sampling that you need might be unobtainably large). Personally, I would solve this by adding a second pull group with a flat-bottom potential that acts only on X and Y to keep the drug within a cylinder that extends along Z away from your protein. However, you would need to modify gromacs source code for this (I have some posts on the uses lists over the years about how to do this). An alternative solution is to define an order parameter that acts in all 3 dimensions, although you may need to compute the free energy of releasing this restraint (a) in the binding pocket and (b) at large separations where your PMF flattens out to zero. Careful also about gen_vel = no (i.e., be sure to load in velocities if you are not generating them) Chris. ________________________________________ From: gromacs.org_gmx-users-boun...@maillist.sys.kth.se <gromacs.org_gmx-users-boun...@maillist.sys.kth.se> on behalf of Arunima Shilpi <writetoas...@gmail.com> Sent: 03 March 2014 09:34 To: gmx-us...@gromacs.org Subject: [gmx-users] md_pull code in umbrella sampling Dear Sir I am trying to calculate the potential mean force (PMF) between protein and Ligand. I have applied position restrain to protein. I have applied reference pulling group to Protein and Pulling force has been applied to ligand. The pul force is applied along Z-axis. I had query as to whether I am proceeding in the right direction. Here I have provided the detail of the content of md_pull.mdp md_pull.mdp title = Umbrella pulling simulation define = -DPOSRES_Protein ; Run parameters integrator = md dt = 0.002 tinit = 0 nsteps = 250000 ; 500 ps nstcomm = 10 ; Output parameters nstxout = 5000 ; every 10 ps nstvout = 5000 nstfout = 500 nstxtcout = 500 ; every 1 ps nstenergy = 500 ; Bond parameters constraint_algorithm = lincs constraints = all-bonds continuation = yes ; continuing from NPT ; Single-range cutoff scheme nstlist = 5 ns_type = grid rlist = 1.4 rcoulomb = 1.4 rvdw = 1.4 ; PME electrostatics parameters coulombtype = PME fourierspacing = 0.12 fourier_nx = 0 fourier_ny = 0 fourier_nz = 0 pme_order = 4 ewald_rtol = 1e-5 optimize_fft = yes ; Berendsen temperature coupling is on in two groups Tcoupl = Nose-Hoover tc_grps = Protein Non-Protein tau_t = 0.5 0.5 ref_t = 310 310 ; Pressure coupling is on Pcoupl = Parrinello-Rahman pcoupltype = isotropic tau_p = 1.0 compressibility = 4.5e-5 ref_p = 1.0 refcoord_scaling = com ; Generate velocities is off gen_vel = no ; Periodic boundary conditions are on in all directions pbc = xyz ; Long-range dispersion correction DispCorr = EnerPres ; Pull code pull = umbrella pull_geometry = distance ; simple distance increase pull_dim = N N Y pull_start = yes ; define initial COM distance > 0 pull_ngroups = 1 pull_group0 = Protein pull_group1 = DRG pull_rate1 = 0.01 ; 0.01 nm per ps = 10 nm per ns pull_k1 = 1000 ; kJ mol^-1 nm^-2 Regards Arunima -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org. -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
