[gmx-users] Problem in calculating static dielectric constant of Lysozymes

[David van der Spoel]

On 31/10/17 20:56, Suman Das wrote:
> Dear Gromacs Users,
> 
> I am trying to calculate static dielectric constant of Lysozymes using
> gromacs version 4.6.7. I tried to follow the protocol used earlier by van
> Gunsteren (J. Phys. Chem. 1993, 97, 2009-2014).  I used PME for the
> electrostatic interaction. First, I converted the trajectory of protein by
> using the command 'trjconv -f traj.xtc -o protein.xtc -pbc whole -e
> 1000.0". Then I removed the translational and rotational motion of the
> protein w.r.t the C-alpha atom by the following command "trjconv -f
> protein.xtc -o protein_fit.xtc -fit rot+trans". I choose C-alpha for the
> least square fit. And finally g_dipoles command to calculate the dielectric
> constant "g_dipoles -f protein_fit.xtc -o Mtot.xvg -eps epsilon.xvg -a
> aver.xvg -d dipdist.xvg -c dipcorr.xvg". However, the dielectric constant
> value (1.04) that I am getting is very less compared to the actual
> dielectric constant (30).

Since the proteins have a net charge you have to translate them to their 
center of charge such that the leading term Q/r in Eqn 1 in the paper 
disappears (in the paper they use center of mass which is not strictly 
the same but a good approximation). The problem is that trjconv does not 
do this and therefore your net dipole moment is incorrect (since it is 
origin dependent). A roundabout way to fix this would be to compute the 
center of mass with g_traj (gmx traj) and multiply that by the total 
charge of the protein and subtract that from the dipole components. The 
compute epsilon yourself.
> 
> 
> Dipole moment (Debye)
> ---------------------
> Average  = 172.7138  Std. Dev. =   7.4113  Error =   0.0524
> 
> The following averages for the complete trajectory have been calculated:
> 
>   Total < M_x > = 38.9789 Debye
>   Total < M_y > = 54.2911 Debye
>   Total < M_z > = -158.896 Debye
> 
>   Total < M_x^2 > = 1575.16 Debye^2
>   Total < M_y^2 > = 3006.62 Debye^2
>   Total < M_z^2 > = 25303.2 Debye^2
> 
>   Total < |M|^2 > = 29885 Debye^2
>   Total |< M >|^2 = 29714.9 Debye^2
> 
>   < |M|^2 > - |< M >|^2 = 170.091 Debye^2
> 
> Finite system Kirkwood g factor G_k = 0.00570201
> Infinite system Kirkwood g factor g_k = 0.00561463
> 
> Epsilon = 1.04819
> 
> Could you help me to find out what probably went wrong with the calculation?
> 
> Thanks,
> Suman
> 


-- 
David van der Spoel, Ph.D., Professor of Biology
Head of Department, Cell & Molecular Biology, Uppsala University.
Box 596, SE-75124 Uppsala, Sweden. Phone: +46184714205.
http://www.icm.uu.se