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 ]
     N    opls_238   -0.500     1
     H    opls_241    0.300     1
    CA    opls_224B   0.140     1
    HA    opls_140    0.060     1
    CB    opls_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_type    mass    charge   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 list    gmx-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

Reply via email to