[gmx-users] analysis of md run results
Justin A. Lemkul
jalemkul at vt.edu
Wed May 26 21:18:19 CEST 2010
Moeed wrote:
> Hello Justin,
>
> The energy values I get differ mainly for LJ and coulomb(SR) terms. I
> read through chapter 4 of manual but I did not find the difference
> between LJ/coulomb 1-4 and SR.
>
> 1-Where can I read about the differences of these?
For starters:
http://www.gromacs.org/Documentation/Gromacs_Utilities/g_energy
Also realize that 1-4 interactions are intramolecular, SR are intermolecular.
The manual does explain these interactions.
For additional reading, any basic MD textbook will undoubtedly discuss potential
energy, its components, and how they're calculated. The Gromacs manual should
provide some references as well.
> 2-Do you think the difference in energies as a result of excluding
> intramolecualr nonbonded interactions are reasonable?
I have no basis to assess that.
> 3-I can not figure out why pressure in the first case in negative and it
> is very large when exclusions are included.Also in the second case
> kinetic energy is zero. :(
For an NVT ensemble, pressure is basically meaningless. Even in the presence of
pressure coupling, it will fluctuate quite a bit.
http://www.gromacs.org/Documentation/Terminology/Pressure
Since pressure is related to the virial, and the virial is affected by
exclusions, you will see a difference.
http://www.gromacs.org/Documentation/Terminology/Virial
The temperature will be zero in the rerun energy file because the .xtc file does
not contain velocities. If you rerun from a .trr file that has velocities, you
should get a temperature term.
> 4-If I want to perform NPT simulation, the only thing I need to do is to
> switch on Pcoupl option in mdp file?
>
Yes.
> Thank you for your help.
>
> No exclusions: mdrun -s Hexane-Stack125_md.tpr -o
> Hexane-Stack125_md.tpr* _-x trajectoy.xtc_ *-c Hexane-Stack125_after_md
> -v >& output.mdrun_md
>
> g_energy -f ener.edr -o energy.xvg:
>
> tatistics over 5001 steps [ 0.0000 thru 10.0000 ps ], 13 data sets
> The term 'Cons. rmsd ()' is averaged over 501 frames
> All other averages are exact over 5001 steps
>
10 ps may or may not be sufficient to equilibrate your system. By the looks of
it, you haven't even reached your desired temperature yet.
> Energy Average RMSD Fluct. Drift
> Tot-Drift
> -------------------------------------------------------------------------------
> Angle 4400.14 259.533 194.602 59.4737
> 594.856
> Ryckaert-Bell. 1000.11 115.318 88.7683 25.494
> 254.991
> LJ-14 648.7 37.6872 36.5595 3.16904
> 31.6967
> Coulomb-14 -281.387 32.3948 12.468 10.3554
> 103.574
> *LJ (SR) -3322.54 46.1011 40.9609 7.32669
> 73.2816
> Coulomb (SR) 667.867 26.5968 9.4795 -8.60663
> -86.0835*
> Coul. recip. 921.791 26.6407 8.48573 -8.74618
> -87.4793
> *Potential 4034.68 386.493 290.055 88.466
> 884.837*
> Kinetic En. 6336.76 250.507 237.567 27.5244
> 275.299
> Total Energy 10371.4 561.13 450.23 115.99
> 1160.14
> Temperature 297.592 11.7645 11.1568 1.29262
> 12.9288
> P*ressure (bar) -37.0417 1165.58 1161.54 -33.5839
> -335.906*
> Cons. rmsd () 3.74741e-06 2.92095e-07 2.34721e-07 6.02133e-08
> 6.02253e-07
> Heat Capacity Cv: 12.5011 J/mol K (factor = 0.00156281)
>
>
> run with exclusions in top file: mdrun -*rerun trajectoy.xtc* -s
> Hexane-Stack125_md.tpr -o Hexane-Stack125_md.tpr -c
> Hexane-Stack125_after_md -v >& output.mdrun_md
>
> Statistics over 5001 steps [ 0.0000 thru 10.0000 ps ], 13 data sets
> All averages are over 501 frames
>
> Energy Average RMSD Fluct. Drift
> Tot-Drift
> -------------------------------------------------------------------------------
> Angle 4400.09 284.259 237.268 54.2184
> 542.293
> Ryckaert-Bell. 1009.62 124.468 94.0379 28.2417
> 282.474
> LJ-14 648.411 40.6844 39.926 2.70793
> 27.0848
> Coulomb-14 -278.654 34.4227 11.8071 11.1987
> 112.01
> *LJ (SR) -3010.16 47.968 42.5388
> 7.67736 76.789
> Coulomb (SR) 74.7846 2.5186 2.47007 0.170417
> 1.70451*
> Coul. recip. -78.0436 2.58606 1.54725 -0.717662
> -7.17805
> Potential 2766.05 454.268 342.141 103.497
> 1035.18
> Kinetic En. 0 0 0
> 0 0
> Total Energy 2766.05 454.268 342.141 103.497
> 1035.18
> *Temperature 0 0 0
> 0 0
> Pressure (bar) 1973.01 177.068 172.803 13.3792
> 133.819*
> Cons. rmsd () 0 0 0
> 0 0
> Heat Capacity Cv: nan J/mol K (factor = nan)
>
> with: mdp file
>
> Electrostatics/VdW
> coulombtype = PME
> vdw-type = cut-off
> ; Cut-offs
> rlist = 1.0
> rcoulomb = 1.0
> rvdw = 1.0
>
> ; Temperature coupling Berendsen temperature coupling is on in
> two groups
> Tcoupl = berendsen
> tc-grps = HEX ;sol
> tau_t = 0.1 ;0.1
> ref_t = 300 ;300
>
> ; Pressure coupling: Pressure coupling is not on
> Pcoupl = no
> tau_p = 0.5
> compressibility = 4.5e-5
> ref_p = 1.0
>
> ; Velocity generation Generate velocites is on at 300 K.
> Manual p155
> gen_vel = yes
> gen_temp = 300.0
> gen_seed = 173529
>
> ; Bonds
> constraints = all-bonds
> constraint-algorithm = lincs
>
> pbc=xyz
>
Is this the entirety of your .mdp file? If so, you are leaving a large number
of parameters at their default values, rolling the dice and hoping they are
appropriate. If this isn't your whole .mdp file, please *always* post the
entire thing, or else I start rambling about unimportant things :)
-Justin
--
========================================
Justin A. Lemkul
Ph.D. Candidate
ICTAS Doctoral Scholar
MILES-IGERT Trainee
Department of Biochemistry
Virginia Tech
Blacksburg, VA
jalemkul[at]vt.edu | (540) 231-9080
http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin
========================================
More information about the gromacs.org_gmx-users
mailing list