[gmx-users] H-bonding autocorrelation function

[fyyan at andrew.cmu.edu]

Hi,

I realized I made a mistake, if I take the simplest case, 1-1-1-1-1-1,
c(t) = <h(0)h(t)> / <h>,
then at t = 0.0,     c(0) = 6 / 6 = 1;
            t = 0.006, c(0.006) = 5 / 5 = 1;
            t = 0.012, c(0.012) = 4 /4 = 1;
            t = 0.018, c(0.018) = 3 / 3 = 1;
            t = 0.024, c(0.024) = 2 / 2 = 1;
            t = 0.003, c(0.003) = 1 / 1 = 1;
            t = 0.036, c(0.036) ~;
now go back to 1-1-1-1-1-1-0-1,
it should be (intermittent bonds, bonds allows break and reform)
            t = 0.0,     c(0) = 7 / 7 = 1;
            t = 0.006, c(0.006) = 5. / 6.
            t = 0.012, c(0.012) =  5. 6.;
            t = 0.018, c(0.018) = 4. / 5.;
            t = 0.024, c(0.024) = 3. /4.;
            t = 0.03,   c(0.03) = 2. /3.;
            t = 0.036, c(0.036) = 1. / 2.;
            t = 0.042, c(0.042) = 1. / 1.;

Since gromacs only give the half of the time frames, so the c(t) would be

t = 0.0,   c(0) = 1;
t = 0.006, c(0.006) = 5./6.;
t = 0.012, c(0.012) = 5./6.;
t= 0.018, c(0.018) = 4./5.;

Is my calculation correct? Could anyone give me some feedback? Thanks!

Fangyong


> Hi,
>
> I have a question about the H-bonding autocorrelation function, if I have
> a trajectory, say
> 1-1-1-1-1-1-0-1 (1 is H-bonding, 0 is no H-bonding), the time step is
> 0.006 ps, the correlation is equal to <h(0)h(t)> / <h(0)*h(0)>, my answer
> is
>
> at t = 0.006, correlation = 7/7 = 1;
> at t = 0.012, correlation = 5/7 = 0.714;
> at t = 0.018, correlation = 5/7 = 0.714;
> at t = 0.024, correlation = 4/7 = 0.571;
> at t = 0.03, correlation = 3/7 = 0.429;
> at t = 0.036, correlation = 2/7. = 0.286,
> at t = 0.042, correlation = 1./7. = 0.143,
> at t = 0.048, correlation = 1./7. = 0.143,
> at t = 0.054, correlation = 0.
>
> However, when i check gromacs, I found it gives a different answer, here
> is the command I use "g_hbond -s md.tpr -f new.gro -ac -dist -temp 400 -b
> 72.444 -e 72.486 -dt 0.006", it uses only the time between 72.444 and
> 72.486 ps with a timestep of 0.006 ps, and the following is the h-bond
> number file obtained from gromacs,
>     time (ps)       H-bond ( H-bond without angle requirement)
>     72.444           1           1
>      72.45            1           1
>     72.456           1           1
>     72.462           1           1
>     72.468           1           1
>     72.474           1           1
>      72.48            0           2
>     72.486           1           1
>
> and the hbac gives ,
> time (ps)            C(t)
> 0                         1
> 0.0059967       0.972222
> 0.012001          0.933333
>
> the following is from hbac.xvg,
>
> @ s0 legend "Ac\sfin sys\v{}\z{}(t)"
> @ s1 legend "Ac(t)"
> @ s2 legend "Cc\scontact,hb\v{}\z{}(t)"
> @ s3 legend "-dAc\sfs\v{}\z{}/dt"
>          0           1           1           0     235.274
>  0.0059967    0.411764    0.972222  -9.93411e-09     117.637
>   0.012001   -0.411764    0.933333    0.166667          -0
>
> Am my calculation wrong? Can anyone help me to resolve this? Thanks!
>
>
> Thanks!
>
> Fangyong
>
>
> --
> Gromacs Users mailing list
>
> * Please search the archive at
> http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before
> posting!
>
> * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
>
> * For (un)subscribe requests visit
> https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send
> a mail to gmx-users-request at gromacs.org.
>