Re: [gmx-users] Charges for Coulomb potential
Alright, I understand now. Many thanks for this! Antoine Le 29/08/2012 18:57, Mark Abraham a écrit : On 30/08/2012 3:18 AM, Delmotte, Antoine wrote: Many thanks for your quick response, Mark. So epsilon_r is actually different from the actual relative permittivity (which would be like 80 for water), is that right? You can have an atomistic solvent or a continuum solvent, but not both at once. If so, how does one usually chooses this value of epsilon_r? Is the default value of 1 applicable for most cases (like for simulations in water), or should one obtain it from the literature somehow? The force field defines it (normally to 1). Mark Thank you, Antoine On 29/08/12 14:47, Mark Abraham wrote: On 08/29/2012 11:08 PM, Delmotte, Antoine wrote: Dear Gromacs users, I would like to know which charges are used by Gromacs in the calculation of electrostatic interactions in the standard coulomb potential: E_electrostatics = 138.935 * q1 * q2 / (epsilon_r * r_12); Are q1 and q2 the values tabulated in the file aminoacids.rtp of the force field used? If you use pdb2gmx to generate your .top, yes. More generally, GROMACS uses the charges found in the .top. Could you also please redirect me towards the relevant literature for the calculation of epsilon_r in Gromacs? (or does Gromacs uses a fixed epsilon_r?) You choose it. See parts of manual 4.1 and 7.3 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] Charges for Coulomb potential
Dear Gromacs users, I would like to know which charges are used by Gromacs in the calculation of electrostatic interactions in the standard coulomb potential: E_electrostatics = 138.935 * q1 * q2 / (epsilon_r * r_12); Are q1 and q2 the values tabulated in the file aminoacids.rtp of the force field used? Could you also please redirect me towards the relevant literature for the calculation of epsilon_r in Gromacs? (or does Gromacs uses a fixed epsilon_r?) Thank you very much in advance, Antoine -- 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] Charges for Coulomb potential
Many thanks for your quick response, Mark. So epsilon_r is actually different from the actual relative permittivity (which would be like 80 for water), is that right? If so, how does one usually chooses this value of epsilon_r? Is the default value of 1 applicable for most cases (like for simulations in water), or should one obtain it from the literature somehow? Thank you, Antoine On 29/08/12 14:47, Mark Abraham wrote: On 08/29/2012 11:08 PM, Delmotte, Antoine wrote: Dear Gromacs users, I would like to know which charges are used by Gromacs in the calculation of electrostatic interactions in the standard coulomb potential: E_electrostatics = 138.935 * q1 * q2 / (epsilon_r * r_12); Are q1 and q2 the values tabulated in the file aminoacids.rtp of the force field used? If you use pdb2gmx to generate your .top, yes. More generally, GROMACS uses the charges found in the .top. Could you also please redirect me towards the relevant literature for the calculation of epsilon_r in Gromacs? (or does Gromacs uses a fixed epsilon_r?) You choose it. See parts of manual 4.1 and 7.3 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] 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] Pi stacking potential and parameters in Gromacs?
Dear Gromacs users, I am currently trying to get an expression of the energy of pi stacking interactions (in DNA) and I would like to know how it is calculated in Gromacs. I guess that Gromacs is estimating the energy of the pi stacking as the sum of the van Der Waals and electrostatics contributions using quardrupoles. Is that correct or am I missing something? Secondly, would anyone know where I could find the exact expression and the values parameters used? I imagine that I would need in particular: - The parameters of the van Der Waals term - The charges associated with each atom - The position of the charges of the quadrupole relative to the nucleus of the atom I have no preference in the force field used really, as long as it can give me a reasonable value of the pi stacking interaction. The calculations I am doing are for DNA, and I would like to calculate the energy of the pi stacking interactions from the pdb file directly, and without having to run Gromacs... Any help on this issue will be greatly appreciated! 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
Re: [gmx-users] Re: Empirical potential for Pi stacking interactions?
Many thanks for the references, Andreea. This is very helpful. This is the type of things I was looking for, although if I could find something even simpler than Hunter's potential (which still seem to require an estimation of the charge distribution as I understand it), that would be even better. I will check the citing articles of Hunter and Sanders' paper, and put the reference if I find something in case someone would be looking for the same thing in the future. Thanks a lot. Antoine Le 11/08/2012 08:45, Andreea a écrit : There is a short subchapter in Leach's Molecular modeling (chapter 4), Using charge schemes to study aromatic-aromatic interactions. Maybe you get your hands on the book. It refers to Hunter and Sanders 1990 article: http://pubs.acs.org/doi/abs/10.1021/ja00170a016 There are more than 1000 articles citing the above one. Skim through and maybe your formula is in one of them. Andreea -- View this message in context: http://gromacs.5086.n6.nabble.com/Empirical-potential-for-Pi-stacking-interactions-tp564p583.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] Empirical potential for Pi stacking interactions?
Dear Gromacs users, I am currently looking for an empirical potential for pi stacking interactions. Something like the Leonard Jones potential, or the Mayo potential for hydrogen bonds. I was wondering if something like that exists for pi stacking interactions? I am not really looking for something that would be super accurate, I am mainly looking for an expression which would give me a good enough approximation of the pi stacking interaction energy, even if it is a relatively crude approximation, in the case of DNA bases in particular. Basically, something I could calculate relatively easily from the coordinates of the bases directly, without having to run my PDB file through a full MD package. I know that this is not a question directly related to Gromacs, but this mailing list seems to be the ideal place to get an answer to this type question. If I am misusing the list in any way, feel free to redirect me to a more appropriate forum or mailing list. Many thanks in advance for your help, 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] Atom types for phosphate group in OPLS-AA force field
Dear Gromacs users, I am here requesting your help with regards to the editing of the OPLS-AA force field, and more specifically, the choice of the atom types. I have added an entry for the molecule CABP to the .rtp and .hdb files. I can run pdb2gmx, genbox and editconf without getting any error. However, when I run grompp, I receive the following errors: ERROR 1 [file 1RBO_Protein_chain_A.itp, line 15054]: No default Bond types ERROR 78 [file 1RBO_Protein_chain_M.itp, line 47276]: No default Angle types ERROR 82 [file 1RBO_Protein_chain_M.itp, line 66639]: No default Ryckaert-Bell. types (there are 160 of those). From what I could see in this mailing list, it seems the error comes from the choice of the atom types for the phosphate groups. I have tried several of them which seemed to make sense, but none of them allowed to remove all the errors. The end of the chain is something like: O | (...)--CH2--O--P--O--H | O (details of the whole molecule here: http://www.rcsb.org/pdb/ligand/ligandsummary.do?hetId=CAPsid=1RBO http://www.rcsb.org/pdb/ligand/ligandsummary.do?hetId=CAPsid=1RBO ) The problem is about the atom type of the -OH group at the end. Here are the entries from the .atp file which I thought were relevant: opls_445 30.97376 ; P in MeOPO3--, MeOPO3H2 opls_446 15.99940 ; O= in MeOPO3--, MeOPO3H2 opls_447 15.99940 ; OMe in MeOPO3--, MeOPO3H2 methyl phosphate opls_448 12.01100 ; C in MeOPO3--, MeOPO3H2 6-31+G* CHELPG opls_449 1.00800 ; H in MeOPO3--, MeOPO3H2 I chose opls_445 for the P, opls_446 for the 3 oxygens, and opls_449 for the hydrogen. Apparently, I got it wrong for the oxygen and the hydrogen of the OH group, at least. Unless the error actually comes from something else... So my question is the following: Have I done something wrong about the atom types? If so, does anyone know which atom type I should choose? Details of my rtp entry: [ CAP ] [ atoms ] P1 opls_445 1.539429 1 O1P opls_446 −0.818751 1 O2P opls_446 −0.899293 1 HO2P opls_449 0.429028 1 O3P opls_446 −0.895147 1 O1 opls_447 −0.772438 2 C1 opls_448 0.022084 3 HC11 opls_449 0.176920 3 HC12 opls_449 0.139342 3 C2 opls_159 0.146383 4 O2 opls_154 −0.883774 4 HO2 opls_155 0.527657 4 C opls_271 0.815849 5 O6 opls_271 −0.766517 5 O7 opls_272 −0.771422 5 C3 opls_158 0.142583 6 HC3 opls_140 0.107158 6 O3 opls_154 −0.817159 6 HO3 opls_155 0.471262 6 C4 opls_158 0.184473 7 HC4 opls_140 0.164192 7 O4 opls_154 −0.822792 7 HO4 opls_155 0.498509 7 C5 opls_448 0.038563 8 HC51 opls_449 0.143027 8 HC52 opls_449 0.135804 8 O5 opls_447 −0.730025 9 P2 opls_445 1.526672 10 O4P opls_446 −0.866625 10 O5P opls_446 −0.810148 10 HO5P opls_449 0.521464 10 O6P opls_446 −0.876311 10 [ bonds ] HO2P O2P O1P P1 O2P P2 O3P P1 P1 O1 O1 C1 C1 HC11 C1 HC12 C1 C2 C2 O2 O2 HO2 C2 C C O6 C O7 C2 C3 C3 O3 C3 HC3 O3 HO3 C3 C4 C4 HC4 C4 O4 O4 HO4 C4 C5 C5 HC51 C5 HC52 C5 O5 O5 P2 P2 O4P P2 O5P P2 O6P O5P HO5P Thank you very much in advance! Antoine -- 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] Atom types for phosphate group in OPLS-AA force field
Dear Justin, Many thanks for your response. I managed to reduce the number of errors by trying other atom types, by looking at what seemed most likely to be right from the ffnonbonded.itp and ffbonded.itp files. Unfortunately, I did not manage to find a combination that both made some sense, and allowed to remove the errors. So I guess these things need to be parametrized, as you said. Could you possibly let me know if there is a relatively easy and straightforward way to find these parameters? (I guess the answer is probably no, but it's always worth asking, I guess...) I also noticed that there is a option -zero in grompp which Sets parameters for bonded interactions without defaults to zero instead of generating an error. What does it mean exactly? Is it just taking the angles and distances inferred from the structure as the ones at equilibrium? Is it safe to use it? I am just a beginner with MD, and I just want to do a quick run to test something. I have no need for super high accuracy, and the run I want to try will be very small (maybe 1 or 2 ns, or something like that), and for a very large system (~80 000 atoms). Would that be unreasonable to use that -zero option (it is tempting, I must say...)? Thanks, Antoine On 26/04/12 12:53, Justin A. Lemkul wrote: On 4/26/12 7:45 AM, Delmotte, Antoine wrote: Dear Gromacs users, I am here requesting your help with regards to the editing of the OPLS-AA force field, and more specifically, the choice of the atom types. I have added an entry for the molecule CABP to the .rtp and .hdb files. I can run pdb2gmx, genbox and editconf without getting any error. However, when I run grompp, I receive the following errors: ERROR 1 [file 1RBO_Protein_chain_A.itp, line 15054]: No default Bond types ERROR 78 [file 1RBO_Protein_chain_M.itp, line 47276]: No default Angle types ERROR 82 [file 1RBO_Protein_chain_M.itp, line 66639]: No default Ryckaert-Bell. types (there are 160 of those). From what I could see in this mailing list, it seems the error comes from the choice of the atom types for the phosphate groups. I have tried several of them which seemed to make sense, but none of them allowed to remove all the errors. The end of the chain is something like: O | (...)--CH2--O--P--O--H | O (details of the whole molecule here: http://www.rcsb.org/pdb/ligand/ligandsummary.do?hetId=CAPsid=1RBO http://www.rcsb.org/pdb/ligand/ligandsummary.do?hetId=CAPsid=1RBO ) The problem is about the atom type of the -OH group at the end. Here are the entries from the .atp file which I thought were relevant: opls_445 30.97376 ; P in MeOPO3--, MeOPO3H2 opls_446 15.99940 ; O= in MeOPO3--, MeOPO3H2 opls_447 15.99940 ; OMe in MeOPO3--, MeOPO3H2 methyl phosphate opls_448 12.01100 ; C in MeOPO3--, MeOPO3H2 6-31+G* CHELPG opls_449 1.00800 ; H in MeOPO3--, MeOPO3H2 I chose opls_445 for the P, opls_446 for the 3 oxygens, and opls_449 for the hydrogen. Apparently, I got it wrong for the oxygen and the hydrogen of the OH group, at least. Unless the error actually comes from something else... So my question is the following: Have I done something wrong about the atom types? If so, does anyone know which atom type I should choose? The problem isn't necessarily the atom types themselves (which seem to be reasonable enough from the information provided), but rather that there are no default bonded parameters for these interactions. If you look into ffnonbonded.itp, the second column is entitled bond_type and thus refers to the names used by grompp when interpreting how atoms are connected. The actual bond, angle, and dihedral parameters are in ffbonded.itp, using these names. You would have to map out what the bond_type names are for each of these and find/derive new parameters for each of the interactions that are missing. -Justin Details of my rtp entry: [ CAP ] [ atoms ] P1 opls_445 1.539429 1 O1P opls_446 −0.818751 1 O2P opls_446 −0.899293 1 HO2P opls_449 0.429028 1 O3P opls_446 −0.895147 1 O1 opls_447 −0.772438 2 C1 opls_448 0.022084 3 HC11 opls_449 0.176920 3 HC12 opls_449 0.139342 3 C2 opls_159 0.146383 4 O2 opls_154 −0.883774 4 HO2 opls_155 0.527657 4 C opls_271 0.815849 5 O6 opls_271 −0.766517 5 O7 opls_272 −0.771422 5 C3 opls_158 0.142583 6 HC3 opls_140 0.107158 6 O3 opls_154 −0.817159 6 HO3 opls_155 0.471262 6 C4 opls_158 0.184473 7 HC4 opls_140 0.164192 7 O4 opls_154 −0.822792 7 HO4 opls_155 0.498509 7 C5 opls_448 0.038563 8 HC51 opls_449 0.143027 8 HC52 opls_449 0.135804 8 O5 opls_447 −0.730025 9 P2 opls_445 1.526672 10 O4P opls_446 −0.866625 10 O5P opls_446 −0.810148 10 HO5P opls_449 0.521464 10 O6P opls_446 −0.876311 10 [ bonds ] HO2P O2P O1P P1 O2P P2 O3P P1 P1 O1 O1 C1 C1 HC11 C1 HC12 C1 C2 C2 O2 O2 HO2 C2 C C O6 C O7 C2 C3 C3 O3 C3 HC3 O3 HO3 C3 C4 C4 HC4 C4 O4 O4 HO4 C4 C5 C5 HC51 C5 HC52 C5 O5 O5 P2 P2 O4P P2 O5P P2 O6P O5P HO5P Thank you very much in advance! Antoine -- gmx-users mailing list
Re: [gmx-users] Atom types for phosphate group in OPLS-AA force field
Thank you very much for your answer. This has been very helpful. I think I will follow your advice and try to find similar molecules to find the parameters. Could I just ask you if you know of a database where I could look for the opls parameters which have calculated for other molecules? Best regards, Antoine On 26/04/12 16:55, Justin A. Lemkul wrote: On 4/26/12 11:37 AM, Delmotte, Antoine wrote: Dear Justin, Many thanks for your response. I managed to reduce the number of errors by trying other atom types, by looking at what seemed most likely to be right from the ffnonbonded.itp and ffbonded.itp files. Unfortunately, I did not manage to find a combination that both made some sense, and allowed to remove the errors. So I guess these things need to be parametrized, as you said. Could you possibly let me know if there is a relatively easy and straightforward way to find these parameters? (I guess the answer is probably no, but it's always worth asking, I guess...) You can make a first effort by assigning parameters based on existing, similar, parameters. They may or may not exist. I also noticed that there is a option -zero in grompp which Sets parameters for bonded interactions without defaults to zero instead of generating an error. What does it mean exactly? Is it just taking the angles and distances inferred from the structure as the ones at equilibrium? Is it safe to use it? What it means is ignore the problem and simply assign the value of zero for equilibrium lengths, force constants, etc. I am just a beginner with MD, and I just want to do a quick run to test something. I have no need for super high accuracy, and the run I want to try will be very small (maybe 1 or 2 ns, or something like that), and for a very large system (~80 000 atoms). Would that be unreasonable to use that -zero option (it is tempting, I must say...)? I say that using -zero is a bad idea. It ignores the problem and whatever results you obtain would be questionable, at best. Parameterization of novel species is an advanced topic. Most people will spend several months developing good parameters for a new residue or molecule. -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
Re: [gmx-users] Error from residues added to rtp file
Dear Gromacs users, I am once again requesting your help for the editing of the opls force field .rtp and .hdb files. I have inserted the parameters for a new residue, N-methyl methionine (MME). In this molecule, the methyl group is attached to the nitrogen atom from the peptide bond, which I think is the origin of my problem. pdb2gmx works fine and generates the topology files without any error. Unfortunately, when I looked at the .gro file generated in PyMol (after conversion to pdb with g_editconf), I realized that 3 hydrogens have been added to the nitrogen atom, so I effectively get R2-N-H4: H1 \ / CH3 H2 - N - CA - (rest of the aminoa acid)... H3 / \ H instead of R2-N-H: CH3 \ N - CA - (rest of the aminoa acid)... H / My guess is that Gromacs automatically adds those 3 hydrogens because MME is a terminal residue, but I don't know how to prevent pdb2gmx from doing so. I have tried different modifications in the .rtp file, including using different atom types, removing the peptide bond (the line N -C ) in the list of bonds, changing the hdb file, but nothing seems to have had any effect. I guess I could also simply remove these hydrogens from the .gro and .itp files after using pdb2gmx, but I would have a better confidence in what I changed in the force field files if pdb2gmx was giving me the right answer directly. See below my .rtp and .hdb entries: [ MME ] [ atoms ] Nopls_238 -0.500 1 CAopls_224B 0.140 1 HAopls_1400.060 1 CBopls_136 -0.120 2 HB1opls_1400.060 2 HB2opls_1400.060 2 CGopls_2100.048 3 HG1opls_1400.060 3 HG2opls_1400.060 3 SDopls_202 -0.335 4 CEopls_209 -0.013 5 HE1opls_1400.060 5 HE2opls_1400.060 5 HE3opls_1400.060 5 Copls_2350.500 6 Oopls_236 -0.500 6 CMopls_244 -0.110 7 HM1opls_1400.037 7 HM2opls_1400.037 7 HM3opls_1400.037 7 [ bonds ] NCA NCM CM HM1 CM HM2 CM HM3 CAHA CACB CA C CB HB1 CB HB2 CBCG CG HG1 CG HG2 CGSD SDCE CE HE1 CE HE2 CE HE3 C O MME 6 1 1 H N CM CA 1 5 HA CA N C CB 2 6 HB CB CG CA 2 6 HG CG SD CB 3 4 HE CE SD CG 3 4 HM CM N CA Once again, any idea about this problem would be greatly appreciated. Thanks to all in advance. Regards, Antoine On 08/26/2011 05:42 PM, Delmotte, Antoine wrote: Oh, thank you so much! That was indeed the error. It's amazing how these little things can sometimes drive you mad Thanks a lot, Antoine On 08/26/2011 05:27 PM, Thomas Piggot wrote: Hi, I think the problem is that you have a dash rather than a minus symbol for the sign of the charge on the OD atom. Cheers Tom Delmotte, Antoine wrote: Dear Gromacs users, I am currently trying to run an MD simulation with the OPLS-AA force field on a protein having different non standard residues and a ligand. I found the charges for the OPLS force field for these residues in the literature and I am now trying to add them in the OPLS force field parameter files. I have edited the aminoacids.rtp and the aminoacids.hdb files for the OPLS-AA force field, as well as the residuetypes.dat file. Here is an example for one of the amino acids, hydroxyproline: [ HYP ] [ atoms ] N opls_239 -0.140 1 CA opls_246 0.010 1 HA opls_140 0.060 1 CB opls_136 -0.120 2 HB1 opls_140 0.060 2 HB2 opls_140 0.060 2 CG opls_137 -0.120 3 HG1 opls_140 0.060 3 OD opls_167 −0.683 3 HD opls_168 0.743 3 CD opls_245 -0.050 4 HD1 opls_140 0.060 4 HD2 opls_140 0.060 4 C opls_235 0.500 5 O opls_236 -0.500 5 [ bonds ] N CA CA HA CA CB CA C CB HB1 CB HB2 CB CG CG HG1 CG OD OD HD CG CD CD HD1 CD HD2 CD N C O -C N [ impropers ] -C CA N CD improper_Z_N_X_Y CA +N C O improper_O_C_X_Y When I run pdb2gmx, I get the following error, which is not very informative: All occupancies are one Opening force field file /usr/share/gromacs/top/oplsaa.ff/atomtypes.atp Atomtype 1 Reading residue database... (oplsaa) Opening force field file /usr/share/gromacs/top/oplsaa.ff/aminoacids.rtp Residue 58 --- Program g_pdb2gmx, VERSION 4.5.3 Source code file: /builddir/build/BUILD/gromacs-4.5.3/src/kernel/resall.c, line: 389 Fatal error: in .rtp file in residue HYP at line: OD opls_167 −0.683 3 I would be grateful if anyone could shed some light on the origin of this error, and on what I can do to correct it. I am using Gromacs 4.5.3. Thanks a lot in advance, Antoine -- gmx-users mailing listgmx-users
Re: [gmx-users] Error from residues added to rtp file
This solved my problem, although I am not totally sure to have perfectly understood how the .tdb had to be edited. I added the following entry: [ MME-NH ] [ delete ] H Anyway, it solved the problem I had. Many thanks for your help! Antoine On 09/05/2011 09:54 AM, Mark Abraham wrote: On 5/09/2011 6:41 PM, Delmotte, Antoine wrote: Dear Gromacs users, I am once again requesting your help for the editing of the opls force field .rtp and .hdb files. I have inserted the parameters for a new residue, N-methyl methionine (MME). In this molecule, the methyl group is attached to the nitrogen atom from the peptide bond, which I think is the origin of my problem. pdb2gmx works fine and generates the topology files without any error. Unfortunately, when I looked at the .gro file generated in PyMol (after conversion to pdb with g_editconf), I realized that 3 hydrogens have been added to the nitrogen atom, so I effectively get R2-N-H4: H1 \ / CH3 H2 - N - CA - (rest of the aminoa acid)... H3 / \ H instead of R2-N-H: CH3 \ N - CA - (rest of the aminoa acid)... H / My guess is that Gromacs automatically adds those 3 hydrogens because MME is a terminal residue, but I don't know how to prevent pdb2gmx from doing so. That uses the .tdb terminal database system. You will need a custom entry to generate your two hydrogen atoms, and to use pdb2gmx -ter to get a chance to select it. Separately, if this residue is ever to be a non-terminal one, you'll need to define one H atom bound to the N atom in the .rtp entry, so do that now. Mark I have tried different modifications in the .rtp file, including using different atom types, removing the peptide bond (the line N -C ) in the list of bonds, changing the hdb file, but nothing seems to have had any effect. I guess I could also simply remove these hydrogens from the .gro and .itp files after using pdb2gmx, but I would have a better confidence in what I changed in the force field files if pdb2gmx was giving me the right answer directly. See below my .rtp and .hdb entries: [ MME ] [ atoms ] Nopls_238 -0.500 1 CAopls_224B 0.140 1 HAopls_1400.060 1 CBopls_136 -0.120 2 HB1opls_1400.060 2 HB2opls_1400.060 2 CGopls_2100.048 3 HG1opls_1400.060 3 HG2opls_1400.060 3 SDopls_202 -0.335 4 CEopls_209 -0.013 5 HE1opls_1400.060 5 HE2opls_1400.060 5 HE3opls_1400.060 5 Copls_2350.500 6 Oopls_236 -0.500 6 CMopls_244 -0.110 7 HM1opls_1400.037 7 HM2opls_1400.037 7 HM3opls_1400.037 7 [ bonds ] NCA NCM CM HM1 CM HM2 CM HM3 CAHA CACB CA C CB HB1 CB HB2 CBCG CG HG1 CG HG2 CGSD SDCE CE HE1 CE HE2 CE HE3 C O MME 6 1 1 H N CM CA 1 5 HA CA N C CB 2 6 HB CB CG CA 2 6 HG CG SD CB 3 4 HE CE SD CG 3 4 HM CM N CA Once again, any idea about this problem would be greatly appreciated. Thanks to all in advance. Regards, Antoine On 08/26/2011 05:42 PM, Delmotte, Antoine wrote: Oh, thank you so much! That was indeed the error. It's amazing how these little things can sometimes drive you mad Thanks a lot, Antoine On 08/26/2011 05:27 PM, Thomas Piggot wrote: Hi, I think the problem is that you have a dash rather than a minus symbol for the sign of the charge on the OD atom. Cheers Tom Delmotte, Antoine wrote: Dear Gromacs users, I am currently trying to run an MD simulation with the OPLS-AA force field on a protein having different non standard residues and a ligand. I found the charges for the OPLS force field for these residues in the literature and I am now trying to add them in the OPLS force field parameter files. I have edited the aminoacids.rtp and the aminoacids.hdb files for the OPLS-AA force field, as well as the residuetypes.dat file. Here is an example for one of the amino acids, hydroxyproline: [ HYP ] [ atoms ] N opls_239 -0.140 1 CA opls_246 0.010 1 HA opls_140 0.060 1 CB opls_136 -0.120 2 HB1 opls_140 0.060 2 HB2 opls_140 0.060 2 CG opls_137 -0.120 3 HG1 opls_140 0.060 3 OD opls_167 −0.683 3 HD opls_168 0.743 3 CD opls_245 -0.050 4 HD1 opls_140 0.060 4 HD2 opls_140 0.060 4 C opls_235 0.500 5 O opls_236 -0.500 5 [ bonds ] N CA CA HA CA CB CA C CB HB1 CB HB2 CB CG CG HG1 CG OD OD HD CG CD CD HD1 CD HD2 CD N C O -C N [ impropers ] -C CA N CD improper_Z_N_X_Y CA +N C O improper_O_C_X_Y When I run pdb2gmx, I get the following error, which is not very informative: All occupancies are one Opening force field file
[gmx-users] Error from residues added to rtp file
Dear Gromacs users, I am currently trying to run an MD simulation with the OPLS-AA force field on a protein having different non standard residues and a ligand. I found the charges for the OPLS force field for these residues in the literature and I am now trying to add them in the OPLS force field parameter files. I have edited the aminoacids.rtp and the aminoacids.hdb files for the OPLS-AA force field, as well as the residuetypes.dat file. Here is an example for one of the amino acids, hydroxyproline: [ HYP ] [ atoms ] N opls_239 -0.140 1 CA opls_246 0.010 1 HA opls_140 0.060 1 CB opls_136 -0.120 2 HB1 opls_140 0.060 2 HB2 opls_140 0.060 2 CG opls_137 -0.120 3 HG1 opls_140 0.060 3 OD opls_167 −0.683 3 HD opls_168 0.743 3 CD opls_245 -0.050 4 HD1 opls_140 0.060 4 HD2 opls_140 0.060 4 C opls_235 0.500 5 O opls_236 -0.500 5 [ bonds ] N CA CA HA CA CB CA C CB HB1 CB HB2 CB CG CG HG1 CG OD OD HD CG CD CD HD1 CD HD2 CD N C O -C N [ impropers ] -C CA N CD improper_Z_N_X_Y CA +N C O improper_O_C_X_Y When I run pdb2gmx, I get the following error, which is not very informative: All occupancies are one Opening force field file /usr/share/gromacs/top/oplsaa.ff/atomtypes.atp Atomtype 1 Reading residue database... (oplsaa) Opening force field file /usr/share/gromacs/top/oplsaa.ff/aminoacids.rtp Residue 58 --- Program g_pdb2gmx, VERSION 4.5.3 Source code file: /builddir/build/BUILD/gromacs-4.5.3/src/kernel/resall.c, line: 389 Fatal error: in .rtp file in residue HYP at line: OD opls_167 −0.683 3 I would be grateful if anyone could shed some light on the origin of this error, and on what I can do to correct it. I am using Gromacs 4.5.3. Thanks a lot in advance, Antoine -- 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] Error from residues added to rtp file
Oh, thank you so much! That was indeed the error. It's amazing how these little things can sometimes drive you mad Thanks a lot, Antoine On 08/26/2011 05:27 PM, Thomas Piggot wrote: Hi, I think the problem is that you have a dash rather than a minus symbol for the sign of the charge on the OD atom. Cheers Tom Delmotte, Antoine wrote: Dear Gromacs users, I am currently trying to run an MD simulation with the OPLS-AA force field on a protein having different non standard residues and a ligand. I found the charges for the OPLS force field for these residues in the literature and I am now trying to add them in the OPLS force field parameter files. I have edited the aminoacids.rtp and the aminoacids.hdb files for the OPLS-AA force field, as well as the residuetypes.dat file. Here is an example for one of the amino acids, hydroxyproline: [ HYP ] [ atoms ] N opls_239 -0.140 1 CA opls_246 0.010 1 HA opls_140 0.060 1 CB opls_136 -0.120 2 HB1 opls_140 0.060 2 HB2 opls_140 0.060 2 CG opls_137 -0.120 3 HG1 opls_140 0.060 3 OD opls_167 −0.683 3 HD opls_168 0.743 3 CD opls_245 -0.050 4 HD1 opls_140 0.060 4 HD2 opls_140 0.060 4 C opls_235 0.500 5 O opls_236 -0.500 5 [ bonds ] N CA CA HA CA CB CA C CB HB1 CB HB2 CB CG CG HG1 CG OD OD HD CG CD CD HD1 CD HD2 CD N C O -C N [ impropers ] -C CA N CD improper_Z_N_X_Y CA +N C O improper_O_C_X_Y When I run pdb2gmx, I get the following error, which is not very informative: All occupancies are one Opening force field file /usr/share/gromacs/top/oplsaa.ff/atomtypes.atp Atomtype 1 Reading residue database... (oplsaa) Opening force field file /usr/share/gromacs/top/oplsaa.ff/aminoacids.rtp Residue 58 --- Program g_pdb2gmx, VERSION 4.5.3 Source code file: /builddir/build/BUILD/gromacs-4.5.3/src/kernel/resall.c, line: 389 Fatal error: in .rtp file in residue HYP at line: OD opls_167 −0.683 3 I would be grateful if anyone could shed some light on the origin of this error, and on what I can do to correct it. I am using Gromacs 4.5.3. Thanks a lot in advance, Antoine -- 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
