Re: [gmx-users] add a group to an amino acid
Dear Mark, I think that i couldn't explain my problem clearly, however i have changed my .rtp but i don't know how should be my .hdb file. a) The structure of *BOC* is: *-CO-CH2-CH(C6H5Cl)-CH2-NH-* b) This is the new charges: [ BOC ] [ atoms ] Nopls_238 -0.500 1 Hopls_2410.300 1 CAopls_224B 0.080 1 HA1opls_1400.060 1 HA2opls_1400.060 1 CBopls_1490.055 2 HBopls_1400.060 2 CG1opls_145 -0.115 2 CD1opls_145 -0.115 3 HD1opls_1460.115 3 CD2opls_145 -0.115 4 HD2opls_1460.115 4 CE1opls_145 -0.115 5 HE1opls_1460.115 5 CE2opls_145 -0.115 6 HE2opls_1460.115 6 CZopls_1451.000 7 Clopls_264 -1.000 7 CG2opls_071 -0.120 8 HG1opls_1400.060 8 HG2opls_1400.060 8 Copls_2350.700 9 Oopls_236 -0.700 9 c) charge for CG1 in PHE (in OPLS .rtp) is negative, too. d) Would you explain for me that what is wrong in my atom types? e) As you see in the error, at first there are 63 residues (like my pdb) but when when it adds COO- and NH3+ there are 69 residues! it should add only one O (for COO-) and two Hs for (NH3+) i think. There are 1 chains and 0 blocks of water and 3 residues with 63 atoms chain #res #atoms 1 ' ' 3 63 All occupancy fields zero. This is probably not an X-Ray structure Opening library file /usr/share/gromacs/top/ffoplsaa.atp Atomtype 1 Reading residue database... (ffoplsaa) Opening library file /usr/share/gromacs/top/ffoplsaa.rtp Residue 57 Sorting it all out... Opening library file /usr/share/gromacs/top/ffoplsaa.hdb Opening library file /usr/share/gromacs/top/ffoplsaa-n.tdb Opening library file /usr/share/gromacs/top/ffoplsaa-c.tdb Back Off! I just backed up topol.top to ./#topol.top.3# Processing chain 1 (63 atoms, 3 residues) There are 3 donors and 3 acceptors There are 4 hydrogen bonds Checking for duplicate atoms Opening library file /usr/share/gromacs/top/specbond.dat 7 out of 7 lines of specbond.dat converted succesfully N-terminus: NH3+ C-terminus: COO- Now there are 3 residues with 69 atoms Making bonds... Warning: Long Bond (53-54 = 0.334847 nm) Warning: Long Bond (53-56 = 0.334912 nm) Warning: Long Bond (64-65 = 0.271173 nm) Warning: Long Bond (64-66 = 0.347808 nm) Warning: Long Bond (64-67 = 0.31288 nm) Opening library file /usr/share/gromacs/top/aminoacids.dat --- Program pdb2gmx, VERSION 4.0.7 Source code file: ../../../../src/kernel/add_par.c, line: 233 Fatal error: Atom HB1 not found in rtp database in residue BOC, it looks a bit like HB -- 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] add a group to an amino acid
On 14/11/2010 7:26 PM, hengame fallah wrote: Dear Mark, I think that i couldn't explain my problem clearly, however i have changed my .rtp but i don't know how should be my .hdb file. a) The structure of *_BOC_* is: *-CO-CH2-CH(C6H5Cl)-CH2-NH-* Great. b) This is the new charges: [ BOC ] [ atoms ] Nopls_238 -0.500 1 Hopls_2410.300 1 CAopls_224B 0.080 1 HA1opls_1400.060 1 HA2opls_1400.060 1 CBopls_1490.055 2 HBopls_1400.060 2 CG1opls_145 -0.115 2 CD1opls_145 -0.115 3 HD1opls_1460.115 3 CD2opls_145 -0.115 4 HD2opls_1460.115 4 CE1opls_145 -0.115 5 HE1opls_1460.115 5 CE2opls_145 -0.115 6 HE2opls_1460.115 6 CZopls_1451.000 7 Clopls_264 -1.000 7 CG2opls_071 -0.120 8 HG1opls_1400.060 8 HG2opls_1400.060 8 Copls_2350.700 9 Oopls_236 -0.700 9 This looks much better. c) charge for CG1 in PHE (in OPLS .rtp) is negative, too. OK d) Would you explain for me that what is wrong in my atom types? No. My point last email was that you may need to consider what atom types are suitable for your backbone methylene C. I don't know what is suitable. You will have to look at the .atp file and perhaps read OPLS/AA papers. e) As you see in the error, at first there are 63 residues (like my pdb) but when when it adds COO- and NH3+ there are 69 residues! it should add only one O (for COO-) and two Hs for (NH3+) i think. Well, what do the .tdb files say? You should model these for BOC upon the PHE ones, of course. You will need to have read the parts of chapter 5 that pertain to the .hdb and .tdb files. There are 1 chains and 0 blocks of water and 3 residues with 63 atoms chain #res #atoms 1 ' ' 3 63 All occupancy fields zero. This is probably not an X-Ray structure Opening library file /usr/share/gromacs/top/ffoplsaa.atp Atomtype 1 Reading residue database... (ffoplsaa) Opening library file /usr/share/gromacs/top/ffoplsaa.rtp Residue 57 Sorting it all out... Opening library file /usr/share/gromacs/top/ffoplsaa.hdb Opening library file /usr/share/gromacs/top/ffoplsaa-n.tdb Opening library file /usr/share/gromacs/top/ffoplsaa-c.tdb Back Off! I just backed up topol.top to ./#topol.top.3# Processing chain 1 (63 atoms, 3 residues) There are 3 donors and 3 acceptors There are 4 hydrogen bonds Checking for duplicate atoms Opening library file /usr/share/gromacs/top/specbond.dat 7 out of 7 lines of specbond.dat converted succesfully N-terminus: NH3+ C-terminus: COO- Now there are 3 residues with 69 atoms Making bonds... Warning: Long Bond (53-54 = 0.334847 nm) Warning: Long Bond (53-56 = 0.334912 nm) Warning: Long Bond (64-65 = 0.271173 nm) Warning: Long Bond (64-66 = 0.347808 nm) Warning: Long Bond (64-67 = 0.31288 nm) These warnings are still suggestive that your atom coordinates are labelled inappropriately in your coordinate file. Opening library file /usr/share/gromacs/top/aminoacids.dat --- Program pdb2gmx, VERSION 4.0.7 Source code file: ../../../../src/kernel/add_par.c, line: 233 Fatal error: Atom HB1 not found in rtp database in residue BOC, it looks a bit like HB I think your .hdb is still trying to add two hydrogens to CB. Mark -- 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] water molecules (interfacial) - IIRC
Dear Mark I have a question as same as atila petrosian, but my system is protein-dna. you said [IIRC g_select has some better documentation available once you run it and ask for help, somewhat like make_ndx]. please explain IIRC more. Is IIRC a program or a list of archives? any help will highly appreciated. -- 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] water molecules (interfacial) - IIRC
On 14/11/2010 7:57 PM, leila karami wrote: Dear Mark I have a question as same as atila petrosian, but my system is protein-dna. you said [IIRC g_select has some better documentation available once you run it and ask for help, somewhat like make_ndx]. please explain IIRC more. Is IIRC a program or a list of archives? IIRC is a standard internet acronym for if I recall correctly In fact, g_select -select help gives access to its online help. Mark -- 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] add a group to an amino acid
I read about .hdb file and then i editted my .hdb file: ... BOC 9 11H N -CCA 26HA CA N CB 15HB CB CA CG1CG2 11HD1CD1CG1CE1 11HD2CD2CG1CE2 11HE1CE1CD1CZ 11HE2CE2CD2CZ 26HG CG2CB C CYS23 11H N-C CA 15HACAN CCB 26HBCBSGCA ... PHE 8 11H N -C CA 15HA CA N CCB 26HB CB CGCA 11HD1CD1CGCE1 11HD2CD2CGCE2 11HE1CE1CD1CZ 11HE2CE2CD2CZ 11HZ CZ CE1CE2 ... but i got this error: ... There are 1 chains and 0 blocks of water and 3 residues with 63 atoms chain #res #atoms 1 ' ' 3 63 All occupancy fields zero. This is probably not an X-Ray structure Opening library file /usr/share/gromacs/top/ffoplsaa.atp Atomtype 1 Reading residue database... (ffoplsaa) Opening library file /usr/share/gromacs/top/ffoplsaa.rtp Residue 57 Sorting it all out... Opening library file /usr/share/gromacs/top/ffoplsaa.hdb --- Program pdb2gmx, VERSION 4.0.7 Source code file: ../../../../src/kernel/h_db.c, line: 87 Fatal error: wrong format in input file ffoplsaa.hdb on line CYS23 -- 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] add a group to an amino acid
On 14/11/2010 8:13 PM, hengame fallah wrote: I read about .hdb file and then i editted my .hdb file: ... BOC 9 Then you'll know from your reading what this number should be :-) It's the number of types of hydrogen atoms that might be added. For each of them, there's one line below. You needed to have left this on 8. Because you've changed it to 9, pdb2gmx tries to make sense of a ninth line, which is CYS2 3 and gives you the error you saw. Mark 11H N -CCA 26HA CA N CB 15HB CB CA CG1CG2 11HD1CD1CG1CE1 11HD2CD2CG1CE2 11HE1CE1CD1CZ 11HE2CE2CD2CZ 26HG CG2CB C CYS23 11H N-C CA 15HACAN CCB 26HBCBSGCA ... PHE 8 11H N -C CA 15HA CA N CCB 26HB CB CGCA 11HD1CD1CGCE1 11HD2CD2CGCE2 11HE1CE1CD1CZ 11HE2CE2CD2CZ 11HZ CZ CE1CE2 ... but i got this error: ... There are 1 chains and 0 blocks of water and 3 residues with 63 atoms chain #res #atoms 1 ' ' 3 63 All occupancy fields zero. This is probably not an X-Ray structure Opening library file /usr/share/gromacs/top/ffoplsaa.atp Atomtype 1 Reading residue database... (ffoplsaa) Opening library file /usr/share/gromacs/top/ffoplsaa.rtp Residue 57 Sorting it all out... Opening library file /usr/share/gromacs/top/ffoplsaa.hdb --- Program pdb2gmx, VERSION 4.0.7 Source code file: ../../../../src/kernel/h_db.c, line: 87 Fatal error: wrong format in input file ffoplsaa.hdb on line CYS23 -- 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] 1-4 interaction, a distance greater than table size
Amin Arabbagheri wrote: Dear All, I'm simulating a DNA duplex using amber99p. starting the simulation(step 0) i face an error which tells: starting mdrun 'Protein in water' 50 steps,500.0 ps. step 0Warning: 1-4 interaction between 136 and 179 at distance 96871502.536 which is larger than the 1-4 table size 2.000 nm These are ignored for the rest of the simulation This usually means your system is exploding, if not, you should increase table-extension in your mdp file or with user tables increase the table size Segmentation fault I've checked the distance between atoms 136 and 179, its quite normal, about 3 angstrom (not 96871502 !!). Im also attaching the output (md.log) which may help finding a reason: Started mdrun on node 0 Sun Nov 14 10:17:10 2010 Step Time Lambda 00.00.0 Grid: 9 x 9 x 9 cells Energies (kJ/mol) Bond AngleProper Dih. Ryckaert-Bell. LJ-14 1.77478e+026.08613e+021.86464e+021.58123e+036.34721e+02 Coulomb-14LJ (SR) Coulomb (SR) Coul. recip. Potential -1.03292e+044.16184e+04 -3.45397e+05 -3.54920e+04 -3.46412e+05 Kinetic En. Total EnergyTemperature Pressure (bar) 5.24527e+01 -3.46359e+052.90932e-01 -6.00103e+03 Maybe its necessary to say that, I'm coupling my system with both thermostat and barostat, using a reasonable temperature of 300K. A complete .mdp file would be significantly more useful. Did you do energy minimization? If so, what was the outcome? http://www.gromacs.org/Documentation/Terminology/Blowing_Up -Justin Thanks in advance for any instruction, Amin -- Justin A. Lemkul Ph.D. Candidate ICTAS Doctoral Scholar MILES-IGERT Trainee 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 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] water molecules (interfacial)
Dear gromacs users I was trying to count waters involved in the interface between Protein_A and Protein_B. I made a selection.dat file as follows: waterO = group SOL and name OW; heavy1 = group Protein_A and group Protein-H; heavy2 = group Protein_B and group Protein-H; inter = waterO and within 0.35 of heavy1 and within 0.35 of heavy2; is my selection.dat file true? I used g_select -f .xtc -s .tpr -n .ndx -os -oc -oi -om -on -sf selection.dat but gromacs: Nothing selected, finishing up. when I used g_select -f .xtc -s .tpr -n .ndx -os -oc -oi -om -on -sf selection.dat -select, what should be put in front of -select to avoid Nothing selected, finishing up? -- 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] Spurious results with pdb2gmx -chargegrp yes option
Dear All, I would like to get little feedback of my findings/experience with recent simulations I did with the CHARMM force field in GROMACS and the -chargegroup options with pdb2gmx. I have performed two simulations with a helical peptide (25 AA) in a cubic box filled with explicit TIP3P water. I used the GMX 4.5.3. The first simulation was performed with -chargegroup option yes and the second with option no (default). For the first simulation, I used the old version of CHARMM - gromacs forcefield conversion given with GMX 4.0.5. I have done this MD because i would like to compare some suspect results obtained with the same peptide performed 6 months ago with the GMX 4.0.5 version and the old CHARMM-gromacs files available in this distribution. In case of the second simulation, I used the most recent charmm27.ff given in GMX 4.5.3 with a single atom charge groups. The two simulations were performed during 24 ns with the same protocol (i.e. with same mdp files, below and the same starting files). Briefly the protocol i used to equilibrate my system : EM (1000 steps with steep) - NVT at 298 K (with berendsen thermostat and restraints) (100 ps) - (with nose hoover thermostat and Parinello-Rahman barostat) 400 ps. Below md.mdp for the two runs ( i am aware that the Coulomb and electrostatic parameters in the mdp are not otimals for the CHARMM ff, but for sake of comparison with others simulations, i didn't change them). title = KTM17-TIP3 MD constraints = all-bonds integrator = md nsteps = 1200 ; 24000ps ou 24ns dt = 0.002 nstlist = 10 nstcalcenergy = 10 nstcomm = 10 continuation= no ; Restarting after NPT vdw-type= cut-off rvdw= 1.0 rlist = 0.9 coulombtype = PME rcoulomb = 0.9 fourierspacing = 0.12 fourier_nx = 0 fourier_ny = 0 fourier_nz = 0 pme_order= 4 ewald_rtol = 1e-05 optimize_fft= yes nstvout = 5 nstxout = 5 nstenergy = 2 nstlog = 5000 ; update log file every 2 ps nstxtcout = 1000 ; frequency to write coordinates to xtc trajectory every 2 ps Tcoupl = nose-hoover: tc-grps = Protein Non-Protein tau-t = 0.4 0.4 ref-t = 298 298 ; Pressure coupling is on Pcoupl = Parrinello-Rahman pcoupltype = isotropic tau_p = 3.0 compressibility = 4.5e-5 ref_p = 1.0135 gen_vel = no I found that the use of the two -chargegrp options influence drastically the MD results especially in my case for the structural properties of the peptides. Indeed, for the first simulation, the peptide keep in its initial helical with a RMSD around 2.0 A along all 24 ns, whereas in the second simulation, the peptide adopt quickly (~10 ns) a unfold conformation with a RMSD value around 5-6.0 A after only 10 ns. Circular dichroism experiments have shown that this peptide adopt an unfold structure in water. So I am sure that the first simulation is not correct (even if the mdp parameters are not optimal) and that my results confirm that chargegrp yes option should not be used in MD with CHARMM force field in GROMACS, as it was discussed recently in this mailing list (for example http://lists.gromacs.org/pipermail/gmx-users/2010-September/054106.html). SA -- 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] Converting between cartesian coordinates and torsion angles
Hi Mark, Thanks (again, again) for your input! On Fri, Nov 12, 2010 at 5:35 PM, Mark Abraham mark.abra...@anu.edu.auwrote: On 13/11/2010 12:49 AM, Martin Kamp Jensen wrote: Hi, As far as I understand, a topology (a .top file) and a conformation (e.g., a .gro file) contain enough information to calculate the torsion angles of that specific conformation. Table 5.5 (page 124) in the GROMACS manual[1] describes possible interactions (which are contained in the topology) between different molecules while the conformation contains the cartesian coordinates. I did not immediately find a way to convert between the cartesian coordinates and the torsion angles. Can GROMACS do it or do I need to understand (or just find) all the functions/formulas that are referenced in Table 5.5? Section 4.2 has the relevant definitions. Table 5.5 pertains to the definition of force field elements that act upon (for example) such internal coordinates, which is not what you're looking for. I have included screenshots of Table 5.5[2] and the relevant part of some example .top file[3]. (Also, it seems that I can use the read_tpx method defined in include/tpxio.h to read in a topology from a .tpr file. This would then, after converting the cartesian coordinates of some conformation, enable me to work with the torsion angles in my own program before writing cartesian coordinates back for use with GROMACS.) Either a) write something that post-processes the result of grompp -pp in concert with the same coordinate file to get the internal coordinates, or Okay, I see that grompp -pp generates a .top file with a lot of parameters (and maybe even some angles), but for some reason atom numbers have been exchanged with atom names, but of course I can change that back. Then I would need to apply the math from Section 4.2 together with a specific conformation to get the angles (and then do the math to get back to cartesian coordinates after having played with the angles... hmm). Unless my contacts tell me that only a few of those interactions will be needed for our purposes, this seems to be unnecessary extra work. b) use a hacked version of mdrun that writes internal coordinates from within src/gmxlib/bondfree.c (probably used as mdrun -rerun) This could be fine since I will only need to go from cartesian coordinates to angles once. However, I will need to go from angles to cartesian coordinates many, many times. And I would need to write inverse methods since the methods in bondfree.c calculate angles from cartesian coordinates to use them for force and energy calculations. Or am I missing something - is it possible to let GROMACS work on angles instead of cartesian coordinates? (I mean give the input not as, e.g., a .gro file with cartesian coordinates, but as another file with torsion angles - which I will first get via bondfree.c.) to create input for your procedure. Mark Regards, Martin. [1] http://www.gromacs.org/@api/deki/files/126/=gromacs_manual-4.5.pdf [2] http://imada.sdu.dk/~mkjens04/gromacs/intra-molecular_interactions_definitions.pnghttp://imada.sdu.dk/%7Emkjens04/gromacs/intra-molecular_interactions_definitions.png [3] http://imada.sdu.dk/~mkjens04/gromacs/part_of_top_file.pnghttp://imada.sdu.dk/%7Emkjens04/gromacs/part_of_top_file.png -- 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 Regards, Martin. -- 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] Fwd: gmx-users Digest, Vol 79, Issue 102
We have a visualization tool that allows us to visualize the changes in the energy during energy minimization. The areas with more intense color are those where the atoms contribute the most to the total energy value. I just wrote a plug-in to use Gromacs for energy minimizations (or MD simulations) and I need the per-atom energies to use our tool's energy visualization feature. Thanks again, Silvia Message: 3 Date: Sun, 14 Nov 2010 06:31:56 +1100 From: Mark Abraham mark.abra...@anu.edu.au Subject: Re: [gmx-users] how to calculate per-atom energies To: Discussion list for GROMACS users gmx-users@gromacs.org Message-ID: 4cdee7ac.9020...@anu.edu.au Content-Type: text/plain; charset=ISO-8859-1; format=flowed On 14/11/2010 5:58 AM, Silvia Crivelli wrote: Hello, When I run energy minimization, I'd like to obtain the energy per atom in addition to the (potential) energy value for the entire protein. Is there a way to do this? You can get the nonbonded energy per energy group, which you could define to be a single atom, but you are limited to 256 of such groups, and to not using PME. Otherwise, there are approaches involving a lot of hacking about with .top or .tpr files and using mdrun -rerun that could do this for both nonbonded and bonded. However, a more important issue before doing such work is to be sure about what you expect such a decomposition to tell you. A high or low energy for a given atom doesn't necessarily mean anything, and even if it did, the force field wasn't necessarily parametrized to produce reliable per-atom energies (though it probably does do so). Mark -- -- 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] Converting between cartesian coordinates and torsion angles
On 15/11/2010 4:13 AM, Martin Kamp Jensen wrote: Hi Mark, Thanks (again, again) for your input! On Fri, Nov 12, 2010 at 5:35 PM, Mark Abraham mark.abra...@anu.edu.au mailto:mark.abra...@anu.edu.au wrote: On 13/11/2010 12:49 AM, Martin Kamp Jensen wrote: Hi, As far as I understand, a topology (a .top file) and a conformation (e.g., a .gro file) contain enough information to calculate the torsion angles of that specific conformation. Table 5.5 (page 124) in the GROMACS manual[1] describes possible interactions (which are contained in the topology) between different molecules while the conformation contains the cartesian coordinates. I did not immediately find a way to convert between the cartesian coordinates and the torsion angles. Can GROMACS do it or do I need to understand (or just find) all the functions/formulas that are referenced in Table 5.5? Section 4.2 has the relevant definitions. Table 5.5 pertains to the definition of force field elements that act upon (for example) such internal coordinates, which is not what you're looking for. I have included screenshots of Table 5.5[2] and the relevant part of some example .top file[3]. (Also, it seems that I can use the read_tpx method defined in include/tpxio.h to read in a topology from a .tpr file. This would then, after converting the cartesian coordinates of some conformation, enable me to work with the torsion angles in my own program before writing cartesian coordinates back for use with GROMACS.) Either a) write something that post-processes the result of grompp -pp in concert with the same coordinate file to get the internal coordinates, or Okay, I see that grompp -pp generates a .top file with a lot of parameters (and maybe even some angles), but for some reason atom numbers have been exchanged with atom names, but of course I can change that back. Then I would need to apply the math from Section 4.2 together with a specific conformation to get the angles (and then do the math to get back to cartesian coordinates after having played with the angles... hmm). Unless my contacts tell me that only a few of those interactions will be needed for our purposes, this seems to be unnecessary extra work. b) use a hacked version of mdrun that writes internal coordinates from within src/gmxlib/bondfree.c (probably used as mdrun -rerun) This could be fine since I will only need to go from cartesian coordinates to angles once. However, I will need to go from angles to cartesian coordinates many, many times. And I would need to write inverse methods since the methods in bondfree.c calculate angles from cartesian coordinates to use them for force and energy calculations. It seems to me that a better way of thinking/describing about what you want to do is to convert something to an internal coordinate representation, then do some operations on that, and then rebuild the Cartesian coordinates for GROMACS. All of that is best done by pretty much anything you can think of, rather than GROMACS (not least because you may have to rebuild the solvent and whatever else goes along with the changed internal coordinates). There is a body of literature and software that deals with this problem, which is of interest to many areas of computational chemistry (e.g. http://onlinelibrary.wiley.com/doi/10.1002/jcc.20237/abstract, and many Google hits). Or am I missing something - is it possible to let GROMACS work on angles instead of cartesian coordinates? (I mean give the input not as, e.g., a .gro file with cartesian coordinates, but as another file with torsion angles - which I will first get via bondfree.c.) No, GROMACS will only accept Cartesian input. While it is well known that geometry optimization is best done in (redundant) internal coordinates when the evaluation of energy and/or force is expensive (e.g. quantum chemistry on small systems), this is not true for the large majority of EM problems on biomolecular systems. There, the problem is dominated by the presence of multiple minima, the result is not of great interest if you're preparing a system for MD, and there is no payoff for the development time. Thus, GROMACS only computes an internal coordinate immediately before using it to evaluate its contribution to bonded energy and/or force, and then discards it. Mark -- 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
