[gmx-users] Heavy water H-H radial distribution function

[Justin Lemkul]

On 6/18/18 3:47 AM, Haelee Hyun wrote:
>
> Dear GROMACS users,
>
> I'm wondering how can I correctly describe H-H radial distribution 
> function of heavy water.
>
> Please check atteched file HH_rdf.PNG which is a calculated result of 
> H-H radial distribution from my simulation.
>

The list does not accept attachments.

> The first peak is due to the intramolecular interaction of water 
> molecules.
>
> It shows almolst 7 at 0.13 nm but when comparing this result with an 
> experimental data, the experimental data shows just almost 2 at the 
> first peak.
>
> I have tried many times of simulations but i couldn't find why this 
> huge difference is caused.
>
> I used tip4p/2005f water model and used potential is atteched below.
>

You need to use the -excl flag to read topology exclusions. Doing so 
will get rid of this intramolecular peak.

-Justin

> [ defaults ]
> ; nbfunc comb-rule gen-pairs fudgeLJ fudgeQQ
> 1  3  yes  0.5 0.5
>
> [ moleculetype ]
> ; molname nrexcl
> SOL  2
>
> [ atoms ]
> ; id at type res nr  residu name at name cg nr charge
> 1       opls_113        1       SOL      OW     1       0.0
> 2       opls_114        1       SOL     HW1     1       0.5564
> 3       opls_114        1       SOL     HW2     1       0.5564
> 4       opls_115        1       SOL      MW     1      -1.1128
>
> ;[nonbond_params]
> ; i j funct q           V            W
> ;1 2 1 0.5564      3.16440e-01  7.74907e-01
> ;1 3 1       0.5564      3.16440e-01  7.74907e-01
>
>
> #ifndef FLEXIBLE
> [ settles ]
> ; OW    funct   doh        dhh
> 1       1       0.09664    0.15555
>
> #else
>
> [ bonds ]
> ; i j funct length   D       beta
> 1 2      3    0.09419   432.581   22.87   ; For TIP4P/2005f Water b0, 
> D, beta
> 1 3      3    0.09419   432.581   22.87   ; For TIP4P/2005f Water b0, 
> D, beta
>
> [ angles ]
> ; i j k funct angle force.c.
> 2 1 3 1 107.4 367.81
> #endif
>
> [ exclusions ]
> 1 2 3 4
> 2 1 3 4
> 3 1 2 4
> 4 1 2 3
>
> ; The position of the virtual site is computed as follows:
> ;
> ; const = distance (OD) / [ cos (angle(DOH))  * distance (OH) ]
> ;   0.015 nm / [ cos (52.26 deg) * 0.09572 nm ]
>
> ; Vsite pos x4 = x1 + a*(x2-x1) + b*(x3-x1)
>
> [ virtual_sites3 ]
> ; Vsite from   funct a  b
> 4 1 2 3 1 0.13288      0.13288
>
>
> I used -DFLEXIBLE option and energy minimization, NVT, NPT 
> equlibration and NVE production run.
>
> If someone would find some wrong things, please let me know.
>
> Thank you.
>
> Haelee Hyun
>
>
>
>
>
>

-- 
==================================================

Justin A. Lemkul, Ph.D.
Assistant Professor
Virginia Tech Department of Biochemistry

303 Engel Hall
340 West Campus Dr.
Blacksburg, VA 24061

jalemkul at vt.edu | (540) 231-3129
http://www.thelemkullab.com

==================================================