Hi Oliver, I think the eigenvalues in NMA are not the same (there used to be a factor of 2PI and the mass weighting). Maybe you can try my script from the user contributions and see if you get something more reasonable (use to flag -n to indicate that your eigenvalues are from NMA).
Ran. oliver.k...@uni-duisburg-essen.de wrote: > Dear Gromacs Users, > I'm trying to calculate entropies from a md trajectory using g_anaeig. > There are two ways to go (question at bottom ;-): > > 1. NMA and quasi-harmonic approximation: Use a bunch of snapshots (maybe > 5-20), minimize each of them to very low maximum forces, calculate the > hessian matrix, diagonalize and use g_anaeig to calculate the entropy from > the resulting eigenvector-matrix assuring that there are no negative > eigenvalues in the eigenvectors 7 to N (first six eigenvectors will not be > part of the calculation). - as follows: > > # Energy Minimization > grompp_d -f em_nma.mdp -t md.fitted.trr -time $t -c md.gro -p protein.top > -o $t.em.tpr > mdrun_d -v -deffnm $t.em -table table6-12_4r_doublePrecision.xvg -tablep > table6-12_4r_doublePrecision.xvg > > # Hessian Matrix > grompp_d -f nma.mdp -t $t.em.trr -c md.gro -p protein.top -o $t.hessian.tpr > mdrun_d -v -deffnm $t.hessian -table table6-12_4r_doublePrecision.xvg > -tablep table6-12_4r_doublePrecision.xvg > > # Diagonalization of the Hessian > g_nmeig_d -f $t.hessian.mtx -s $t.hessian.tpr -first 1 -last 10000 -v > $t.eigenvec.trr > > # Entropy calculation - vibrational (without first 6 modes) > g_anaeig_d -v $t.eigenvec.trr -f $t.em.trr -s $t.hessian.tpr -temp 298.15 > -nevskip 6 -entropy 2>&1 | tee out.anaeig.Svib.$t > > grep 'The Entropy due to the Quasi Harmonic approximation is' > out.anaeig.Svib.$t | awk '{print $10}' >> result/Svib.nma > > I use distance-dependent dielectric e=4r, but that doesn't make much > difference. > > 2. Schlitter approximation based on covariance: Use all snapshots of the > md trajectory, calculate the covariance matrix (g_covar), - diagonalized > matrix will be returned -, and subsequently calculate the entropy with > g_anaeig. - as follows: > > # covariance matrix as time average over configurations > g_covar_d -f md$i.trr -s md.gro -v md$i.eigenvec.trr -mwa -av > average.$i.pdb -ascii covar.$i -xpm covar.$i -xpma covara.$i -l covar.$i > -o md$i.eigenval.xvg <<- EOF > 0 > 0 > EOF > > # Analysis of the principal components (and entropy calculation) > g_anaeig -v md$i.eigenvec.trr -f md$i.trr -s md.gro -first 1 -last -1 > -entropy > out.anaeig.schlitter.$i > > grep 'The Entropy due to the Schlitter formula is' out.anaeig.schlitter.$i > | awk '{print $9}' >> result/Svib.schlitter > > Somebody before mentioned, he would like to have the undiagonalized > covariance matrix as input for the entropy calculation, I think, that > doesn't make a difference, am I right? > > So, practically, I tried to reproduce entropy from Schlitter 1993. A > simulation of a deca-alanine-helix in vacuo in the old gmx force-field > with vdw-cut-off etc. and I could reproduce the value of ca. 700 > kJoule/mol K with the Schlitter approximation. > And now the question, why don't I get the same range of values when doing > normal-mode analysis (as described above)? > > values of the Schlitter-approximation (for different simulation lengths): > 667.365 > 685.594 > 681.259 > 680.269 > values given by the quasi-harmonic approximation when calculating from > covariance: > 582.731 > 596.97 > 590.71 > 589.07 > > values from NMA and quasi-harmonic approximation (for 3 snapshots): > 21662.9 > 21674.9 > 21662.9 > > So, there is a factor of round 30 between hessian- and covariance-based > entropy!? > > I'm totally stuck with this. > If anybody has experience with this phenomenon, any help is appreciated. > Thanx in advance > > Oliver Kuhn > > > > > > > > _______________________________________________ > gmx-users mailing list gmx-users@gromacs.org > http://www.gromacs.org/mailman/listinfo/gmx-users > Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php > > -- ------------------------------------------------------ Ran Friedman Postdoctoral Fellow Computational Structural Biology Group (A. Caflisch) Department of Biochemistry University of Zurich Winterthurerstrasse 190 CH-8057 Zurich, Switzerland Tel. +41-44-6355593 Email: r.fried...@bioc.unizh.ch Skype: ran.friedman ------------------------------------------------------ _______________________________________________ gmx-users mailing list gmx-users@gromacs.org http://www.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php