Hello all, For some time now I have been unable to find a solution to my problem. I am attempting to do umbrella sampling on a protein ligand system. To generate the starting configurations, I use the following mdp file:
title = Protein-ligand complex MD Simulation define = -DPOSRES_SMD ;keeps from moving ; Run parameters integrator = md ; leap-frog integrator nsteps = 250000 ; 2 * 250000 fs = 500 ps (0.5 ns) dt = 0.002 ; 2 fs ; Output control nstxout = 0 ; suppress .trr output nstvout = 0 ; suppress .trr output nstenergy = 500 ; save energies every 1 ps nstlog = 500 ; update log file every 1 ps nstxtcout = 500 ; write .xtc trajectory every 1 ps energygrps = Protein_CA MOL ; Bond parameters continuation = yes ; first dynamics run constraint_algorithm = lincs ; holonomic constraints constraints = all-bonds ; all bonds (even heavy atom-H bonds) constrained lincs_iter = 1 ; accuracy of LINCS lincs_order = 4 ; also related to accuracy ; Neighborsearching ns_type = grid ; search neighboring grid cells nstlist = 5 ; 10 fs rlist = 0.9 ; short-range neighborlist cutoff (in nm) rcoulomb = 0.9 ; short-range electrostatic cutoff (in nm) rvdw = 1.4 ; short-range van der Waals cutoff (in nm) ; Electrostatics coulombtype = PME ; Particle Mesh Ewald for long-range electrostatics pme_order = 4 ; cubic interpolation fourierspacing = 0.16 ; grid spacing for FFT ;COM Pulling pull = umbrella pull_geometry = direction pull_start = yes pull_dim = Y Y Y pull_group0 = Protein_CA pull_group1 = MOL pull_vec1 = -2.4 -1.0 -0.4 pull_rate1 = 0.008 ; will pull 4nm during the 500ps run pull_k1 = 500 ; Temperature coupling is on tcoupl = nose-hoover tc-grps = Protein_MOL Water_and_ions ; two coupling groups - more accurate tau_t = 0.2 0.2 ; time constant, in ps ref_t = 300 300 ; reference temperature, one for each group, in K ; Pressure coupling is off pcoupl = Parrinello-Rahman ; pressure coupling is on for NPT pcoupltype = isotropic ; uniform scaling of box vectors tau_p = 2.0 ; time constant, in ps ref_p = 1.0 ; reference pressure, in bar compressibility = 4.5e-5 ; isothermal compressibility of water, bar^-1 ; Periodic boundary conditions pbc = xyz ; 3-D PBC ; Dispersion correction DispCorr = EnerPres ; account for cut-off vdW scheme ; Velocity generation gen_vel = yes ; assign velocities from Maxwell distribution gen_temp = 300 ; temperature for Maxwell distribution gen_seed = -1 ; generate a random seed A distance vs time graph of this 500ps simulation is here: http://dl.dropbox.com/u/32738905/Dist_SMD500_upped50.png I then run umbrella sampling on 21 of the 500 generated starting configurations to use as input for g_wham. The mdp file I use for the umbrella sampling is as follows: title = MD Run Umbrella define = -DPOSRES_SMD ; protein does not move ; Run parameters integrator = md ; leap-frog integrator nsteps = 1000000 ; 2 * 1000000 = 2000 ps (2ns) dt = 0.002 ; 2 fs tinit = 0 ; Output control nstxout = 50000 ; save coordinates every 10 ps nstvout = 50000 ; save velocities every 10 ps nstfout = 5000 ; every 10 ps nstxtcout = 5000 ; every 10 ps nstenergy = 5000 ; save energies every 10 ps nstlog = 5000 ; update log file every 10 ps energygrps = Protein_CA MOL ; Bond parameters continuation = yes ; first dynamics run constraint_algorithm = lincs ; holonomic constraints constraints = all-bonds ; all bonds (even heavy atom-H bonds) constrained lincs_iter = 1 ; accuracy of LINCS lincs_order = 4 ; also related to accuracy ; Neighborsearching ns_type = grid ; search neighboring grid cells nstlist = 5 ; 10 fs rlist = 0.9 ; short-range neighborlist cutoff (in nm) rcoulomb = 0.9 ; short-range electrostatic cutoff (in nm) rvdw = 1.4 ; short-range van der Waals cutoff (in nm) ; Electrostatics coulombtype = PME ; Particle Mesh Ewald for long-range electrostatics pme_order = 4 ; cubic interpolation fourierspacing = 0.16 ; grid spacing for FFT ; Temperature coupling is on tcoupl = nose-hoover tc-grps = Protein_MOL Water_and_ions ; two coupling groups - more accurate tau_t = 0.2 0.2 ; time constant, in ps ref_t = 300 300 ; reference temperature, one for each group, in K ; Pressure coupling is off pcoupl = Parrinello-Rahman ; pressure coupling is on for NPT pcoupltype = isotropic ; uniform scaling of box vectors tau_p = 2.0 ; time constant, in ps ref_p = 1.0 ; reference pressure, in bar compressibility = 4.5e-5 ; isothermal compressibility of water, bar^-1 ; Periodic boundary conditions pbc = xyz ; 3-D PBC ; Dispersion correction DispCorr = EnerPres ; account for cut-off vdW scheme ; Velocity generation gen_vel = no ; do not assign velocities ; Pull Code pull = umbrella pull_geometry = distance pull_start = yes pull_dim = Y Y Y pull_group0 = Protein_CA pull_group1 = MOL pull_rate1 = 0.0 pull_k1 = 25000 pull_nstxout = 1000 ; every 2 ps pull_nstfout = 1000 ; every 2 ps Once this is finished, I run g_wham and I get the following output files: Histogram: http://dl.dropbox.com/u/32738905/Histogram_deletable2.png Profile: http://dl.dropbox.com/u/32738905/Profile_deletable2.png Both graphs begin at 1.28nm which is the starting distance between the protein and the ligand but there is something very wrong with them. I was expecting something similar the output of Justin's umbrella sampling tutorial, since I based what I did on his tutorial. All that is different is I am using the AMBERsb99 forcefield and I am doing a protein-ligand system not a peptide-peptide system. If somebody could recommend something for me to do about this problem it would be greatly appreciated. Also, if you would like additional information about the problem feel free to email me. Thanks for your time, Eric LeGresley P.S. I'm sorry that what I have written may seem a little simple to some people, but I am 15 years old and I only started using GROMACS a few months ago. -- gmx-users mailing list [email protected] http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/Search before posting! Please don't post (un)subscribe requests to the list. Use the www interface or send it to [email protected]. Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
