[gmx-developers] What is dH/dl functional for Coulomb interaction?
Berk Hess
hessb at mpip-mainz.mpg.de
Wed Jun 25 10:32:19 CEST 2008
Hi,
All the equations for free-energy Coulomb interaction in the manual are
incorrect.
The correct linear Coulomb interpolation is:
V=1/(4 pi eps0 rij) ((1-lambda) qAi qAj + lamda qBi qBj)
dV/dl = 1/(4 pi eps0 rij) (-qAi qAj + qBi qBj)
So there is no explicit lambda dependence for dV/dl
(unless you use soft-core interactions).
Equation 4.120 gives the correct dV/dl for soft-core interactions
and also without soft-core when you set alpha=0.
We will correct the equations for the 4.0 manual.
Berk.
Leontyev Igor wrote:
> Dear gmx-developers,
> I am confused by the fact that dhdl-values calculated with the same
> molecular structure for different values of the Coulomb interaction
> coupling
> parameter lambda (L) does not depend on the lambda. According to the
> Gromacs
> 3.3 Manual the L-dependence of Coulomb part of Hamiltonian has form
> given by
> Eq.(4.103) such as dH(L)/dL should explicitly depend on lambda.
> However, my
> tests show opposite. I have tried to calculate a charging free energy of
> QM-system (solute) for two different systems:
> 1) QM-system (solute) in protein. The QM-system consist of metal center
> alone with imidozole rings of ligating residues (35 atoms) and the
> protein
> is represented by 11938 atoms. The QM-system has non-zero atomic
> charges in
> both the initial and final states. Total charge difference of the system
> between final and initial state is -1.
> 2) Solute alone (no solvent). The model QM-system is constituted by two
> point charges separated by 1nm . Point charges have values (0,0) in the
> initial state and (-1,1) in the final state.
>
> The independence of dhdl-values on the parameter lambda (L) was
> observed for
> both systems. For the model system (2 point charges), dhdl-values can be
> calculated manually. According to Eq. (4.103), H(L)=f/(1nm)*L*(-1)*L*(1),
> thus:
> H(L)=-138.935*L*L (1)
> and
> dH(L)/dL=-2*138.935*L
> However, results obtained by Gromacs 3.3.2 are not described by the above
> functional (1). The results are:
> L Potential energy dhdl
> 0.00000 0.00000e+00 -1.38935e+02
> 0.50000 -6.94677e+01 -1.38935e+02
> 1.00000 -1.38935e+02 -1.38935e+02
>
> It is interesting to notice that these results are precisely described by
> the conventional functional form for Coulomb coupling:
> H(L)=HA+L*(HB-HA), (2)
>
> where HB and HA are Hamiltonians in the final and initial states,
> respectively. Thus, could someone confirm that the Coulomb coupling
> actually
> implemented in Gromacs is given by functional (2) instead of fuctional
> given
> by Eqs (4.103)-(4.106) of the Manual ?
> Thanks,
> Igor
>
>
> _______________________________________________
> gmx-developers mailing list
> gmx-developers at gromacs.org
> http://www.gromacs.org/mailman/listinfo/gmx-developers
> Please don't post (un)subscribe requests to the list. Use the www
> interface or send it to gmx-developers-request at gromacs.org.
>
> This email was Anti Virus checked by Astaro Security Gateway.
> http://www.astaro.com
>
More information about the gromacs.org_gmx-developers
mailing list