Re: [gmx-users] Re: Energy Drift in Gromacs 4.0
Ilya Chorny wrote: Ok, I got it. Set electrostatics to PME-switch. Made rcoulomb-switch .1 less than rcoulomb and set rlist accordingly and now I am conserving energy. Yipee!!! It did not make sense to me that the system had to be fully relaxed to conserve energy. It should conserve energy energy as long as it has been minimized and has no major clashes. This is correct, provided T and P coupling are turned off. Thanks, Ilya On Wed, Mar 11, 2009 at 7:38 PM, Ilya Chorny icho...@gmail.com mailto:icho...@gmail.com wrote: If I plot the K.E(t) and P.E.(t) they are both going down. In an NVE simulation they should counter balance. I am leaking K.E. somewhere. I reran with a relaxed P.E(t) and saw the same behavior. This is happening both at 2 and 4 fs so its not a time step issue. Any thoughts? I thought it might be COM removal so I turned it off to no avail. Any chance one of the Gromacs people can give me the mdp file from the DHFR simulation in the 4.0 paper. On Tue, Mar 10, 2009 at 10:21 PM, Mark Abraham mark.abra...@anu.edu.au mailto:mark.abra...@anu.edu.au wrote: Ilya Chorny wrote: Forgot column labels on energies. Average RMSD Fluct. Drift Tot-Drift 2 fs time step Total Energy46497.2 36445.6 2630.7 -148.833-148833 4 fs time step Total Energy41580.8 37693.42856.08 -130.198-130198 On Tue, Mar 10, 2009 at 7:17 PM, Ilya Chorny icho...@gmail.com mailto:icho...@gmail.com mailto:icho...@gmail.com mailto:icho...@gmail.com wrote: Hello All, I am trying to run some calibration calculations with 2/4 fs time steps. I am trying to reproduce the results in the Gromacs 4.0 paper on a protein/water (not the same protein) system with ~100K atoms. I ran 1 ns simulation in the NVE ensemble. My mdp params are shown below. Topolgy was created using -vsite h. My system has 218441 degrees of freedom according to GMXDUMP. My P.E. energy is ~ -1E5 and has not full relaxed after 1ns but is close. 2 fs time step Total Energy46497.2 36445.62630.7 -148.833-148833 4 fs time step Total Energy41580.8 37693.4 2856.08 -130.198-130198 Drift (kt/ns) 2fs/4fs = .34/.29 The drift seems to be much larger for both the 2fs and 4fs time steps then in the paper. Do I have a clear bug in my params? Should I wait till the system fully relaxes before doing the drift calculation? You should definitely equilibrate the system before trying to measure energy drift. I suggest doing so in NPT to fix your density and temperature, then switching to NVE, wait a bit, and only then start to collect your data. coulombtype = PME rcoulomb = 1.0 fourierspacing = 0.12 ; FFT grid size, when a value is 0 fourierspacing will be used fourier_nx = 0 fourier_ny = 0 fourier_nz = 0 ; EWALD/PME/PPPM parameters pme_order= 4 ewald_rtol = 1.e-04 I would think the ewald_rtol value used in the paper was 1e-5, but it is unfortunately not stated. The different numbers will lead to different costs and accuracies for the direct and FFT parts of the PME algorithm. rcoulomb, fourier_whatever, pme_order and ewald_rtol are all relevant if you're trying to reproduce results. Mark ___ gmx-users mailing listgmx-users@gromacs.org mailto:gmx-users@gromacs.org http://www.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/search before posting! Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org mailto:gmx-users-requ...@gromacs.org. Can't post? Read http://www.gromacs.org/mailing_lists/users.php -- Ilya Chorny Ph.D. -- Ilya Chorny Ph.D. ___ gmx-users mailing listgmx-users@gromacs.org http://www.gromacs.org/mailman/listinfo/gmx-users
[gmx-users] Re: Energy Drift in Gromacs 4.0
Even if a system is not well equilibrated (but not so far from equilibrium) ,energy should be conserved in NVE ensemble. Try to lower the nstlist from 10 to 1 i.e. renew every time the neighbourlist. In systems I have studied, that is the key that improve much more energy conservation. It works even without switching of VdW potential. Luca Dre Gaetani Message: 1 Date: Thu, 12 Mar 2009 13:59:11 +1100 From: Mark Abraham mark.abra...@anu.edu.au Subject: Re: [gmx-users] Re: Energy Drift in Gromacs 4.0 To: Discussion list for GROMACS users gmx-users@gromacs.org Message-ID: 49b87a7f.30...@anu.edu.au Content-Type: text/plain; charset=ISO-8859-1; format=flowed Ilya Chorny wrote: If I plot the K.E(t) and P.E.(t) they are both going down. In an NVE simulation they should counter balance. I am leaking K.E. somewhere. I reran with a relaxed P.E(t) and saw the same behavior. This is happening both at 2 and 4 fs so its not a time step issue. Any thoughts? I thought it might be COM removal so I turned it off to no avail. going down is seriously unhelpful. I told you you needed to equilibrate before measuring such things. If you want people to take your problem seriously, you'll need to announce that you have equilibrated and for how long and in what ensembles, etc. If you're having a problem, you need to demonstrate what areas you feel are not causing it. We're not your family doctor who will ask a pile of questions to find out the problem :-) Mark ___ gmx-users mailing listgmx-users@gromacs.org http://www.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/search before posting! Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. Can't post? Read http://www.gromacs.org/mailing_lists/users.php
RE: [gmx-users] Re: Energy Drift in Gromacs 4.0
Hi, Everything depends completely on how much drift you are talking about. There will always be energy drift, it is just a question how much. In general shifted potentials will be cheaper (and lead to less energy drift) than setting nstlist to 1, but for most purposes a bit of energy drift is no problem. Berk Date: Thu, 12 Mar 2009 10:55:47 +0100 From: degaet...@dcci.unipi.it To: gmx-users@gromacs.org Subject: [gmx-users] Re: Energy Drift in Gromacs 4.0 Even if a system is not well equilibrated (but not so far from equilibrium) ,energy should be conserved in NVE ensemble. Try to lower the nstlist from 10 to 1 i.e. renew every time the neighbourlist. In systems I have studied, that is the key that improve much more energy conservation. It works even without switching of VdW potential. Luca Dre Gaetani Message: 1 Date: Thu, 12 Mar 2009 13:59:11 +1100 From: Mark Abraham mark.abra...@anu.edu.au Subject: Re: [gmx-users] Re: Energy Drift in Gromacs 4.0 To: Discussion list for GROMACS users gmx-users@gromacs.org Message-ID: 49b87a7f.30...@anu.edu.au Content-Type: text/plain; charset=ISO-8859-1; format=flowed Ilya Chorny wrote: If I plot the K.E(t) and P.E.(t) they are both going down. In an NVE simulation they should counter balance. I am leaking K.E. somewhere. I reran with a relaxed P.E(t) and saw the same behavior. This is happening both at 2 and 4 fs so its not a time step issue. Any thoughts? I thought it might be COM removal so I turned it off to no avail. going down is seriously unhelpful. I told you you needed to equilibrate before measuring such things. If you want people to take your problem seriously, you'll need to announce that you have equilibrated and for how long and in what ensembles, etc. If you're having a problem, you need to demonstrate what areas you feel are not causing it. We're not your family doctor who will ask a pile of questions to find out the problem :-) Mark ___ gmx-users mailing listgmx-users@gromacs.org http://www.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/search before posting! Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. Can't post? Read http://www.gromacs.org/mailing_lists/users.php _ Express yourself instantly with MSN Messenger! Download today it's FREE! http://messenger.msn.click-url.com/go/onm00200471ave/direct/01/___ gmx-users mailing listgmx-users@gromacs.org http://www.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/search before posting! Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. Can't post? Read http://www.gromacs.org/mailing_lists/users.php
Re: [gmx-users] Re: Energy Drift in Gromacs 4.0
If I plot the K.E(t) and P.E.(t) they are both going down. In an NVE simulation they should counter balance. I am leaking K.E. somewhere. I reran with a relaxed P.E(t) and saw the same behavior. This is happening both at 2 and 4 fs so its not a time step issue. Any thoughts? I thought it might be COM removal so I turned it off to no avail. Any chance one of the Gromacs people can give me the mdp file from the DHFR simulation in the 4.0 paper. On Tue, Mar 10, 2009 at 10:21 PM, Mark Abraham mark.abra...@anu.edu.auwrote: Ilya Chorny wrote: Forgot column labels on energies. Average RMSD Fluct. Drift Tot-Drift 2 fs time step Total Energy46497.236445.6 2630.7 -148.833-148833 4 fs time step Total Energy41580.837693.42856.08 -130.198-130198 On Tue, Mar 10, 2009 at 7:17 PM, Ilya Chorny icho...@gmail.com mailto: icho...@gmail.com wrote: Hello All, I am trying to run some calibration calculations with 2/4 fs time steps. I am trying to reproduce the results in the Gromacs 4.0 paper on a protein/water (not the same protein) system with ~100K atoms. I ran 1 ns simulation in the NVE ensemble. My mdp params are shown below. Topolgy was created using -vsite h. My system has 218441 degrees of freedom according to GMXDUMP. My P.E. energy is ~ -1E5 and has not full relaxed after 1ns but is close. 2 fs time step Total Energy46497.236445.6 2630.7 -148.833-148833 4 fs time step Total Energy41580.837693.4 2856.08 -130.198-130198 Drift (kt/ns) 2fs/4fs = .34/.29 The drift seems to be much larger for both the 2fs and 4fs time steps then in the paper. Do I have a clear bug in my params? Should I wait till the system fully relaxes before doing the drift calculation? You should definitely equilibrate the system before trying to measure energy drift. I suggest doing so in NPT to fix your density and temperature, then switching to NVE, wait a bit, and only then start to collect your data. coulombtype = PME rcoulomb = 1.0 fourierspacing = 0.12 ; FFT grid size, when a value is 0 fourierspacing will be used fourier_nx = 0 fourier_ny = 0 fourier_nz = 0 ; EWALD/PME/PPPM parameters pme_order= 4 ewald_rtol = 1.e-04 I would think the ewald_rtol value used in the paper was 1e-5, but it is unfortunately not stated. The different numbers will lead to different costs and accuracies for the direct and FFT parts of the PME algorithm. rcoulomb, fourier_whatever, pme_order and ewald_rtol are all relevant if you're trying to reproduce results. Mark ___ gmx-users mailing listgmx-users@gromacs.org http://www.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/search before posting! Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. Can't post? Read http://www.gromacs.org/mailing_lists/users.php -- Ilya Chorny Ph.D. ___ gmx-users mailing listgmx-users@gromacs.org http://www.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/search before posting! Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. Can't post? Read http://www.gromacs.org/mailing_lists/users.php
Re: [gmx-users] Re: Energy Drift in Gromacs 4.0
Ok, I got it. Set electrostatics to PME-switch. Made rcoulomb-switch .1 less than rcoulomb and set rlist accordingly and now I am conserving energy. Yipee!!! It did not make sense to me that the system had to be fully relaxed to conserve energy. It should conserve energy energy as long as it has been minimized and has no major clashes. Thanks, Ilya On Wed, Mar 11, 2009 at 7:38 PM, Ilya Chorny icho...@gmail.com wrote: If I plot the K.E(t) and P.E.(t) they are both going down. In an NVE simulation they should counter balance. I am leaking K.E. somewhere. I reran with a relaxed P.E(t) and saw the same behavior. This is happening both at 2 and 4 fs so its not a time step issue. Any thoughts? I thought it might be COM removal so I turned it off to no avail. Any chance one of the Gromacs people can give me the mdp file from the DHFR simulation in the 4.0 paper. On Tue, Mar 10, 2009 at 10:21 PM, Mark Abraham mark.abra...@anu.edu.auwrote: Ilya Chorny wrote: Forgot column labels on energies. Average RMSD Fluct. Drift Tot-Drift 2 fs time step Total Energy46497.236445.6 2630.7 -148.833-148833 4 fs time step Total Energy41580.837693.42856.08 -130.198-130198 On Tue, Mar 10, 2009 at 7:17 PM, Ilya Chorny icho...@gmail.com mailto: icho...@gmail.com wrote: Hello All, I am trying to run some calibration calculations with 2/4 fs time steps. I am trying to reproduce the results in the Gromacs 4.0 paper on a protein/water (not the same protein) system with ~100K atoms. I ran 1 ns simulation in the NVE ensemble. My mdp params are shown below. Topolgy was created using -vsite h. My system has 218441 degrees of freedom according to GMXDUMP. My P.E. energy is ~ -1E5 and has not full relaxed after 1ns but is close. 2 fs time step Total Energy46497.236445.6 2630.7 -148.833-148833 4 fs time step Total Energy41580.837693.4 2856.08 -130.198-130198 Drift (kt/ns) 2fs/4fs = .34/.29 The drift seems to be much larger for both the 2fs and 4fs time steps then in the paper. Do I have a clear bug in my params? Should I wait till the system fully relaxes before doing the drift calculation? You should definitely equilibrate the system before trying to measure energy drift. I suggest doing so in NPT to fix your density and temperature, then switching to NVE, wait a bit, and only then start to collect your data. coulombtype = PME rcoulomb = 1.0 fourierspacing = 0.12 ; FFT grid size, when a value is 0 fourierspacing will be used fourier_nx = 0 fourier_ny = 0 fourier_nz = 0 ; EWALD/PME/PPPM parameters pme_order= 4 ewald_rtol = 1.e-04 I would think the ewald_rtol value used in the paper was 1e-5, but it is unfortunately not stated. The different numbers will lead to different costs and accuracies for the direct and FFT parts of the PME algorithm. rcoulomb, fourier_whatever, pme_order and ewald_rtol are all relevant if you're trying to reproduce results. Mark ___ gmx-users mailing listgmx-users@gromacs.org http://www.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/search before posting! Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. Can't post? Read http://www.gromacs.org/mailing_lists/users.php -- Ilya Chorny Ph.D. -- Ilya Chorny Ph.D. ___ gmx-users mailing listgmx-users@gromacs.org http://www.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/search before posting! Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. Can't post? Read http://www.gromacs.org/mailing_lists/users.php
Re: [gmx-users] Re: Energy Drift in Gromacs 4.0
Ilya Chorny wrote: If I plot the K.E(t) and P.E.(t) they are both going down. In an NVE simulation they should counter balance. I am leaking K.E. somewhere. I reran with a relaxed P.E(t) and saw the same behavior. This is happening both at 2 and 4 fs so its not a time step issue. Any thoughts? I thought it might be COM removal so I turned it off to no avail. going down is seriously unhelpful. I told you you needed to equilibrate before measuring such things. If you want people to take your problem seriously, you'll need to announce that you have equilibrated and for how long and in what ensembles, etc. If you're having a problem, you need to demonstrate what areas you feel are not causing it. We're not your family doctor who will ask a pile of questions to find out the problem :-) Mark ___ gmx-users mailing listgmx-users@gromacs.org http://www.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/search before posting! Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. Can't post? Read http://www.gromacs.org/mailing_lists/users.php
[gmx-users] Re: Energy Drift in Gromacs 4.0
Forgot column labels on energies. Average RMSD Fluct. Drift Tot-Drift2 fs time step Total Energy 46497.236445.6 2630.7 -148.833-148833 4 fs time step Total Energy41580.837693.42856.08 -130.198-130198 On Tue, Mar 10, 2009 at 7:17 PM, Ilya Chorny icho...@gmail.com wrote: Hello All, I am trying to run some calibration calculations with 2/4 fs time steps. I am trying to reproduce the results in the Gromacs 4.0 paper on a protein/water (not the same protein) system with ~100K atoms. I ran 1 ns simulation in the NVE ensemble. My mdp params are shown below. Topolgy was created using -vsite h. My system has 218441 degrees of freedom according to GMXDUMP. My P.E. energy is ~ -1E5 and has not full relaxed after 1ns but is close. 2 fs time step Total Energy46497.236445.6 2630.7 -148.833-148833 4 fs time step Total Energy41580.837693.42856.08 -130.198-130198 Drift (kt/ns) 2fs/4fs = .34/.29 The drift seems to be much larger for both the 2fs and 4fs time steps then in the paper. Do I have a clear bug in my params? Should I wait till the system fully relaxes before doing the drift calculation? integrator = md ; Start time and timestep in ps tinit= 0 dt = 0.002/0.004 nsteps = 50/25 ; For exact run continuation or redoing part of a run ; Part index is updated automatically on checkpointing (keeps files separate) simulation_part = 1 init_step= 0 ; mode for center of mass motion removal comm_mode= linear ; number of steps for center of mass motion removal nstcomm = 1 ; group(s) for center of mass motion removal comm_grps= system ; NEIGHBORSEARCHING PARAMETERS ; nblist update frequency nstlist = 10/5 ; ns algorithm (simple or grid) ns_type = grid ; Periodic boundary conditions: xyz, no, xy pbc = xyz periodic_molecules = no ; nblist cut-off rlist= 1.0 ; OPTIONS FOR ELECTROSTATICS AND VDW ; Method for doing electrostatics coulombtype = PME rcoulomb-switch = 0 rcoulomb = 1.0 ; Relative dielectric constant for the medium and the reaction field epsilon-r= 80 epsilon_rf = 1 ; Method for doing Van der Waals vdw-type = Switch ; cut-off lengths rvdw-switch = .9 rvdw = 1.0 ; Apply long range dispersion corrections for Energy and Pressure DispCorr = EnerPres ; Extension of the potential lookup tables beyond the cut-off table-extension = 1 ; Seperate tables between energy group pairs energygrp_table = ; Spacing for the PME/PPPM FFT grid fourierspacing = 0.12 ; FFT grid size, when a value is 0 fourierspacing will be used fourier_nx = 0 fourier_ny = 0 fourier_nz = 0 ; EWALD/PME/PPPM parameters pme_order= 4 ewald_rtol = 1.e-04 ewald_geometry = 3d epsilon_surface = 0 optimize_fft = no ; OPTIONS FOR WEAK COUPLING ALGORITHMS ; Temperature coupling Tcoupl = no ; Groups to couple separately tc-grps = System ; Time constant (ps) and reference temperature (K) tau_t= 0.1 ref_t= 298.0 ; Pressure coupling Pcoupl = No Pcoupltype = Isotropic ; Time constant (ps), compressibility (1/bar) and reference P (bar) tau_p= 1.67 compressibility = 4.5e-5 ref_p= 1.01325 ; Scaling of reference coordinates, No, All or COM refcoord_scaling = No ; Random seed for Andersen thermostat andersen_seed= 815131 ; GENERATE VELOCITIES FOR STARTUP RUN gen_vel = yes gen_temp = 298.0 gen-seed = 173529 ; OPTIONS FOR BONDS constraints = hbonds/all-bonds ; Type of constraint algorithm constraint-algorithm = lincs ; Do not constrain the start configuration continuation = no ; Use successive overrelaxation to reduce the number of shake iterations Shake-SOR= no ; Relative tolerance of shake shake-tol= 0.0001 ; Highest order in the expansion of the constraint coupling matrix lincs-order = 6 ; Number of iterations in the final step of LINCS. 1 is fine for ; normal simulations, but use 2 to conserve energy in NVE runs. ; For energy minimization with constraints it should be 4 to 8. lincs-iter = 2 ; Lincs will write a warning to the stderr if in one step a bond ; rotates over more degrees than
Re: [gmx-users] Re: Energy Drift in Gromacs 4.0
Ilya Chorny wrote: Forgot column labels on energies. Average RMSD Fluct. Drift Tot-Drift 2 fs time step Total Energy46497.236445.6 2630.7 -148.833-148833 4 fs time step Total Energy41580.837693.42856.08 -130.198-130198 On Tue, Mar 10, 2009 at 7:17 PM, Ilya Chorny icho...@gmail.com mailto:icho...@gmail.com wrote: Hello All, I am trying to run some calibration calculations with 2/4 fs time steps. I am trying to reproduce the results in the Gromacs 4.0 paper on a protein/water (not the same protein) system with ~100K atoms. I ran 1 ns simulation in the NVE ensemble. My mdp params are shown below. Topolgy was created using -vsite h. My system has 218441 degrees of freedom according to GMXDUMP. My P.E. energy is ~ -1E5 and has not full relaxed after 1ns but is close. 2 fs time step Total Energy46497.236445.6 2630.7 -148.833-148833 4 fs time step Total Energy41580.837693.4 2856.08 -130.198-130198 Drift (kt/ns) 2fs/4fs = .34/.29 The drift seems to be much larger for both the 2fs and 4fs time steps then in the paper. Do I have a clear bug in my params? Should I wait till the system fully relaxes before doing the drift calculation? You should definitely equilibrate the system before trying to measure energy drift. I suggest doing so in NPT to fix your density and temperature, then switching to NVE, wait a bit, and only then start to collect your data. coulombtype = PME rcoulomb = 1.0 fourierspacing = 0.12 ; FFT grid size, when a value is 0 fourierspacing will be used fourier_nx = 0 fourier_ny = 0 fourier_nz = 0 ; EWALD/PME/PPPM parameters pme_order= 4 ewald_rtol = 1.e-04 I would think the ewald_rtol value used in the paper was 1e-5, but it is unfortunately not stated. The different numbers will lead to different costs and accuracies for the direct and FFT parts of the PME algorithm. rcoulomb, fourier_whatever, pme_order and ewald_rtol are all relevant if you're trying to reproduce results. Mark ___ gmx-users mailing listgmx-users@gromacs.org http://www.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/search before posting! Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. Can't post? Read http://www.gromacs.org/mailing_lists/users.php