[gmx-users] Re: Cross-correlation maps
up :) It's appeared two additional questions. 1) In addition to the pca's cross-correlation maps I wounder to know about possibility of calculation of such cross-correlation's from the trajectories indirectly without calculation of the covariance matrices. 2) is there any way to calculate degree of fluctuations of side-chains ( degree of torsion's dynamics) from the different trajectories and to compare it ? E.g I have one protein in two different (apo and holo) forms. I've calculated two trajectories for both structures and observed different degree of dynamics in case of each structure ( e.g fluctuations in case of apo form were more frequent than in case of liganded form). IS there any way to direct mesure and comprison of such side-chain's dynamics for two trajectories? James 2012/8/15, James Starlight jmsstarli...@gmail.com: Dear Gromacs users! I want to obtain Cross-correlation maps ( for indication of the cross-correlated fluctuations of the residues). The example of such maps can be found here http://pubs.acs.org/doi/abs/10.1021/ja076046a I found that modificied version of the G_covar from users contributions can do such things. But because of the older version of that program (3.3.3) I've obtained the below error using it with 4.5.5 gromacs Reading file md_GO.tpr, VERSION 4.5.5 (single precision) --- Program g_covar_mod, VERSION 3.3.3 Source code file: tpxio.c, line: 1192 Fatal error: reading tpx file (md_GO.tpr) version 73 with version 40 program --- Is there newest versions of the G_covar for such things or alternativelly any others ways to calculate such correlations maps ? Thanks for help James -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Only plain text messages are allowed! * 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 gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
[gmx-users] Gromacs-4.5.5 on GPU : Floating point exception
Gromacs-4.5.5 is built on GPU with OpenMM-4.1.1 binary package. But the standard gromacs gpu benchmark - dhfr failed with Floating point exception. The benchmark name is : dhfr/GPU/dhfr-impl-1nm.bench Error Output: Reading file topol.tpr, VERSION 4.5.1-dev-20100917-b1d66 (single precision) WARNING: OpenMM does not support leap-frog, will use velocity-verlet integrator. WARNING: OpenMM supports only Andersen thermostat with the md/md-vv/md-vv-avek integrators. WARNING: OpenMM provides contraints as a combination of SHAKE, SETTLE and CCMA. Accuracy is based on the SHAKE tolerance set by the shake_tol option. WARNING: The OBC scale factors alpha, beta and gamma are hardcoded in OpenMM with the default Gromacs values. Floating point exception What's wrong here? The standard installation steps given at http://www.gromacs.org/Documentation/Installation_Instructions/GPUs are followed. -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Only plain text messages are allowed! * 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 gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
[gmx-users] what's the difference between gen_seed and ld_seed?
Am 23.08.2012 23:11, schrieb gmx-users-requ...@gromacs.org: hello : I am a little bit confused about the difference between gen_seed and ld_seed. I checked the manual, it is said: gen_seed used to initialize random generator for random velocities, when gen_seed is set to -1, the seed is calculated from the process ID number. This is often used coupled with gen_vel which is Generate velocities in grompp according to a Maxwell distribution at temperaturegen_temp [K], with random seed gen_seed. This is only meaningful with integratormd. As indicated in Jonhn E.Kerrigan's tutorial, gen_seed=-1 is always turned on in NVT, NPT and MD production step. However, in Justin's Lysozyme in Water tutorial, this option is only turned on in NVT and the following NPT and MD production were turned off. Yes because you undo your equilibration when you reassign initial velocities using gen_seed in NPT and MD production (NPT) (md/md-vv integrator). I am just wondering which option would be better for our simulations? How about ld_seed? here is the statement from manual: ld_seed: (1993) [integer] used to initialize random generator for thermal noise for stochastic and Brownian dynamics. When ld_seed is set to -1, the seed is calculated from the process ID. When running BD or SD on multiple processors, each processor uses a seed equal to ld_seed plus the processor number. when should we turn this option on? As stated, you use this option for use with the bd or sd integrators. Also the 'v-rescale' thermostat uses the ld_seed. -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Only plain text messages are allowed! * 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 gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
[gmx-users] Re: equilibrium stage runs for too long
Dear Vitaly, Thank you so much for your suggestion!! I've run the simultion using one core as you suggested, and I can see two outputs on step 0 and step 100 now ( shown on the bottom), but it still took a while to only generate these two outputs. Maybe it is indeed the problem of scaling as mentioned by Justin? There are: 34001 Atoms Charge group distribution at step 0: 1953 1973 2103 2032 2044 2072 2063 1967 Grid: 3 x 7 x 18 cells Constraining the starting coordinates (step 0) Constraining the coordinates at t0-dt (step 0) RMS relative constraint deviation after constraining: 1.99e-05 Initial temperature: 300.057 K Started mdrun on node 0 Fri Aug 24 10:37:20 2012 Step Time Lambda 00.00.0 Energies (kJ/mol) AngleProper Dih. Ryckaert-Bell. Improper Dih. LJ-14 6.66008e+033.06081e+037.94034e+022.55031e+028.61710e+02 Coulomb-14LJ (SR) Disper. corr. Coulomb (SR) Coul. recip. 2.40538e+035.11091e+04 -3.00143e+03 -4.49637e+05 -8.11103e+04 Position Rest. PotentialKinetic En. Total Energy Conserved En. 5.53187e+01 -4.68547e+058.57672e+04 -3.82780e+05 -3.82780e+05 Temperature Pres. DC (bar) Pressure (bar) Constr. rmsd 3.00450e+02 -9.73436e+01 -2.97627e+033.74608e-05 DD step 4 load imb.: force 6.0% DD step 99 load imb.: force 6.1% Step Time Lambda 1000.20.0 Energies (kJ/mol) AngleProper Dih. Ryckaert-Bell. Improper Dih. LJ-14 1.10274e+042.82478e+031.23525e+032.24466e+021.54991e+03 Coulomb-14LJ (SR) Disper. corr. Coulomb (SR) Coul. recip. 4.67051e+036.34091e+04 -3.00143e+03 -4.31792e+05 -8.19081e+04 Position Rest. PotentialKinetic En. Total Energy Conserved En. 1.44619e+04 -4.17298e+057.13854e+04 -3.45913e+05 -3.82303e+05 Temperature Pres. DC (bar) Pressure (bar) Constr. rmsd 2.50069e+02 -9.73436e+01 -3.79866e+023.45085e-06 -- View this message in context: http://gromacs.5086.n6.nabble.com/equilibrium-stage-runs-for-too-long-tp5000461p5000483.html Sent from the GROMACS Users Forum mailing list archive at Nabble.com. -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Only plain text messages are allowed! * 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 gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
Re: [gmx-users] Re: equilibrium stage runs for too long
On 8/24/12 6:01 AM, Clare wrote: Dear Vitaly, Thank you so much for your suggestion!! I've run the simultion using one core as you suggested, and I can see two outputs on step 0 and step 100 now ( shown on the bottom), but it still took a while to only generate these two outputs. Maybe it is indeed the problem of scaling as mentioned by Justin? There are: 34001 Atoms The reason that it is slow now is that you have 34001 atoms being calculated on one core. What that demonstrates is that your simulation is stable and likely did not crash when run in parallel. My previous comments about scaling issues were simply a guess based on incomplete information. Your log file from the MPI run indicated the job never even started. That's not an issue of scaling; it's an issue of functionality. It looks to me like either mdrun_mpi or your MPI implementation is broken or not being invoked properly. Try the test case I suggested before, something that doesn't involve Gromacs. -Justin Charge group distribution at step 0: 1953 1973 2103 2032 2044 2072 2063 1967 Grid: 3 x 7 x 18 cells Constraining the starting coordinates (step 0) Constraining the coordinates at t0-dt (step 0) RMS relative constraint deviation after constraining: 1.99e-05 Initial temperature: 300.057 K Started mdrun on node 0 Fri Aug 24 10:37:20 2012 Step Time Lambda 00.00.0 Energies (kJ/mol) AngleProper Dih. Ryckaert-Bell. Improper Dih. LJ-14 6.66008e+033.06081e+037.94034e+022.55031e+028.61710e+02 Coulomb-14LJ (SR) Disper. corr. Coulomb (SR) Coul. recip. 2.40538e+035.11091e+04 -3.00143e+03 -4.49637e+05 -8.11103e+04 Position Rest. PotentialKinetic En. Total Energy Conserved En. 5.53187e+01 -4.68547e+058.57672e+04 -3.82780e+05 -3.82780e+05 Temperature Pres. DC (bar) Pressure (bar) Constr. rmsd 3.00450e+02 -9.73436e+01 -2.97627e+033.74608e-05 DD step 4 load imb.: force 6.0% DD step 99 load imb.: force 6.1% Step Time Lambda 1000.20.0 Energies (kJ/mol) AngleProper Dih. Ryckaert-Bell. Improper Dih. LJ-14 1.10274e+042.82478e+031.23525e+032.24466e+021.54991e+03 Coulomb-14LJ (SR) Disper. corr. Coulomb (SR) Coul. recip. 4.67051e+036.34091e+04 -3.00143e+03 -4.31792e+05 -8.19081e+04 Position Rest. PotentialKinetic En. Total Energy Conserved En. 1.44619e+04 -4.17298e+057.13854e+04 -3.45913e+05 -3.82303e+05 Temperature Pres. DC (bar) Pressure (bar) Constr. rmsd 2.50069e+02 -9.73436e+01 -3.79866e+023.45085e-06 -- View this message in context: http://gromacs.5086.n6.nabble.com/equilibrium-stage-runs-for-too-long-tp5000461p5000483.html Sent from the GROMACS Users Forum mailing list archive at Nabble.com. -- Justin A. Lemkul, Ph.D. Research Scientist Department of Biochemistry Virginia Tech Blacksburg, VA jalemkul[at]vt.edu | (540) 231-9080 http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Only plain text messages are allowed! * 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 gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
[gmx-users] Re: equilibrium stage runs for too long
Dear Justin, I just tried running it using one node with 12 cores, and it's been completed successfully within a short time! Yes you are right, the system is under reparation and this probably causes those problems. I'm really grateful for your and Vitaly's help. Thank you very much!!! Clare -- View this message in context: http://gromacs.5086.n6.nabble.com/equilibrium-stage-runs-for-too-long-tp5000461p5000488.html Sent from the GROMACS Users Forum mailing list archive at Nabble.com. -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Only plain text messages are allowed! * 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 gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
Re: [gmx-users] deltaG from PMF
I did g_wham calculation as you suggested for: 10-30 ns 30-50 ns 10-50 ns In some cases of PMF curves the drift is 0.3 kcal/mol but in some cases it is over 1.5 kcal/mol so huge drifts. In terms of huge drifts what would you suggest ? Steven On Tue, Aug 21, 2012 at 5:48 PM, Justin Lemkul jalem...@vt.edu wrote: On 8/21/12 12:46 PM, Steven Neumann wrote: Thanks Thomas. Justin, could you please comment on this? I agree with everything Thomas has said. There's not really anything to say. -Justin Steven On Tue, Aug 21, 2012 at 5:37 PM, Thomas Schlesier schl...@uni-mainz.de wrote: Am 21.08.2012 18:22, schrieb gmx-users-requ...@gromacs.org: On Tue, Aug 21, 2012 at 4:49 PM, Thomas Schlesierschl...@uni-mainz.de wrote: Since your simulations of the individual windows are about 50 ns, you could first dismiss the first 10 ns for equilibration, and then perform the WHAM analysis for 10-30 ns and 30-50 ns. If everything is fine, you should see no drift. If you want to have more data for the analysis you could also use 5ns ; 5-27.5ns and 27.5-50ns. From the PMF it seems that the equilibrium state should be around 0.6 nm. To be sure, you can perform a normal simulation (without any restraints) from you initial starting window (~0.4nm) and a window near the minima (~0.6nm). Then after the equilibration phase, look at the distribution of the distance along the reaction coordinate. If in both cases the maximum is at ~0.6nm, this should be the 'true' equilibrium state of the system (instead of the first window of the PMF calculation) and i would measure \Delta G from this point. Greetings Thomas Thanks Thomas for this but finally I realised that my first configuration corresponds to 0.6 nm which is the minima so I take the free energy difference based on this value and plateau. I want also to calculate error bars. Would you do this: Final PMF curve for 10-50 ns Error bars from: g_wham -b 3 -e 4 g_wham -b 5 -e 6 Think this approach would be good to see if you have any drifts. But for error bars there is something implemented in 'g_wham'. But i never used it, since for my system umbrella sampling is not really applicable, only TI. So i can't comment on it, if there is anything one should be aware of, or similar. But 'g_wham -h' prints some info about how to use the error analysis Greetings Thomas Steven Am 21.08.2012 17:25, schriebgmx-users-requ...@gromacs.org: On Tue, Aug 21, 2012 at 4:21 PM, Justin Lemkuljalem...@vt.edu wrote: On 8/21/12 11:18 AM, Steven Neumann wrote: On Tue, Aug 21, 2012 at 4:13 PM, Justin Lemkuljalem...@vt.edu wrote: On 8/21/12 11:09 AM, Steven Neumann wrote: On Tue, Aug 21, 2012 at 3:48 PM, Justin Lemkuljalem...@vt.edu wrote: On 8/21/12 10:42 AM, Steven Neumann wrote: Please see the example of the plot from US simulations and WHAM: http://speedy.sh/Ecr3A/PMF.JPG First grompp of frame 0 corresponds to 0.8 nm - this is what was shown by grompp at the end. The mdp file: ; Run parameters define = -DPOSRES_T integrator = md; leap-frog integrator nsteps = 2500 ; 100ns dt = 0.002 ; 2 fs tinit = 0 nstcomm = 10 ; Output control nstxout = 0 ; save coordinates every 100 ps nstvout = 0 ; save velocities every nstxtcout = 5; every 10 ps nstenergy = 1000 ; Bond parameters continuation= no ; first dynamics run constraint_algorithm = lincs; holonomic constraints constraints = all-bonds ; all bonds (even heavy atom-H bonds) constrained ; Neighborsearching ns_type = grid ; search neighboring grid cells nstlist = 5 ; 10 fs vdwtype = Switch rvdw-switch = 1.0 rlist = 1.4 ; short-range neighborlist cutoff (in nm) rcoulomb= 1.4 ; short-range electrostatic cutoff (in nm) rvdw= 1.2 ; short-range van der Waals cutoff (in nm) ewald_rtol = 1e-5 ; relative strenght of the Ewald-shifted potential rcoulomb ; Electrostatics coulombtype = PME ; Particle Mesh Ewald for long-range electrostatics pme_order = 4 ; cubic interpolation fourierspacing = 0.12 ; grid spacing for FFT fourier_nx = 0 fourier_ny = 0 fourier_nz = 0 optimize_fft= yes ; Temperature coupling is on tcoupl = V-rescale ; modified Berendsen thermostat tc_grps = Protein LIG_Water_and_ions ; two coupling groups - more accurate tau_t = 0.1 0.1 ; time constant, in ps ref_t = 318 318 ; reference temperature, one for each group, in K ; Pressure coupling is on 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,
[gmx-users] BAR / g_bar problems
Hi, we have terrible problems to get reasonable results from BAR free energy calculations. We basically follow Justin's tutorial for solvation of methanol, ethanol, diethyl-ether and neopentane in water using OPLS but get very strange results. The free energy curves are here: http://folding.bmc.uu.se/images/koko.jpg Molecule Energy Exper Methanol -8-21 Ethanol -9-21 Diethylether -11-7 Neopentane -10 +11 any clue? -- David van der Spoel, Ph.D., Professor of Biology Dept. of Cell Molec. Biol., Uppsala University. Box 596, 75124 Uppsala, Sweden. Phone: +46184714205. sp...@xray.bmc.uu.sehttp://folding.bmc.uu.se -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Only plain text messages are allowed! * 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 gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
[gmx-users] CHARMM and FEP calculations
Hi I asked a similar question a year ago about implementing CMAP in FEP calculations for gromacs 4.5.4. At the time I could not get the simulation to run with CMAP as it came up with: Fatal error: Function type CMAP Dih. not implemented in ip_pert Helpfully, Mark suggested a modification to the code to get round this issue, but I am unclear as to whether in the newer version of gromacs (4.5.5) CMAP terms can actually be scaled from state A to state B or whether for all values of lambda this part of the potential remains in state A even if I were to manually write state B CMAP terms in the topology file and the code no longer gives a warning. In addition, for CHARMM, and other FFs, there are some dihedrals which have multiple potentials, using the dihedral function type 9. I have manually specified the state B parameters for these (and all other dihedrals) in my topology file but get the error message: Cannot automatically perturb a torsion with multiple terms to different form. Please specify perturbed parameters manually for this torsion in your topology! Is there a way to implement scaling of these multiply defined dihedrals in FEP calculations? Many thanks for your help Louise -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Only plain text messages are allowed! * 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 gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
[gmx-users] g_rms and g_rmsdist on initial structure
Hi, For example, I have a A.pdb as a initial structure file. And I just used pdb2gmx on it to generate another B.pdb file with GROMOS96 43a1 as its force filed. Then I select C-alpha atoms to calculate RMSD. echo 3 | g_rms -f B.pdb -s A.pdb I suppose the RMSD value should be 0, but the value is high to about 0.5nm. Can someone explain for me? Sincerely yours, Hsin-Lin -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Only plain text messages are allowed! * 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 gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
Re: [gmx-users] g_rms and g_rmsdist on initial structure
Hi, Since gromos-forcefields are not strictly all-atom forcefields, there might be a mismatch between atoms in the two structures. Best, Erik 24 aug 2012 kl. 15.35 skrev Hsin-Lin Chiang: Hi, For example, I have a A.pdb as a initial structure file. And I just used pdb2gmx on it to generate another B.pdb file with GROMOS96 43a1 as its force filed. Then I select C-alpha atoms to calculate RMSD. echo 3 | g_rms -f B.pdb -s A.pdb I suppose the RMSD value should be 0, but the value is high to about 0.5nm. Can someone explain for me? Sincerely yours, Hsin-Lin -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Only plain text messages are allowed! * 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 gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists --- Erik Marklund, PhD Dept. of Cell and Molecular Biology, Uppsala University. Husargatan 3, Box 596,75124 Uppsala, Sweden phone:+46 18 471 6688fax: +46 18 511 755 er...@xray.bmc.uu.se http://www2.icm.uu.se/molbio/elflab/index.html -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Only plain text messages are allowed! * 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 gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
[gmx-users] deltaG from PMF
As a first step, i would shift all curves so, that the energy of the minium is for all plots at the same (aribarity) value. The minimum should be the point which has sampled the best. If you shift then all values, it should be easier to spot differences between the plots. And probably make the analysis for every 10ns slices 10-20, 20-30, ... then it's easier to see if you have a drift or fluctuations. If you identify problematic regions you could do there additional umbrella simulations, probably with higher force constants, in order to sample that region better. In theory the PMF should converge if you wait long enough, till the system equilibrates under the external restraints (how long this takes? nobody knows). On could probably wait a long time, big question here is what is the biggest error you would tolerate. On the other hand one can't simulate forever... Look for what other people used for simulation-times for umbrella sampling (i have the impression that 50ns is rather long, but i could fool me here) and what they did to estimate if the calculations are converged. Then decide, what to do. Think nobody here would / will tell you this or this error is ok ... but if you do something reasonable, which is in accord what others do / did it should be ok. Am 24.08.2012 14:59, schrieb gmx-users-requ...@gromacs.org: I did g_wham calculation as you suggested for: 10-30 ns 30-50 ns 10-50 ns In some cases of PMF curves the drift is 0.3 kcal/mol but in some cases it is over 1.5 kcal/mol so huge drifts. In terms of huge drifts what would you suggest ? Steven On Tue, Aug 21, 2012 at 5:48 PM, Justin Lemkul jalem...@vt.edu wrote: On 8/21/12 12:46 PM, Steven Neumann wrote: Thanks Thomas. Justin, could you please comment on this? I agree with everything Thomas has said. There's not really anything to say. -Justin Steven On Tue, Aug 21, 2012 at 5:37 PM, Thomas Schlesierschl...@uni-mainz.de wrote: Am 21.08.2012 18:22, schriebgmx-users-requ...@gromacs.org: On Tue, Aug 21, 2012 at 4:49 PM, Thomas Schlesierschl...@uni-mainz.de wrote: Since your simulations of the individual windows are about 50 ns, you could first dismiss the first 10 ns for equilibration, and then perform the WHAM analysis for 10-30 ns and 30-50 ns. If everything is fine, you should see no drift. If you want to have more data for the analysis you could also use 5ns ; 5-27.5ns and 27.5-50ns. From the PMF it seems that the equilibrium state should be around 0.6 nm. To be sure, you can perform a normal simulation (without any restraints) from you initial starting window (~0.4nm) and a window near the minima (~0.6nm). Then after the equilibration phase, look at the distribution of the distance along the reaction coordinate. If in both cases the maximum is at ~0.6nm, this should be the 'true' equilibrium state of the system (instead of the first window of the PMF calculation) and i would measure \Delta G from this point. Greetings Thomas Thanks Thomas for this but finally I realised that my first configuration corresponds to 0.6 nm which is the minima so I take the free energy difference based on this value and plateau. I want also to calculate error bars. Would you do this: Final PMF curve for 10-50 ns Error bars from: g_wham -b 3 -e 4 g_wham -b 5 -e 6 Think this approach would be good to see if you have any drifts. But for error bars there is something implemented in 'g_wham'. But i never used it, since for my system umbrella sampling is not really applicable, only TI. So i can't comment on it, if there is anything one should be aware of, or similar. But 'g_wham -h' prints some info about how to use the error analysis Greetings Thomas Steven Am 21.08.2012 17:25,schriebgmx-users-requ...@gromacs.org: On Tue, Aug 21, 2012 at 4:21 PM, Justin Lemkuljalem...@vt.edu wrote: On 8/21/12 11:18 AM, Steven Neumann wrote: On Tue, Aug 21, 2012 at 4:13 PM, Justin Lemkuljalem...@vt.edu wrote: On 8/21/12 11:09 AM, Steven Neumann wrote: On Tue, Aug 21, 2012 at 3:48 PM, Justin Lemkuljalem...@vt.edu wrote: On 8/21/12 10:42 AM, Steven Neumann wrote: Please see the example of the plot from US simulations and WHAM: http://speedy.sh/Ecr3A/PMF.JPG First grompp of frame 0 corresponds to 0.8 nm - this is what was shown by grompp at the end. The mdp file: ; Run parameters define = -DPOSRES_T integrator = md; leap-frog integrator nsteps = 2500 ; 100ns dt = 0.002 ; 2 fs tinit = 0 nstcomm = 10 ; Output control nstxout = 0 ; save coordinates every 100 ps nstvout = 0 ; save velocities every nstxtcout = 5; every 10 ps nstenergy = 1000 ; Bond parameters continuation= no ; first dynamics run constraint_algorithm = lincs; holonomic constraints constraints = all-bonds ; all bonds (even heavy atom-H bonds) constrained ;
[gmx-users] Re: g_rms and g_rmsdist on initial structure
-- View this message in context: http://gromacs.5086.n6.nabble.com/g-rms-and-g-rmsdist-on-initial-structure-tp5000492p5000495.html Sent from the GROMACS Users Forum mailing list archive at Nabble.com. -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Only plain text messages are allowed! * 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 gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
[gmx-users] Re: g_rms and g_rmsdist on initial structure
Thank you for your reply. I checked the coordinates of all C-alpha atoms between two .pdb files. They are exactly the same. So I don't understand why g_rms(g_rmsdist) selecting group C-alpha give me high values. Sincerely yours, Hsin-Lin Erik Marklund wrote Hi, Since gromos-forcefields are not strictly all-atom forcefields, there might be a mismatch between atoms in the two structures. Best, Erik 24 aug 2012 kl. 15.35 skrev Hsin-Lin Chiang: Hi, For example, I have a A.pdb as a initial structure file. And I just used pdb2gmx on it to generate another B.pdb file with GROMOS96 43a1 as its force filed. Then I select C-alpha atoms to calculate RMSD. echo 3 | g_rms -f B.pdb -s A.pdb I suppose the RMSD value should be 0, but the value is high to about 0.5nm. Can someone explain for me? Sincerely yours, Hsin-Lin -- gmx-users mailing listgmx-users@ http://lists.gromacs.org/mailman/listinfo/gmx-users * Only plain text messages are allowed! * 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 gmx-users-request@. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists --- Erik Marklund, PhD Dept. of Cell and Molecular Biology, Uppsala University. Husargatan 3, Box 596,75124 Uppsala, Sweden phone:+46 18 471 6688fax: +46 18 511 755 erikm@.uu http://www2.icm.uu.se/molbio/elflab/index.html -- gmx-users mailing listgmx-users@ http://lists.gromacs.org/mailman/listinfo/gmx-users * Only plain text messages are allowed! * 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 gmx-users-request@. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists -- View this message in context: http://gromacs.5086.n6.nabble.com/g-rms-and-g-rmsdist-on-initial-structure-tp5000492p5000496.html Sent from the GROMACS Users Forum mailing list archive at Nabble.com. -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Only plain text messages are allowed! * 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 gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
[gmx-users] KALP15-DPPC Tutorial
Hi, I am doing the simulation of kalp15 in DPPC following the Justin's tutorial. In the 2. Pack the lipids around the protein step, to get the correct area per lipid, I enter these commands: #perl inflategro.pl system.gro 4 DPPC 14 system_inflated.gro 5 area.dat Updating the .top file and adding DPPC 126 to it. First minimization: Energy minimization with restrained protein #grompp -f minim.mdp -c system_inflated.gro -p topol.top -o em.tpr #mdrun -deffnm em scale down the lipids by a factor of 0.95 #perl inflategro.pl em.gro 0.95 DPPC 0 system_shrink1.gro 5 area_shrink1.dat Second minimization: #grompp -f minim.mdp -c system_shrink1.gro -p topol.top -o em_shrink1.tpr Then, I get this error: Fatal error: number of coordinates in coordinate file (system_shrink1.gro, 6400) does not match topology (topol.top, 6438) Is this not the correct command to start the second minimization? Please help me. Your suggestions are welcomed. Sincerely, Shima -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Only plain text messages are allowed! * 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 gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
Re: [gmx-users] BAR / g_bar problems
Hi, David- Perhaps we can work with you to compare the internal g_bar implementation with our external BAR and MBAR implementation run from the .dhdl.xvg data output in 4.6. It would probably be easier to run these calculations after the optimization updates are added in the next few days(?), but we can plan now. Note that withe errors this big, 100-200 ps should be enough to see what's going on, so it can be done rapidly. Drop me a line off the list to figure out details? Best, Michael On Fri, Aug 24, 2012 at 9:22 AM, David van der Spoel sp...@xray.bmc.uu.se wrote: Hi, we have terrible problems to get reasonable results from BAR free energy calculations. We basically follow Justin's tutorial for solvation of methanol, ethanol, diethyl-ether and neopentane in water using OPLS but get very strange results. The free energy curves are here: http://folding.bmc.uu.se/images/koko.jpg Molecule Energy Exper Methanol -8-21 Ethanol -9-21 Diethylether -11-7 Neopentane -10 +11 any clue? -- David van der Spoel, Ph.D., Professor of Biology Dept. of Cell Molec. Biol., Uppsala University. Box 596, 75124 Uppsala, Sweden. Phone: +46184714205. sp...@xray.bmc.uu.sehttp://folding.bmc.uu.se -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Only plain text messages are allowed! * 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 gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Only plain text messages are allowed! * 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 gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
Re: [gmx-users] Re: g_rms and g_rmsdist on initial structure
Because it uses the atom index, which is most likely NOT the same as your processed file will not contain certain hydrogens, which are present in A.pdb. 24 aug 2012 kl. 16.33 skrev Hsinlin: Thank you for your reply. I checked the coordinates of all C-alpha atoms between two .pdb files. They are exactly the same. So I don't understand why g_rms(g_rmsdist) selecting group C-alpha give me high values. Sincerely yours, Hsin-Lin Erik Marklund wrote Hi, Since gromos-forcefields are not strictly all-atom forcefields, there might be a mismatch between atoms in the two structures. Best, Erik 24 aug 2012 kl. 15.35 skrev Hsin-Lin Chiang: Hi, For example, I have a A.pdb as a initial structure file. And I just used pdb2gmx on it to generate another B.pdb file with GROMOS96 43a1 as its force filed. Then I select C-alpha atoms to calculate RMSD. echo 3 | g_rms -f B.pdb -s A.pdb I suppose the RMSD value should be 0, but the value is high to about 0.5nm. Can someone explain for me? Sincerely yours, Hsin-Lin -- gmx-users mailing listgmx-users@ http://lists.gromacs.org/mailman/listinfo/gmx-users * Only plain text messages are allowed! * 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 gmx-users-request@. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists --- Erik Marklund, PhD Dept. of Cell and Molecular Biology, Uppsala University. Husargatan 3, Box 596,75124 Uppsala, Sweden phone:+46 18 471 6688fax: +46 18 511 755 erikm@.uu http://www2.icm.uu.se/molbio/elflab/index.html -- gmx-users mailing listgmx-users@ http://lists.gromacs.org/mailman/listinfo/gmx-users * Only plain text messages are allowed! * 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 gmx-users-request@. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists -- View this message in context: http://gromacs.5086.n6.nabble.com/g-rms-and-g-rmsdist-on-initial-structure-tp5000492p5000496.html Sent from the GROMACS Users Forum mailing list archive at Nabble.com. -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Only plain text messages are allowed! * 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 gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists --- Erik Marklund, PhD Dept. of Cell and Molecular Biology, Uppsala University. Husargatan 3, Box 596,75124 Uppsala, Sweden phone:+46 18 471 6688fax: +46 18 511 755 er...@xray.bmc.uu.se http://www2.icm.uu.se/molbio/elflab/index.html -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Only plain text messages are allowed! * 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 gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
Re: [gmx-users] deltaG from PMF
Thanks for this! Steven On Fri, Aug 24, 2012 at 3:26 PM, Thomas Schlesier schl...@uni-mainz.de wrote: As a first step, i would shift all curves so, that the energy of the minium is for all plots at the same (aribarity) value. The minimum should be the point which has sampled the best. If you shift then all values, it should be easier to spot differences between the plots. And probably make the analysis for every 10ns slices 10-20, 20-30, ... then it's easier to see if you have a drift or fluctuations. If you identify problematic regions you could do there additional umbrella simulations, probably with higher force constants, in order to sample that region better. In theory the PMF should converge if you wait long enough, till the system equilibrates under the external restraints (how long this takes? nobody knows). On could probably wait a long time, big question here is what is the biggest error you would tolerate. On the other hand one can't simulate forever... Look for what other people used for simulation-times for umbrella sampling (i have the impression that 50ns is rather long, but i could fool me here) and what they did to estimate if the calculations are converged. Then decide, what to do. Think nobody here would / will tell you this or this error is ok ... but if you do something reasonable, which is in accord what others do / did it should be ok. Am 24.08.2012 14:59, schrieb gmx-users-requ...@gromacs.org: I did g_wham calculation as you suggested for: 10-30 ns 30-50 ns 10-50 ns In some cases of PMF curves the drift is 0.3 kcal/mol but in some cases it is over 1.5 kcal/mol so huge drifts. In terms of huge drifts what would you suggest ? Steven On Tue, Aug 21, 2012 at 5:48 PM, Justin Lemkul jalem...@vt.edu wrote: On 8/21/12 12:46 PM, Steven Neumann wrote: Thanks Thomas. Justin, could you please comment on this? I agree with everything Thomas has said. There's not really anything to say. -Justin Steven On Tue, Aug 21, 2012 at 5:37 PM, Thomas Schlesierschl...@uni-mainz.de wrote: Am 21.08.2012 18:22, schriebgmx-users-requ...@gromacs.org: On Tue, Aug 21, 2012 at 4:49 PM, Thomas Schlesierschl...@uni-mainz.de wrote: Since your simulations of the individual windows are about 50 ns, you could first dismiss the first 10 ns for equilibration, and then perform the WHAM analysis for 10-30 ns and 30-50 ns. If everything is fine, you should see no drift. If you want to have more data for the analysis you could also use 5ns ; 5-27.5ns and 27.5-50ns. From the PMF it seems that the equilibrium state should be around 0.6 nm. To be sure, you can perform a normal simulation (without any restraints) from you initial starting window (~0.4nm) and a window near the minima (~0.6nm). Then after the equilibration phase, look at the distribution of the distance along the reaction coordinate. If in both cases the maximum is at ~0.6nm, this should be the 'true' equilibrium state of the system (instead of the first window of the PMF calculation) and i would measure \Delta G from this point. Greetings Thomas Thanks Thomas for this but finally I realised that my first configuration corresponds to 0.6 nm which is the minima so I take the free energy difference based on this value and plateau. I want also to calculate error bars. Would you do this: Final PMF curve for 10-50 ns Error bars from: g_wham -b 3 -e 4 g_wham -b 5 -e 6 Think this approach would be good to see if you have any drifts. But for error bars there is something implemented in 'g_wham'. But i never used it, since for my system umbrella sampling is not really applicable, only TI. So i can't comment on it, if there is anything one should be aware of, or similar. But 'g_wham -h' prints some info about how to use the error analysis Greetings Thomas Steven Am 21.08.2012 17:25,schriebgmx-users-requ...@gromacs.org: On Tue, Aug 21, 2012 at 4:21 PM, Justin Lemkuljalem...@vt.edu wrote: On 8/21/12 11:18 AM, Steven Neumann wrote: On Tue, Aug 21, 2012 at 4:13 PM, Justin Lemkuljalem...@vt.edu wrote: On 8/21/12 11:09 AM, Steven Neumann wrote: On Tue, Aug 21, 2012 at 3:48 PM, Justin Lemkuljalem...@vt.edu wrote: On 8/21/12 10:42 AM, Steven Neumann wrote: Please see the example of the plot from US simulations and WHAM: http://speedy.sh/Ecr3A/PMF.JPG First grompp of frame 0 corresponds to 0.8 nm - this is what was shown by grompp at the end. The mdp file: ; Run parameters define = -DPOSRES_T integrator = md ; leap-frog integrator nsteps = 2500 ; 100ns dt = 0.002 ; 2 fs tinit = 0 nstcomm = 10 ; Output control nstxout = 0
[gmx-users] Electrostatics and van Der Waals in Gromacs
Dear Gromacs users, I am trying to get a full understanding of how the energy of electrostatics and VDW interactions are calculated in Gromacs, but I am not sure I have the right picture in mind. Could someone let me know if the procedure below is correct? So, let's say I have a PDB file, and I want to calculate the electrostatic interaction between a pair of atoms. Can I do this: 1) Take the charges q_1 and q_2 of the two atoms from those in the aminoacids.rtp file of the force field. Example: From aminoacids.rtp: [ ALA ] [ atoms ] Nopls_238 -0.500 1 Hopls_2410.300 1 CAopls_224B 0.140 1 HAopls_1400.060 1 CBopls_135 -0.180 2 If atom 1 is CA from alanine, I would choose q_1 = 0.140? 2) Choose a value for the dielectric constant = I read that epsilon_r = 4 is a reasonable value for the inside of a protein, is that correct? 3) Calculate the distance r_12 between the two atoms 4) Simply evaluate: E_electrostatics = 138.935 * q1 * q2 / (epsilon_r * r_12); Is that the way Gromacs is calculating the actual energy or am I missing something? Is this at least a good approximation, or not at all? Is the aminoacids.rtp the right place to look for the two charges? Is epsilon_r = 4 reasonable? Similarly, if I want to calculate the van der Waals energy between my two atoms, can I do this: 1) Take the values of epsilon_1, epsilon_2, sigma_1 and sigma_2 from the last two columns of the ffbonded.itp for the line corresponding to my two atoms/ Example: From ffnonbonded.itp: ; name bond_typemasscharge ptype sigma epsilon opls_001 C 6 12.01100 0.500 A 3.75000e-01 4.39320e-01 ; SIG opls_002 O 8 15.99940-0.500 A 2.96000e-01 8.78640e-01 ; SIG (...) if atom 1 is opls_001 and atom 2 opls_002, epsilon_1 = 4.39320e-01, sigma_1 = 3.75000e-01, etc... 2) calculate epsilon = sqrt(epsilon_1*epsilon_2) and sigma = sqrt(sigma_1*sigma_2) 3) Calculate the distance r_12 between the two atoms 4) Simply evaluate E_vdw = 4*epsilon * [ (sigma/r_12)^12 - (sigma/r_12)^6 ] (I am aware that there are other ways to do step 2) and other potentials than L-J, but let's assume they are those my forcefield is using). Same question: Is that how gromacs is calculating E_vdw or am I missing something important here? One last detail: is there a specific term for hydrogen bonds or their energy is simply included in the sum of electrostatics and VDW? My apologies for this long and probably very basic question and many thanks in advance... Antoine -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Only plain text messages are allowed! * 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 gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
[gmx-users] g_tune_pme restart
Dear: I use g_tune_pme for MD production but it clashed before the job finished since it is over cluster walltime limitation. I get the following information from perf.out: mpirun -np 144 mdrun -npme -1 -s tuned.tpr -v -o md.trr -c md.gro -e md.edr -g md.log I am just wondering, is it correct using the following command to append the jobs? mpirun -np 144 mdrun -npme -1 -s tuned.tpr -v -f md.trr -e md.edr -g md.log -o md.trr -cpi -append thank you very much best Albert -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Only plain text messages are allowed! * 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 gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
Re: [gmx-users] KALP15-DPPC Tutorial
On 8/24/12 10:48 AM, Shima Arasteh wrote: Hi, I am doing the simulation of kalp15 in DPPC following the Justin's tutorial. In the 2. Pack the lipids around the protein step, to get the correct area per lipid, I enter these commands: #perl inflategro.pl system.gro 4 DPPC 14 system_inflated.gro 5 area.dat Updating the .top file and adding DPPC 126 to it. First minimization: Energy minimization with restrained protein #grompp -f minim.mdp -c system_inflated.gro -p topol.top -o em.tpr #mdrun -deffnm em scale down the lipids by a factor of 0.95 #perl inflategro.pl em.gro 0.95 DPPC 0 system_shrink1.gro 5 area_shrink1.dat Second minimization: #grompp -f minim.mdp -c system_shrink1.gro -p topol.top -o em_shrink1.tpr Then, I get this error: Fatal error: number of coordinates in coordinate file (system_shrink1.gro, 6400) does not match topology (topol.top, 6438) Is this not the correct command to start the second minimization? After the first inflation step and lipids are removed (and the corresponding change reflected in the topology), the topology is invariant. I can offer no explanation as to why the first invocation of grompp worked and the subsequent one did not, unless you have made other changes or are using the wrong file. -Justin -- Justin A. Lemkul, Ph.D. Research Scientist Department of Biochemistry Virginia Tech Blacksburg, VA jalemkul[at]vt.edu | (540) 231-9080 http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Only plain text messages are allowed! * 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 gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
Re: [gmx-users] g_tune_pme restart
On Fri, Aug 24, 2012 at 12:23 PM, Albert mailmd2...@gmail.com wrote: Dear: I use g_tune_pme for MD production but it clashed before the job finished since it is over cluster walltime limitation. I get the following information from perf.out: mpirun -np 144 mdrun -npme -1 -s tuned.tpr -v -o md.trr -c md.gro -e md.edr -g md.log I am just wondering, is it correct using the following command to append the jobs? mpirun -np 144 mdrun -npme -1 -s tuned.tpr -v -f md.trr -e md.edr -g md.log -o md.trr -cpi -append mdrun doesn't have -f as option. It has to be -o. Otherwise it seems OK. Roland thank you very much best Albert -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Only plain text messages are allowed! * 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 gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists -- ORNL/UT Center for Molecular Biophysics cmb.ornl.gov 865-241-1537, ORNL PO BOX 2008 MS6309 -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Only plain text messages are allowed! * 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 gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists