[gmx-users] Md simulation mdp options error
Hi, I am trying to run a simulation with the following optins in my production mdp file. Previously I have run the simulations with the production run parameters using V-rescale thermostat and berendsen barostat. But now recently I read in discussion that Justin had advised using Parinello-Rahman barostat as it would be much more suitable and provides adequate sampling. Also, I had to use the Verlet cut-off scheme because the cluster where I run my production simulation has updated gomacs to the recent version 4.6.5 which requires one to use the verlet scheme. My mdp file is below title= production MD ; Run parameters integrator= md; leap-frog algorithm nsteps= 5000; 0.002 * 5000 = 10 ps or 100 ns dt = 0.002 ; 2 fs ; Output control nstxout= 0; save coordinates every 2 ps nstvout= 0; save velocities every 2 ps nstxtcout= 1000; xtc compressed trajectory output every 5 ps nstenergy= 1000; save energies every 5 ps nstlog= 1000; update log file every 5 ps ; Bond parameters constraint_algorithm = lincs; holonomic constraints constraints= all-bonds; all bonds (even heavy atom-H bonds) constraine d lincs_iter= 1; accuracy of LINCS lincs_order= 4; also related to accuracy ; Neighborsearching ns_type= grid; search neighboring grid cells nstlist= 5; 25 fs rlist= 1.0; short-range neighborlist cutoff (in nm) rcoulomb= 1.0; short-range electrostatic cutoff (in nm) rvdw= 1.4; short-range van der Waals cutoff (in nm) rlistlong= 1.0; long-range neighborlist cutoff (in nm) cutoff-scheme = Verlet ; Electrostatics coulombtype= PME; Particle Mesh Ewald for long-range electrostat ics pme_order= 4; cubic interpolation fourierspacing= 0.16; grid spacing for FFT nstcomm = 10; remove com every 10 steps ; Temperature coupling is on tcoupl= V-rescale; modified Berendsen thermostat tc-grps= Protein Non-Protein; two coupling groups - more accurate tau_t= 0.10.1; time constant, in ps ref_t= 318 318; reference temperature, one for each group, in K ; Pressure coupling is off pcoupl = Parrinello-Rahman ; pressure coupling is on for NP T 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; Velocity generation is on gen_temp= 318; reference temperature, for protein in K ERROR 1 [file md.mdp]: With Verlet lists rcoulomb!=rvdw is not supported NOTE 1 [file md.mdp]: With Verlet lists the optimal nstlist is = 10, with GPUs = 20. Note that with the Verlet scheme, nstlist has no effect on the accuracy of your simulation. NOTE 2 [file md.mdp]: nstcomm nstcalcenergy defeats the purpose of nstcalcenergy, setting nstcomm to nstcalcenergy WARNING 1 [file md.mdp]: You are generating velocities so I am assuming you are equilibrating a system. You are using Parrinello-Rahman pressure coupling, but this can be unstable for equilibration. If your system crashes, try equilibrating first with Berendsen pressure coupling. If you are not equilibrating the system, you can probably ignore this warning. Velocities were taken from a Maxwell distribution at 318 K Removing all charge groups because cutoff-scheme=Verlet There were 2 notes There was 1 warning --- Program grompp_d, VERSION 4.6.5 Source code file: /home/y07/y07/gmx/4.6.5-phase1/source/src/kernel/grompp.c, line: 1610 Fatal error: There was 1 error in input file(s) For more information and tips for troubleshooting, please check the GROMACS website at http://www.gromacs.org/Documentation/Errors --- So, I guess, NOTES 1 and 2 are harmless and can be ignored, right? Regarding error 1, I think I would have to change my rvdw to 1.0 and that should sort the problem out. I wanted to get advice on the warning too, should I ignore and proceed? Or should I go back to berendsen as I have been using for my previous simulations. For NVT and NPT, I have used berendsen. As this is a very crucial simulation run, I just wanted to confirm since I have not had any previous experience with Parinello-Rahman so, any advice would surely help me. -- Ankita Naithani -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support
Re: [gmx-users] Md simulation mdp options error
Thank you Justin for your reply. Yes indeed this is a production run. So would the following mdp file be okay now? title= production MD ; Run parameters integrator= md; leap-frog algorithm nsteps= 5000; 0.002 * 5000 = 10 ps or 100 ns dt = 0.002 ; 2 fs ; Output control nstxout= 0; save coordinates every 2 ps nstvout= 0; save velocities every 2 ps nstxtcout= 1000; xtc compressed trajectory output every 5 ps nstenergy= 1000; save energies every 5 ps nstlog= 1000; update log file every 5 ps ; Bond parameters constraint_algorithm = lincs; holonomic constraints constraints= all-bonds; all bonds (even heavy atom-H bonds) constraine d lincs_iter= 1; accuracy of LINCS lincs_order= 4; also related to accuracy ; Neighborsearching ns_type= grid; search neighboring grid cells nstlist= 5; 25 fs rlist= 1.0; short-range neighborlist cutoff (in nm) rcoulomb= 1.0; short-range electrostatic cutoff (in nm) rvdw= 1.0; short-range van der Waals cutoff (in nm) rlistlong= 1.0; long-range neighborlist cutoff (in nm) cutoff-scheme = Verlet ; Electrostatics coulombtype= PME; Particle Mesh Ewald for long-range electrostat ics pme_order= 4; cubic interpolation fourierspacing= 0.16; grid spacing for FFT nstcomm = 10; remove com every 10 steps ; Temperature coupling is on tcoupl= V-rescale; modified Berendsen thermostat tc-grps= Protein Non-Protein; two coupling groups - more accurate tau_t= 0.10.1; time constant, in ps ref_t= 318 318; reference temperature, one for each group, in K ; Pressure coupling is off pcoupl = Parrinello-Rahman ; pressure coupling is on for NP T 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; Velocity generation is off gen_temp= 318; reference temperature, for protein in K grompp -f md.mdp -c npt.gro -p complex.top -o md.tpr? Could you please advice? Kind regards, Ankita On Fri, Apr 18, 2014 at 12:58 PM, Justin Lemkul jalem...@vt.edu wrote: On 4/18/14, 4:57 AM, Ankita Naithani wrote: Hi, I am trying to run a simulation with the following optins in my production mdp file. Previously I have run the simulations with the production run parameters using V-rescale thermostat and berendsen barostat. But now recently I read in discussion that Justin had advised using Parinello-Rahman barostat as it would be much more suitable and provides adequate sampling. Also, I had to use the Verlet cut-off scheme because the cluster where I run my production simulation has updated gomacs to the recent version 4.6.5 which requires one to use the verlet scheme. A few things to clarify: 1. Parrinello-Rahman is superior to Berendsen for production runs, not necessarily equilibration. In fact, P-R is often unstable during equilibration (see below). 2. The Verlet scheme is *not* required in version 4.6.5. It is the new default, as the group scheme will eventually be phased out, but that is not the case at the moment. My mdp file is below title= production MD ; Run parameters integrator= md; leap-frog algorithm nsteps= 5000; 0.002 * 5000 = 10 ps or 100 ns dt = 0.002 ; 2 fs ; Output control nstxout= 0; save coordinates every 2 ps nstvout= 0; save velocities every 2 ps nstxtcout= 1000; xtc compressed trajectory output every 5 ps nstenergy= 1000; save energies every 5 ps nstlog= 1000; update log file every 5 ps ; Bond parameters constraint_algorithm = lincs; holonomic constraints constraints= all-bonds; all bonds (even heavy atom-H bonds) constraine d lincs_iter= 1; accuracy of LINCS lincs_order= 4; also related to accuracy ; Neighborsearching ns_type= grid; search neighboring grid cells nstlist= 5; 25 fs rlist= 1.0; short-range neighborlist cutoff (in nm) rcoulomb= 1.0; short-range electrostatic cutoff (in nm) rvdw= 1.4; short-range van der Waals cutoff (in nm) rlistlong= 1.0; long-range neighborlist cutoff (in nm) cutoff-scheme = Verlet ; Electrostatics coulombtype= PME
Re: [gmx-users] Md simulation mdp options error
Also, just one more floating query... In presssure coupling, Previously, I have been using tau_p = 1.0 but for this simulaion, I have changed it to tau_p = 2.0 This was after seeing the mdp file of one of your tutorials. Just wanted to confirm if that's okay. -- 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.
Re: [gmx-users] MPI error in gromacs 4.6
Hi Pavan, Thank you for your response. I am trying to generate the tpr file with the following parameter; ; Neighborsearching ns_type= grid; search neighboring grid cells nstlist= 5; 25 fs rlist= 1.0; short-range neighborlist cutoff (in nm) rcoulomb= 1.0; short-range electrostatic cutoff (in nm) rvdw= 1.0; short-range van der Waals cutoff (in nm) rlistlong= 1.0; long-range neighborlist cutoff (in nm) cutoff-scheme = Verlet But, I get a warning of Unknown left-hand 'cutoff-scheme' in parameter file. On Mon, Mar 24, 2014 at 11:26 AM, Pavan Kumar kumar.pavan...@gmail.comwrote: Hello Ankita You have to just include the following line in your mdp file cutoff-scheme=Verlet And run your grompp with the modfied mdp file to generate tpr file and then mdrun. Hope this doesn't give you the same error On Mon, Mar 24, 2014 at 4:47 PM, Ankita Naithani ankitanaith...@gmail.comwrote: Hi, I am trying to run a simulation of my protein (monomer ~500 residues). I had few questions and erors regarding the same. I have previously run the simulation of the apo form of the same protein using Gromacs 4.5.5 which was available at the cluster facility I was using and also which is installed in my system. However, when I tried to run the holo form, I got error : Fatal error: 11 particles communicated to PME node 106 are more than 2/3 times the cut-off out of the domain decomposition cell of their charge group in dimension y. This usually means that your system is not well equilibrated. For more information and tips for troubleshooting, please check the GROMACS website at http://www.gromacs.org/Documentation/Errors This I figured out could be solved using a lower timestep as my previous timestep was 4fs and now I have reduced it to 3fs which should work fine now. However, after producing the tpr file for production run in my GROMACS 4.5.5, I realised that the grant for the cluster facility is over and the new clusters which I am trying to set up the same protein for support only gromacs 4.6. I am trying to run the code in these clusters and I get he following error: --- Program mdrun_mpi, VERSION 4.6.3 Source code file: /home/gromacs-4.6.3/src/kernel/runner .c, line: 824 Fatal error: OpenMP threads have been requested with cut-off scheme Group, but these are only supported with cut-off scheme Verlet For more information and tips for troubleshooting, please check the GROMACS website at http://www.gromacs.org/Documentation/Errors - 1. I wanted help with my mdp options to make it compatible. 2. Since my pevious calculations were based on gromacs 4.5.5, switching to gromacs 4.6, would that break the continuity of the run or would that bring about differences in the way the trajectories would be analysed? Below, is my mdp file title= production MD ; Run parameters integrator= md; leap-frog algorithm nsteps= ; 0.003 * = 10 ps or 100 n dt= 0.003; 3 fs ; Output control nstxout= 0; save coordinates every 2 ps nstvout= 0; save velocities every 2 ps nstxtcout= 1000; xtc compressed trajectory output every 5 ps nstenergy= 1000; save energies every 5 ps nstlog= 1000; update log file every 5 ps energygrps = Protein ATP ; Bond parameters 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; 25 fs rlist= 1.0; short-range neighborlist cutoff (in nm) rcoulomb= 1.0; short-range electrostatic cutoff (in nm) rvdw= 1.0; short-range van der Waals cutoff (in nm) rlistlong= 1.0; long-range neighborlist 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 nstcomm = 10; remove com every 10 steps ; Temperature coupling is on tcoupl= V-rescale; modified Berendsen thermostat tc-grps= Protein Non-Protein; two coupling groups - more accurate tau_t= 0.10.1; time constant, in ps ref_t= 318 318; reference temperature, one for each group, in K ; Pressure coupling is off pcoupl = berendsen; Berendsen thermostat
[gmx-users] Binning of Free Energy landscapes
Hi, I was wondering as to what exactly happens during the binning of Free Energy Landscapes as generated by g_sham? I generated Free Energy Landscapes using g_sham and my Principal Components but I am trying to understand what does bin account for and how exactly does it function? Kind regards -- Ankita Naithani -- 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.
Re: [gmx-users] Bin files of g_sham
I would like to add that maybe the bins store Frame numbers. Because I seem to have around 21000 frames and I guess I am unable to see a number higher than that. Apologies for the confusion about earlier number. So, in this case if I want a confirmation at a particular frame number, how should I go about that? Could you please help me with that? Kind regards, Ankita On Tue, Dec 17, 2013 at 7:20 AM, Ankita Naithani ankitanaith...@gmail.comwrote: Hi Llyod, Thank you for your reply. I actually wanted to extract frames corresponding to the g_sham output. The input for g_sham in my case was column 1 had time points, column 2 and column 3 had Principal component 1 and Principal component 2 respectively. So, the binning which has been done, i.e. if I say that the bin number 786 has the following few numbers like : 3876, 9870, 22002 etc. These should technically correspond to the time points right which I can then dump using trjconv too maybe. Also, I am a little confused because if the bins have all the time points, I don't seem to find the time points uptil the end of simulation in that case i.e. I am unable to see higher time scale. My simulation was for 48ns and had 3fs time step. I am unable to trace anything beyond 22000 or 23000 number in my bins. I hope you could help me solve this confusion. Kind regards, Ankita On Sun, Dec 15, 2013 at 11:17 PM, lloyd riggs lloyd.ri...@gmx.ch wrote: I will attempt an answer, however there might be a better response. It depends on your input .xvg (what the first column is) and what you are telling it to do (defualt is time), and if these correspond to your time frames (can be timepoint or siply frame # I believe). g_sham seems (to me) to work like a comand line spread sheet manipulator with some built in functions. If you want a particular frame from a trajectory time point, then use g_traj with the first/last set to the frame you want but there's also 2 other tools that do this the same way... Stephan Watkins *Gesendet:* Sonntag, 15. Dezember 2013 um 14:49 Uhr *Von:* Ankita Naithani ankitanaith...@gmail.com *An:* Discussion list for GROMACS users gmx-us...@gromacs.org *Betreff:* [gmx-users] Bin files of g_sham Hi, I had a query regarding the bin files generated by g_sham. Does the bin index contain timepoints? i.e. there are 1024 bins so each bin has a number of time points? And supposes bin 756 has numbers like 3876, 20020, 7896 so these would all correspond to time and we could then dump the representative frame corresponding to these time points? I am a little confused as to what is stored in the bins? Timepoints or frame number and then alternatively how to get the representative frame/snapshot stored at that specific bin number. Kind regards, -- Ankita Naithani -- 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. -- Ankita Naithani -- Ankita Naithani -- 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.
[gmx-users] Bin files of g_sham
Hi, I had a query regarding the bin files generated by g_sham. Does the bin index contain timepoints? i.e. there are 1024 bins so each bin has a number of time points? And supposes bin 756 has numbers like 3876, 20020, 7896 so these would all correspond to time and we could then dump the representative frame corresponding to these time points? I am a little confused as to what is stored in the bins? Timepoints or frame number and then alternatively how to get the representative frame/snapshot stored at that specific bin number. Kind regards, -- Ankita Naithani -- 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.
[gmx-users] Fwd: Help on the vector component
Hi, I wanted some help regarding vector components. I.e. when we use the g_anaeig command with a flag of -comp, it outputs a file of the corresponding eigenvector per atom. So, the x-axis has the atom number and the Y-axis has the total, x, y and z components. I am not sure about what do they mean? Also, what are the units of these components? Like the numbers on Y-axis, what are the corresponding units? I would be really grateful if anyone could kindly help me on this. Apologies for reposting but was wondering if anyone has any thoughts on it? Best wishes, -- Ankita Naithani -- Ankita Naithani -- 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.
