Re: [gmx-users] Huge acceleration needed to reproduce results!
Hi, Yes I realised that gromacs works in ps. I converted my force in kj mol-1 A-1 to acceleration in nm/ps2. I also took into account that the msd.xvg is plotted in nm and ps-2 and the calculated gradient printed at the top of the msd.xvg file is in cm2/s. One strange thing that I do get is the message ?There were 228 inconsistent shifts. Check your topology? when I carry out the g_msd with the ?mol option but not when I don?t use -mol. Why is this? I also came across a forum post (http://www.mail-archive.com/gmx-users@gromacs.org/msg5.html) that said ?If the distance between two atoms is close to half the box, the force may arbitrarily change sign. This is an ill-defined situation for which there is no obvious solution.? Could this somehow be affecting my simulations? Below are the relevant parts of my .mdp file and other files. If you see something suspicious please let me know because I?m stuck, Thanks ; RUN CONTROL PARAMETERS integrator = md ; Start time and timestep in ps tinit= 0 dt = 0.001 nsteps = 100 ; 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= ; OUTPUT CONTROL OPTIONS ; Output frequency for coords (x), velocities (v) and forces (f) nstxout = 1000 nstvout = 1000 nstfout = 0 ; Output frequency for energies to log file and energy file nstlog = 1000 nstenergy= 1000 ; Output frequency and precision for xtc file nstxtcout= 1000 xtc-precision= 1000 ; This selects the subset of atoms for the xtc file. You can ; select multiple groups. By default all atoms will be written. xtc-grps = ; Selection of energy groups energygrps = ; NEIGHBORSEARCHING PARAMETERS ; nblist update frequency nstlist = ; ns algorithm (simple or grid) ns_type = grid ; Periodic boundary conditions: xyz, no, xy pbc = xyz periodic_molecules = yes ; nblist cut-off rlist= 0.9 ; OPTIONS FOR ELECTROSTATICS AND VDW ; Method for doing electrostatics coulombtype = PME rcoulomb-switch = 0 rcoulomb = 0.9 ; Relative dielectric constant for the medium and the reaction field epsilon_r= epsilon_rf = ; Method for doing Van der Waals vdw-type = Cut-off ; cut-off lengths rvdw-switch = 0 rvdw = 0.9 ; Apply long range dispersion corrections for Energy and Pressure DispCorr = No ; Extension of the potential lookup tables beyond the cut-off table-extension = ; 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= ewald_rtol = 1e-05 ewald_geometry = 3d epsilon_surface = 0 optimize_fft = yes ; OPTIONS FOR WEAK COUPLING ALGORITHMS ; Temperature coupling tcoupl = nose-hoover ; Groups to couple separately tc-grps = System ; Time constant (ps) and reference temperature (K) tau_t= 0.1 ref_t= 150 ; Pressure coupling Pcoupl = No Pcoupltype = ; Time constant (ps), compressibility (1/bar) and reference P (bar) tau-p= compressibility = ref-p= ; Scaling of reference coordinates, No, All or COM refcoord_scaling = no ; Random seed for Andersen thermostat andersen_seed= ; GENERATE VELOCITIES FOR STARTUP RUN gen_vel = yes gen_temp = 150 gen_seed = 173529 ; OPTIONS FOR BONDS constraints = none ; 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 = 4 ; 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
Re: [gmx-users] Huge acceleration needed to reproduce results!
Title indicates you think that you have CH4 in MCM: [ system ] ; Name CH4 in MCM Actual system is MCM in CH4: [ molecules ] ; Compound#mols MCM1 CH4340 .mdp accelerates the CH4: acc-grps = CH4 Now I'm not sure if this is the source of any problem, equal and opposite being true, but you are certainly trying to accelerate the more massive group and it seems like a strange thing to do. Could this be it? Chris. -- original message -- Quoting Jennifer Williams jennifer.willi...@ed.ac.uk: Hi, Yes I realised that gromacs works in ps. I converted my force in kj mol-1 A-1 to acceleration in nm/ps2. I also took into account that the msd.xvg is plotted in nm and ps-2 and the calculated gradient printed at the top of the msd.xvg file is in cm2/s. One strange thing that I do get is the message ?There were 228 inconsistent shifts. Check your topology? when I carry out the g_msd with the ?mol option but not when I don?t use -mol. Why is this? I also came across a forum post (http://www.mail-archive.com/gmx-users@gromacs.org/msg5.html) that said ?If the distance between two atoms is close to half the box, the force may arbitrarily change sign. This is an ill-defined situation for which there is no obvious solution.? Could this somehow be affecting my simulations? Below are the relevant parts of my .mdp file and other files. If you see something suspicious please let me know because I?m stuck, Thanks ; RUN CONTROL PARAMETERS integrator = md ; Start time and timestep in ps tinit= 0 dt = 0.001 nsteps = 100 ; 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= ; OUTPUT CONTROL OPTIONS ; Output frequency for coords (x), velocities (v) and forces (f) nstxout = 1000 nstvout = 1000 nstfout = 0 ; Output frequency for energies to log file and energy file nstlog = 1000 nstenergy= 1000 ; Output frequency and precision for xtc file nstxtcout= 1000 xtc-precision= 1000 ; This selects the subset of atoms for the xtc file. You can ; select multiple groups. By default all atoms will be written. xtc-grps = ; Selection of energy groups energygrps = ; NEIGHBORSEARCHING PARAMETERS ; nblist update frequency nstlist = ; ns algorithm (simple or grid) ns_type = grid ; Periodic boundary conditions: xyz, no, xy pbc = xyz periodic_molecules = yes ; nblist cut-off rlist= 0.9 ; OPTIONS FOR ELECTROSTATICS AND VDW ; Method for doing electrostatics coulombtype = PME rcoulomb-switch = 0 rcoulomb = 0.9 ; Relative dielectric constant for the medium and the reaction field epsilon_r= epsilon_rf = ; Method for doing Van der Waals vdw-type = Cut-off ; cut-off lengths rvdw-switch = 0 rvdw = 0.9 ; Apply long range dispersion corrections for Energy and Pressure DispCorr = No ; Extension of the potential lookup tables beyond the cut-off table-extension = ; 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= ewald_rtol = 1e-05 ewald_geometry = 3d epsilon_surface = 0 optimize_fft = yes ; OPTIONS FOR WEAK COUPLING ALGORITHMS ; Temperature coupling tcoupl = nose-hoover ; Groups to couple separately tc-grps = System ; Time constant (ps) and reference temperature (K) tau_t= 0.1 ref_t= 150 ; Pressure coupling Pcoupl = No Pcoupltype = ; Time constant (ps), compressibility (1/bar) and reference P (bar) tau-p= compressibility = ref-p= ; Scaling of reference coordinates, No, All or COM refcoord_scaling = no ; Random seed for Andersen thermostat andersen_seed= ; GENERATE VELOCITIES FOR STARTUP RUN gen_vel = yes gen_temp = 150 gen_seed = 173529 ; OPTIONS FOR BONDS constraints = none ; Type of constraint algorithm
Re: [gmx-users] Huge acceleration needed to reproduce results!
Hi Chris, Thanks for the input-its good to get another person's perspective. Let me explain a bit more about my system... MCM stands for MCM-41, a mesoporous silica which consits of 1071 atoms so even though it is only 1 molecule, it is the more massive group. I am purposely accelerating the methane to try and calculate their diffusion thourgh the pore of the MCM. As I understand it, the acceleration I specify should be applied to each and every molecule of CH4. I.e if I specify an acceleration of 0.1 nm ps-2, each of my 340 molecules should experience an acceleration of 0.1. I.e the acceleration is not in any way divided up between the total number of accelerated molecules? Jenny Quoting chris.ne...@utoronto.ca: Title indicates you think that you have CH4 in MCM: [ system ] ; Name CH4 in MCM Actual system is MCM in CH4: [ molecules ] ; Compound#mols MCM1 CH4340 .mdp accelerates the CH4: acc-grps = CH4 Now I'm not sure if this is the source of any problem, equal and opposite being true, but you are certainly trying to accelerate the more massive group and it seems like a strange thing to do. Could this be it? Chris. -- original message -- Quoting Jennifer Williams jennifer.willi...@ed.ac.uk: Hi, Yes I realised that gromacs works in ps. I converted my force in kj mol-1 A-1 to acceleration in nm/ps2. I also took into account that the msd.xvg is plotted in nm and ps-2 and the calculated gradient printed at the top of the msd.xvg file is in cm2/s. One strange thing that I do get is the message ?There were 228 inconsistent shifts. Check your topology? when I carry out the g_msd with the ?mol option but not when I don?t use -mol. Why is this? I also came across a forum post (http://www.mail-archive.com/gmx-users@gromacs.org/msg5.html) that said ?If the distance between two atoms is close to half the box, the force may arbitrarily change sign. This is an ill-defined situation for which there is no obvious solution.? Could this somehow be affecting my simulations? Below are the relevant parts of my .mdp file and other files. If you see something suspicious please let me know because I?m stuck, Thanks ; RUN CONTROL PARAMETERS integrator = md ; Start time and timestep in ps tinit= 0 dt = 0.001 nsteps = 100 ; 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= ; OUTPUT CONTROL OPTIONS ; Output frequency for coords (x), velocities (v) and forces (f) nstxout = 1000 nstvout = 1000 nstfout = 0 ; Output frequency for energies to log file and energy file nstlog = 1000 nstenergy= 1000 ; Output frequency and precision for xtc file nstxtcout= 1000 xtc-precision= 1000 ; This selects the subset of atoms for the xtc file. You can ; select multiple groups. By default all atoms will be written. xtc-grps = ; Selection of energy groups energygrps = ; NEIGHBORSEARCHING PARAMETERS ; nblist update frequency nstlist = ; ns algorithm (simple or grid) ns_type = grid ; Periodic boundary conditions: xyz, no, xy pbc = xyz periodic_molecules = yes ; nblist cut-off rlist= 0.9 ; OPTIONS FOR ELECTROSTATICS AND VDW ; Method for doing electrostatics coulombtype = PME rcoulomb-switch = 0 rcoulomb = 0.9 ; Relative dielectric constant for the medium and the reaction field epsilon_r= epsilon_rf = ; Method for doing Van der Waals vdw-type = Cut-off ; cut-off lengths rvdw-switch = 0 rvdw = 0.9 ; Apply long range dispersion corrections for Energy and Pressure DispCorr = No ; Extension of the potential lookup tables beyond the cut-off table-extension = ; 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= ewald_rtol = 1e-05 ewald_geometry = 3d epsilon_surface = 0 optimize_fft = yes ; OPTIONS FOR WEAK COUPLING ALGORITHMS ; Temperature coupling tcoupl = nose-hoover ; Groups to couple
Re: [gmx-users] Huge acceleration needed to reproduce results!
As for your posted problem, I don't think that I can be of any more assistance. I picked the low-hanging fruit, but probably you need some assistance from somebody who has actually done constant acceleration simulations. I will take this opportunity to note, though, that you will never get diffusion information from a constant acceleration simulation (as far as I can figure). Although I am certain that you require a resolution to your posted problem, I suggest that you concurrently ensure that your approach is even warranted. Chris. Quoting Jennifer Williams jennifer.willi...@ed.ac.uk: Hi Chris, Thanks for the input-its good to get another person's perspective. Let me explain a bit more about my system... MCM stands for MCM-41, a mesoporous silica which consits of 1071 atoms so even though it is only 1 molecule, it is the more massive group. I am purposely accelerating the methane to try and calculate their diffusion thourgh the pore of the MCM. As I understand it, the acceleration I specify should be applied to each and every molecule of CH4. I.e if I specify an acceleration of 0.1 nm ps-2, each of my 340 molecules should experience an acceleration of 0.1. I.e the acceleration is not in any way divided up between the total number of accelerated molecules? Jenny Quoting chris.ne...@utoronto.ca: Title indicates you think that you have CH4 in MCM: [ system ] ; Name CH4 in MCM Actual system is MCM in CH4: [ molecules ] ; Compound#mols MCM1 CH4340 .mdp accelerates the CH4: acc-grps = CH4 Now I'm not sure if this is the source of any problem, equal and opposite being true, but you are certainly trying to accelerate the more massive group and it seems like a strange thing to do. Could this be it? Chris. -- original message -- Quoting Jennifer Williams jennifer.willi...@ed.ac.uk: Hi, Yes I realised that gromacs works in ps. I converted my force in kj mol-1 A-1 to acceleration in nm/ps2. I also took into account that the msd.xvg is plotted in nm and ps-2 and the calculated gradient printed at the top of the msd.xvg file is in cm2/s. One strange thing that I do get is the message ?There were 228 inconsistent shifts. Check your topology? when I carry out the g_msd with the ?mol option but not when I don?t use -mol. Why is this? I also came across a forum post (http://www.mail-archive.com/gmx-users@gromacs.org/msg5.html) that said ?If the distance between two atoms is close to half the box, the force may arbitrarily change sign. This is an ill-defined situation for which there is no obvious solution.? Could this somehow be affecting my simulations? Below are the relevant parts of my .mdp file and other files. If you see something suspicious please let me know because I?m stuck, Thanks ; RUN CONTROL PARAMETERS integrator = md ; Start time and timestep in ps tinit= 0 dt = 0.001 nsteps = 100 ; 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= ; OUTPUT CONTROL OPTIONS ; Output frequency for coords (x), velocities (v) and forces (f) nstxout = 1000 nstvout = 1000 nstfout = 0 ; Output frequency for energies to log file and energy file nstlog = 1000 nstenergy= 1000 ; Output frequency and precision for xtc file nstxtcout= 1000 xtc-precision= 1000 ; This selects the subset of atoms for the xtc file. You can ; select multiple groups. By default all atoms will be written. xtc-grps = ; Selection of energy groups energygrps = ; NEIGHBORSEARCHING PARAMETERS ; nblist update frequency nstlist = ; ns algorithm (simple or grid) ns_type = grid ; Periodic boundary conditions: xyz, no, xy pbc = xyz periodic_molecules = yes ; nblist cut-off rlist= 0.9 ; OPTIONS FOR ELECTROSTATICS AND VDW ; Method for doing electrostatics coulombtype = PME rcoulomb-switch = 0 rcoulomb = 0.9 ; Relative dielectric constant for the medium and the reaction field epsilon_r= epsilon_rf = ; Method for doing Van der Waals vdw-type = Cut-off ; cut-off lengths rvdw-switch = 0 rvdw = 0.9 ; Apply long range dispersion corrections for Energy and Pressure DispCorr = No ; Extension of the potential lookup tables beyond the cut-off table-extension = ;