[gmx-users] Ligand charge issues
Hi Gromacs users, I am studying the protein-ligand interaction using amber99sb-ILDN force field in gromacs 4.6.2. To create the ligand topology (lipid A), I have used online version of ACPYPE/antechamber.http://webapps.ccpn.ac.uk/acpype/ I have read a lot more time that charge on different atoms are not correct when created through various softwares/scripts. My question is whether ACPYPE/antechamber also fall in this category? If yes, can anyone please provide me with any link from where I can get charge values to assign? Thank you.mayaz -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/Search before posting! * Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
Re: [gmx-users] Ligand charge issues
In https://code.google.com/p/acpype/ you can look the wikis and you see the explanations about the partial charges. The best solution, though not straightforward, would be using http://q4md-forcefieldtools.org/REDS/ Alan On 2 September 2013 11:27, Muhammad Ayaz Anwar wrote: > Hi Gromacs users, > I am studying the protein-ligand interaction using amber99sb-ILDN force > field in gromacs 4.6.2. To create the ligand topology (lipid A), I have > used online version of ACPYPE/antechamber. > http://webapps.ccpn.ac.uk/acpype/ > I have read a lot more time that charge on different atoms are not correct > when created through various softwares/scripts. My question is whether > ACPYPE/antechamber also fall in this category? If yes, can anyone please > provide me with any link from where I can get charge values to assign? > Thank you.mayaz > -- > gmx-users mailing listgmx-users@gromacs.org > http://lists.gromacs.org/mailman/listinfo/gmx-users > * Please search the archive at > http://www.gromacs.org/Support/Mailing_Lists/Search before posting! > * Please don't post (un)subscribe requests to the list. Use the > www interface or send it to gmx-users-requ...@gromacs.org. > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists > -- Alan Wilter SOUSA da SILVA, DSc Bioinformatician, UniProt European Bioinformatics Institute (EMBL-EBI) European Molecular Biology Laboratory Wellcome Trust Genome Campus Hinxton Cambridge CB10 1SD United Kingdom Tel: +44 (0)1223 494588 -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/Search before posting! * Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
[gmx-users] adding ETHANOL to CGennFF in gromacs has problem
Hello, I am trying to add ethanol as a new residue to CGenFF force field in gromacs. Basically, I need to bring this information from charmm rtf file to the gromacs format. This is the information taken from rtf file of charmm. RESI ETOH 0.00 ! C2H6O Ethanol, adm jr. GROUP ATOM C1 CG321 0.05 ! H21 H11 H12 ATOM O1 OG311 -0.65 ! \ \ / ATOM HO1 HGP1 0.42 ! H22--C2--C1 ATOM H11 HGA2 0.09 ! / \ ATOM H12 HGA2 0.09 ! H23 O1--HO1 GROUP ATOM C2 CG331 -0.27 ATOM H21 HGA3 0.09 ATOM H22 HGA3 0.09 ATOM H23 HGA3 0.09 BOND C1 C2 C1 O1 C1 H11 C1 H12 O1 HO1 BOND C2 H21 C2 H22 C2 H23 DONO HO1 O1 ACCE O1 ! for ic build IC O1 C1 C2 H21 0. 0. 180. 0. 0. IC C1 H21 *C2 H22 0. 0. 120. 0. 0. IC C1 H21 *C2 H23 0. 0. 240. 0. 0. IC O1 C2 *C1 H11 0. 0. 240. 0. 0. IC O1 C2 *C1 H12 0. 0. 120. 0. 0. IC C2 C1 O1 HO1 0. 0. 180. 0. 0. I have added The ATOM and BOND information to the aminoacid.rtp file. I thought that the rest of the information are optional and pdb2gmx doesn't complain while generating the topology. But when I use the topology and gro file to do a minimization. The molecule bonds break apart and what I see is an exploded molecule. So I thought probably the rest of the information are important too, for example: DONO HO1 O1 ACCE O1 But I dont know where to put them. Thanks G. -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/Search before posting! * Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
Re: [gmx-users] adding ETHANOL to CGennFF in gromacs has problem
On 9/2/13 12:59 PM, Golshan Hejazi wrote: Hello, I am trying to add ethanol as a new residue to CGenFF force field in gromacs. Basically, I need to bring this information from charmm rtf file to the gromacs format. This is the information taken from rtf file of charmm. RESI ETOH 0.00 ! C2H6O Ethanol, adm jr. GROUP ATOM C1 CG3210.05 ! H21 H11 H12 ATOM O1 OG311 -0.65 ! \ \ / ATOM HO1 HGP1 0.42 ! H22--C2--C1 ATOM H11 HGA2 0.09 ! / \ ATOM H12 HGA2 0.09 ! H23O1--HO1 GROUP ATOM C2 CG331 -0.27 ATOM H21 HGA3 0.09 ATOM H22 HGA3 0.09 ATOM H23 HGA3 0.09 BOND C1 C2 C1 O1 C1 H11 C1 H12 O1 HO1 BOND C2 H21 C2 H22 C2 H23 DONO HO1 O1 ACCE O1 ! for ic build IC O1 C1 C2 H21 0. 0. 180. 0. 0. IC C1 H21 *C2 H22 0. 0. 120. 0. 0. IC C1 H21 *C2 H23 0. 0. 240. 0. 0. IC O1 C2 *C1 H11 0. 0. 240. 0. 0. IC O1 C2 *C1 H12 0. 0. 120. 0. 0. IC C2 C1 O1 HO1 0. 0. 180. 0. 0. I have added The ATOM and BOND information to the aminoacid.rtp file. I thought that the rest of the information are optional and pdb2gmx doesn't complain while generating the topology. But when I use the topology and gro file to do a minimization. The molecule bonds break apart and what I see is an exploded molecule. Please post your complete .rtp entry and the contents of the resulting .top file. So I thought probably the rest of the information are important too, for example: DONO HO1 O1 ACCE O1 But I dont know where to put them. These lines are only used by CHARMM during hydrogen bond analysis or in the archaic form that included a hydrogen bonding term to the energy function. -Justin -- == Justin A. Lemkul, Ph.D. Postdoctoral Fellow Department of Pharmaceutical Sciences School of Pharmacy Health Sciences Facility II, Room 601 University of Maryland, Baltimore 20 Penn St. Baltimore, MD 21201 jalem...@outerbanks.umaryland.edu | (410) 706-7441 == -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/Search before posting! * Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
Re: [gmx-users] adding ETHANOL to CGennFF in gromacs has problem
Please make sure to keep the discussion on the mailing list. On 9/2/13 1:27 PM, Golshan Hejazi wrote: This is the gro file that I am using: ETHANOL 9 1ETHAC11 3.133 3.266 0.008 1ETHAO12 3.244 3.350 -0.021 1ETHA HO13 3.324 3.296 -0.015 1ETHA H114 3.143 3.225 0.111 1ETHA H125 3.129 3.182 -0.065 1ETHAC26 3.007 3.351 -0.002 1ETHA H217 3.011 3.435 0.071 1ETHA H228 2.998 3.394 -0.104 1ETHA H239 2.917 3.290 0.019 3.41840 3.44350 3.38360 This is what I added in the rtp: [ ETHA ] [ atoms ] C1 CG321 0.050 O1 OG311 -0.65 1 HO1 HGP10.422 H11 HGA20.093 H12 HGA20.094 C2 CG331 -0.27 5 H21 HGA30.096 H22 HGA30.097 H23 HGA30.098 [ bonds ] C1 C2 C1 O1 C1 H11 C1 H12 O1 HO1 C2 H21 C2 H22 C2 H23 And this is the resulting topology: ; Include forcefield parameters #include "./charmm36.ff/forcefield.itp" [ moleculetype ] ; Namenrexcl drug3 [ atoms ] ; nr type resnr residue atom cgnr charge mass typeB chargeB massB ; residue 1 ETHA rtp ETHA q 0.0 1 CG321 1 ETHA C1 1 0.05 12.011 ; qtot 0.05 2 OG311 1 ETHA O1 2 -0.6515.9994 ; qtot -0.6 3 HGP1 1 ETHAHO1 3 0.42 1.008 ; qtot -0.18 4 HGA2 1 ETHAH11 4 0.09 1.008 ; qtot -0.09 5 HGA2 1 ETHAH12 5 0.09 1.008 ; qtot 0 6 CG331 1 ETHA C2 6 -0.27 12.011 ; qtot -0.27 7 HGA3 1 ETHAH21 7 0.09 1.008 ; qtot -0.18 8 HGA3 1 ETHAH22 8 0.09 1.008 ; qtot -0.09 9 HGA3 1 ETHAH23 9 0.09 1.008 ; qtot 0 [ bonds ] ; aiaj functc0c1c2c3 1 2 1 1 4 1 1 5 1 1 6 1 2 3 1 6 7 1 6 8 1 6 9 1 [ pairs ] ; aiaj functc0c1c2c3 2 7 1 2 8 1 2 9 1 3 4 1 3 5 1 3 6 1 4 7 1 4 8 1 4 9 1 5 7 1 5 8 1 5 9 1 [ angles ] ; aiajak functc0c1c2c3 2 1 4 5 2 1 5 5 2 1 6 5 4 1 5 5 4 1 6 5 5 1 6 5 1 2 3 5 1 6 7 5 1 6 8 5 1 6 9 5 7 6 8 5 7 6 9 5 8 6 9 5 [ dihedrals ] ; aiajakal functc0c1c2 c3c4c5 4 1 2 3 9 5 1 2 3 9 6 1 2 3 9 2 1 6 7 9 2 1 6 8 9 2 1 6 9 9 4 1 6 7 9 4 1 6 8 9 4 1 6 9 9 5 1 6 7 9 5 1 6 8 9 5 1 6 9 9 ; Include Position restraint file #ifdef POSRES #include "posre.itp" #endif [ system ] ; Name ETHANOL [ molecules ] ; Compound#mols drug1 All of this looks perfectly normal. What are you trying to do when it flies apart? Are you dealing with a single molecule only? -Justin -- == Justin A. Lemkul, Ph.D. Postdoctoral Fellow Department of Pharmaceutical Sciences School of Pharmacy Health Sciences Facility II, Room 601 University of Maryland, Baltimore 20 Penn St. Baltimore, MD 21201 jalem...@outerbanks.umaryland.edu | (410) 706-7441 == -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/Search before posting! * Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
Re: [gmx-users] adding ETHANOL to CGennFF in gromacs has problem
Thanks Justin for looking at the files. I am using the gro and top file to generate a tpr with the following mdp file in order to see if everything is all right: title = geometry-optimization cpp = /usr/bin/cpp include = integrator = steep comm_mode = Linear dt = 0.001 nsteps = 100 ; Output frequency for coords (x), velocities (v) and forces (f) = nstxout = 1000 nstvout = 1000 ; Output frequency for energies to log file and energy file = nstlog = 1000 nstenergy = 1000 nstxtcout = 1000 xtc-precision = 10 xtc_grps = System energygrps = System pbc = xyz nstlist = 10 epsilon_r = 1. ns_type = grid ; or grid I don't know coulombtype = pme vdwtype = Cut-off fourierspacing = 0.12 ; EWALD/PME/PPPM parameters pme_order = 4 ewald_rtol = 1e-05 epsilon_surface = 0 optimize_fft = yes rlist = 0.8 rcoulomb = 0.8 rvdw = 0.8 emtol = 1.0 Yes, I am dealing with single molecule of ethanol. I can attach the tpr file and the resulting trr and gro file but it seems that the file size is big. When I visualize the trajectory. in the second step, the molecule is already without any bond ... it is separated atoms. I don't think this is visualization problem since I already tried all trjconv -pbc option. I am very confused. Thank for helping G. From: Justin Lemkul To: Discussion list for GROMACS users Sent: Monday, September 2, 2013 1:34 PM Subject: Re: [gmx-users] adding ETHANOL to CGennFF in gromacs has problem Please make sure to keep the discussion on the mailing list. On 9/2/13 1:27 PM, Golshan Hejazi wrote: > This is the gro file that I am using: > ETHANOL > 9 > 1ETHA C1 1 3.133 3.266 0.008 > 1ETHA O1 2 3.244 3.350 -0.021 > 1ETHA HO1 3 3.324 3.296 -0.015 > 1ETHA H11 4 3.143 3.225 0.111 > 1ETHA H12 5 3.129 3.182 -0.065 > 1ETHA C2 6 3.007 3.351 -0.002 > 1ETHA H21 7 3.011 3.435 0.071 > 1ETHA H22 8 2.998 3.394 -0.104 > 1ETHA H23 9 2.917 3.290 0.019 > 3.41840 3.44350 3.38360 > > This is what I added in the rtp: > > [ ETHA ] > [ atoms ] > C1 CG321 0.05 0 > O1 OG311 -0.65 1 > HO1 HGP1 0.42 2 > H11 HGA2 0.09 3 > H12 HGA2 0.09 4 > C2 CG331 -0.27 5 > H21 HGA3 0.09 6 > H22 HGA3 0.09 7 > H23 HGA3 0.09 8 > [ bonds ] > C1 C2 > C1 O1 > C1 H11 > C1 H12 > O1 HO1 > C2 H21 > C2 H22 > C2 H23 > > And this is the resulting topology: > > ; Include forcefield parameters > #include "./charmm36.ff/forcefield.itp" > > [ moleculetype ] > ; Name nrexcl > drug 3 > > [ atoms ] > ; nr type resnr residue atom cgnr charge mass typeB > chargeB massB > ; residue 1 ETHA rtp ETHA q 0.0 > 1 CG321 1 ETHA C1 1 0.05 12.011 ; qtot >0.05 > 2 OG311 1 ETHA O1 2 -0.65 15.9994 ; qtot >-0.6 > 3 HGP1 1 ETHA HO1 3 0.42 1.008 ; qtot >-0.18 > 4 HGA2 1 ETHA H11 4 0.09 1.008 ; qtot >-0.09 > 5 HGA2 1 ETHA H12 5 0.09 1.008 ; qtot >0 > 6 CG331 1 ETHA C2 6 -0.27 12.011 ; qtot >-0.27 > 7 HGA3 1 ETHA H21 7 0.09 1.008 ; qtot >-0.18 > 8 HGA3 1 ETHA H22 8 0.09 1.008 ; qtot >-0.09 > 9 HGA3 1 ETHA H23 9 0.09 1.008 ; qtot >0 > [ bonds ] > ; ai aj funct c0 c1 c2 c3 > 1 2 1 > 1 4 1 > 1 5 1 > 1 6 1 > 2 3 1 > 6 7 1 > 6 8 1 > 6 9 1 > > [ pairs ] > ; ai aj funct c0 c1 c2 c3 > 2 7 1 > 2 8 1 > 2 9 1 > 3 4 1 > 3 5 1 > 3 6 1 > 4 7 1 > 4 8 1 > 4 9 1 > 5 7 1 > 5 8 1 > 5 9 1 > [ angles ] > ; ai aj ak funct c0 c1 c2 > c3 > 2 1 4 5 > 2 1 5
Re: [gmx-users] adding ETHANOL to CGennFF in gromacs has problem
On 9/2/13 2:02 PM, Golshan Hejazi wrote: Thanks Justin for looking at the files. I am using the gro and top file to generate a tpr with the following mdp file in order to see if everything is all right: title= geometry-optimization cpp = /usr/bin/cpp include = integrator = steep comm_mode= Linear dt = 0.001 nsteps = 100 ; Output frequency for coords (x), velocities (v) and forces (f) = nstxout = 1000 nstvout = 1000 ; Output frequency for energies to log file and energy file = nstlog = 1000 nstenergy= 1000 nstxtcout= 1000 xtc-precision= 10 xtc_grps = System energygrps = System pbc = xyz nstlist = 10 epsilon_r= 1. ns_type = grid ; or grid I don't know coulombtype = pme vdwtype = Cut-off fourierspacing = 0.12 ; EWALD/PME/PPPM parameters pme_order= 4 ewald_rtol = 1e-05 epsilon_surface = 0 optimize_fft = yes rlist= 0.8 rcoulomb = 0.8 rvdw = 0.8 emtol= 1.0 Note that your cutoffs are incorrect for strict adherence to CHARMM settings. See extensive threads from recent days and weeks on this topic. Yes, I am dealing with single molecule of ethanol. I can attach the tpr file and the resulting trr and gro file but it seems that the file size is big. The mailing list does not accept attachments, nor are those files particularly informative (maybe the .gro, but that can easily just be pasted in an email for a single-molecule system). Screenshots linked from some publicly accessible image-sharing site work well. When I visualize the trajectory. in the second step, the molecule is already without any bond ... it is separated atoms. I don't think this is visualization problem since I already tried all trjconv -pbc option. Measurements of intramolecular distances are more informative than whether or not visualization software thinks there's actually a bond there. The topology is definitive, visualization is not. There may certainly be a problem, but be sure you're using quantitative assessments, not qualitative (and error-prone). -Justin -- == Justin A. Lemkul, Ph.D. Postdoctoral Fellow Department of Pharmaceutical Sciences School of Pharmacy Health Sciences Facility II, Room 601 University of Maryland, Baltimore 20 Penn St. Baltimore, MD 21201 jalem...@outerbanks.umaryland.edu | (410) 706-7441 == -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/Search before posting! * Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
Re: [gmx-users] Distance restraints exploding system
Have you tried with even less restraints? I found systems are not always stable with more than the bare minimum of restraints On Sun, Sep 1, 2013 at 7:54 PM, Trayder Thomas wrote: > It still explodes on a 0.1fs timestep, turning off P-R doesn't seem to have > an impact. I've tried being gentle, slowly turning up the force constant > and running for 1 ns for each value but as soon as the force constant > approaches 100 it crashes. > The starting structure was generated with the same restraints, so it is > very close. I have tried using slightly different starting structures as > well. > > I've tried running it with only 3 restraints (1 methyl group to 1 hydrogen > with the restraints extended by 0.2nm so that all hydrogens are within the > restraint distance) and I'm getting a segmentation fault (no LINCS > warnings) at 100 steps. It runs fine with any combination of 2 hydrogens > restrained, but as soon as I restrain the 3rd one I get a segmentation > fault (this occurs with either setting for disre-weighting). So it seems to > fair better with more restraints? > > The more I try to solve this problem the less it makes sense! > > -Trayder > > > On Fri, Aug 30, 2013 at 6:02 PM, Mark Abraham >wrote: > > > On Fri, Aug 30, 2013 at 8:56 AM, Trayder > > wrote: > > > Hello, > > > I am attempting to simulate a protein-ligand complex using distance > > > restraints to match it to NMR data. > > > The system runs stably without restraints. With restraints it tends to > > spit > > > out LINCS angle warnings and blow up under most conditions. > > > > > > I'm attempting to use: > > > ; Restraints > > > disre = simple > > > disre-weighting = conservative > > > disre-fc= 1000 > > > > > > It blows up within 100 steps unless: > > > I run on a single core (+gpu) or > > > disre-fc = <100 or > > > disre-weighting = equal > > > > > > If disre-weighting = conservative is causing extreme forces, then I > > figure > > > it should do the same on 1 core. > > > > Not really. MD is chaotic. Small changes in initial conditions lead to > > different results. > > > > > If domain decomposition is the problem, then I would think > > disre-weighting = > > > equal shouldn't work either. > > > I'm stumped... anyone got any ideas? > > > > http://www.gromacs.org/Documentation/Terminology/Blowing_Up has the > > usual suggestions - don't use P-R yet, try a smaller time step, make > > sure your system is close to the restrained regime (or be extra gentle > > until it is). > > > > Mark > > > > > Thanks in advance, > > > -Trayder > > > > > > Distance restraints excerpt: > > > ; aiaj typeindex type’ low up1 up2 fac > > > ; 2 symmetric hydrogens > > > 1306 1389 1 10 1 0.0 0.548 1.0 1.0 > > > 1306 1396 1 10 1 0.0 0.548 1.0 1.0 > > > ; Diastereotopic methyl groups > > > 1306 1374 1 11 1 0.0 0.654 1.0 1.0 > > > 1306 1375 1 11 1 0.0 0.654 1.0 1.0 > > > 1306 1376 1 11 1 0.0 0.654 1.0 1.0 > > > 1306 1385 1 11 1 0.0 0.654 1.0 1.0 > > > 1306 1386 1 11 1 0.0 0.654 1.0 1.0 > > > 1306 1387 1 11 1 0.0 0.654 1.0 1.0 > > > > > > Full mdp: > > > ; Run Control > > > integrator = md ; simulation algorithm > > > tinit= 0 > > > dt = 0.002 > > > nsteps = 50 > > > ; > > > ; Output Control > > > nstxout = 20; write coordinates to > > .trr > > > nstvout = 20; write velocities to > > .trr > > > nstlog = 1000 ; write energies to > .log > > > nstenergy = 4000 ; write energies to > .edr > > > nstxtcout = 1000 ; write coordinates to > > .xtc > > > ; > > > ; Neighbour Searching > > > nstlist = 10 ; update neighbour list > > > ns_type = grid ; neighbour list method > > > pbc = xyz ; periodic boundary > > > conditions > > > rlist = 0.9 ; cut-off for > short-range > > > neighbour (nm) > > > cutoff-scheme = verlet > > > ; > > > ; Electrostatics and VdW > > > coulombtype = PME ; type of coulomb > > > interaction > > > rcoulomb= 0.9 ; cut-off distance for > > > coulomb > > > epsilon_r = 1; dielectric constant > > > rvdw= 0.9 ; cut-off for vdw > > > fourierspacing = 0.12 ; maximum grid spacing > > for > > > FFT > > > pme_order = 4; interpolation order > f
Re: [gmx-users] Long range Lennard Jones
On Thu, Aug 29, 2013 at 7:18 AM, Gianluca Interlandi wrote: > Justin, > > I respect your opinion on this. However, in the paper indicated below by BR > Brooks they used a cutoff of 10 A on LJ when testing IPS in CHARMM: > > Title: Pressure-based long-range correction for Lennard-Jones interactions > in molecular dynamics simulations: Application to alkanes and interfaces > Author(s): Lague, P; Pastor, RW; Brooks, BR > Source: JOURNAL OF PHYSICAL CHEMISTRY B Volume: 108 Issue: 1 Pages: 363-368 > Published: JAN 8 2004 > > There is also a paper by Piana and Shaw where different cutoffs for > non-bonded are tested with CHARMM22 on Anton: > > http://www.plosone.org/article/info%3Adoi%2F10.1371%2Fjournal.pone.0039918 > > They found some subtle differences, in particular for cutoffs shorter than 9 > A. However, Anton uses abrupt truncation (no switching) and I believe that > the differences they found at cutoffs > 9 A would be much smaller if they > had used a finer mesh (as they show at the 8 A cutoff). I always use > fourierspacing=1.0 > > I agree though that it strongly depends on the system and I have always run > control simulations but never found significant differences in the case of > just proteins. > > Finally, I have not tested it in gromacs, but in NAMD there is a performance > gain of 25% when using the shorter cutoff. This is a factor to consider. I'm not sure about what "the sorter cutoff" refers to, neither about the context Let me pitch in with little bit of clarification on computational costs, I hope somebody can make use of it. Note that I don't want to argue for the use of shorter cut-offs - although revising a "magic simulation recipe" one used just because that's what one was provided may be advantageous anyway. The performance gain from a shorter cut-off depends on several factors (and most of these factors are _not_ GROMACS-specific): 1) Performance of the short-range and (vs) long-range interaction algorithm/implementation; 2) Scaling of short-range and long-range parts of the code (a given implementation on a given hardware, note that not all implementations are equal!); 3) If some sort of task-parallelization is used, e.g. accelerator offloading (GPUs) or "PME nodes" in GROMACS, depending on the implementation, a longer cut-off can sometimes be free lunch or it can even improve performance. For instance, with GROMACS 4.6 if you have a very fast GPU (wrt the CPU) you may get the same or better performance with a longer cut-off. 4) Single or twin cut-off (rcoul = rvdw?). (+probably a few more that I forgot about :) Given the fact that e.g. a 1.2 nm cut-off will result in a 1.73x larger interaction volume than 1.0 nm, this will obviously drive up the cost of short-range interactions. Hence, in this case, 1 nm cut-off will in most cases give faster simulations. As mentioned above, exception can be cases when, due to imbalanced PP-PME load, reducing the PME cost is advantageous e.g. fast GPU case or at high parallelization. -- Szilárd > When I asked for Teragrid supercomputing allocations back in 2006 and 2007 > and I suggested 10/12/14 cutoff, the reviewers always complained and cut my > requested time of 20% with the justification that I must use a shorter > cutoff. > > Gianluca > > > On Wed, 28 Aug 2013, Justin Lemkul wrote: > >> On 8/28/13 7:28 PM, Gianluca Interlandi wrote: >>> >>> Thanks for your replies, Mark. What do you think about the current >>> DispCorr >>> option in gromacs? Is it worth it trying it? Also, I wonder whether using >>> DispCorr for LJ + PME for Cb justifies reducing the cutoff for non-bonded >>> to 1 >>> nm with the CHARMM force field, where 1.2 nm is usually recommended. >>> >> >> This is risky. Current CHARMM development relies on a 1.2-nm cutoff for >> LJ, so that's how we balance all of the forces during parameterization. To >> me, ad hoc changes like these are not worth the tiny (potential) increase in >> performance. As I recently told someone else on this topic, if you're intent >> on fiddling with the typical workings of a force field, especially if you're >> making changes to something so fundamental, be prepared to undertake a >> demonstration that you can recapitulate all of the expected outcomes of the >> force field or improve upon them. My gut feeling, in this case and others, >> is that you won't be able to. You're messing with something that is fairly >> critical to obtaining sensible results. >> >> As for dispersion correction, it is generally helpful, but it assumes a >> homogeneous environment. If you simulate with a membrane, for instance, >> this assumption breaks down, though some literature suggests that use of >> dispersion correction in these cases is still better than nothing. >> >> -Justin >> >> -- >> == >> >> Justin A. Lemkul, Ph.D. >> Postdoctoral Fellow >> >> Department of Pharmaceutical Sciences >> School of Pharmacy >> Health Sciences Facility II, Room 601 >> University of
[gmx-users] Re: MD vs. free energy simulations
Dear Justin, I set the couple_intramol parameter to yes and rerun the free energy simulations. mdrun was able to fully utilize all the cores in the system but there's one issue. The free energy of solvation is vastly different if the parameter is set to 'no' (DG=-408.10 +/- 8.69 kJ/mol) or 'yes' (1890.14 +/- 15.52). The behavior observed in the simulations is more consistent with the latter number. The manual states that 'yes' can be used can be "useful for partitioning free-energies of relatively large molecules, where the intra-molecular non-bonded interactions might lead to kinetically trapped vacuum conformations. The 1-4 pair interactions are not turned off." My molecule has 161 atoms. How large is "relatively large" ? Thanks in advance, Jernej Zidar On Thu, Aug 29, 2013 at 8:00 PM, wrote: > This issue is discussed very frequently, particularly in the context of free > energy calculations. Please refer to the list archive and consider the effect > that couple_intramol setting has on the DD setup. -- Windows: Re-Boot, Linux: Be-Root. -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/Search before posting! * Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
Re: [gmx-users] Re: MD vs. free energy simulations
On 9/2/13 9:43 PM, Jernej Zidar wrote: Dear Justin, I set the couple_intramol parameter to yes and rerun the free energy simulations. mdrun was able to fully utilize all the cores in the system but there's one issue. The free energy of solvation is vastly different if the parameter is set to 'no' (DG=-408.10 +/- 8.69 kJ/mol) or 'yes' (1890.14 +/- 15.52). The behavior observed in the simulations is more consistent with the latter number. The latter number, of course, is not an actual final answer. You have to do a corresponding gas-phase calculation to determine the real DG of solvation. The manual states that 'yes' can be used can be "useful for partitioning free-energies of relatively large molecules, where the intra-molecular non-bonded interactions might lead to kinetically trapped vacuum conformations. The 1-4 pair interactions are not turned off." My molecule has 161 atoms. How large is "relatively large" ? I think a bigger factor is flexibility. For a relatively rigid, conjugated system, it may not be a big issue. Gut feeling? Yes, I'd say that something with 161 atoms is relatively large. -Justin -- == Justin A. Lemkul, Ph.D. Postdoctoral Fellow Department of Pharmaceutical Sciences School of Pharmacy Health Sciences Facility II, Room 601 University of Maryland, Baltimore 20 Penn St. Baltimore, MD 21201 jalem...@outerbanks.umaryland.edu | (410) 706-7441 == -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/Search before posting! * Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
