[gmx-users] NPT for compression
Justin A. Lemkul
jalemkul at vt.edu
Wed Jul 14 00:18:22 CEST 2010
Moeed wrote:
> Hello,
>
> I followed the hints. I did several NPT runs and I found out pressure
> above 100 eg 120 bar gives final boxsize close to volume of box I want.
> But I dont know why pressure given by g_energy is totally different.
>
If you're compressing the system, it likely has not had time to reach
equilibrium. You have a negative (average) pressure, which to me indicates your
system has collapsed in on itself, which sounds correct. If you want to
equilibrate at a higher pressure for some reason, run a longer simulation.
Without seeing your actual .mdp file, it's hard to comment further.
> 1- I looked at values in xvg file and I notoced pressre starts froma
> very low value, also takes negative and postive values. If I am doing
> NPT why P is changing ?
The terms "constant" and "fixed" are not equivalent. A barostat will produce a
distribution of pressures such that (after sufficient time) the average should
be correct. For more:
http://www.gromacs.org/Documentation/Terminology/Pressure
Since you have only one molecule with a few hundred atoms, I would expect you'd
see very large fluctuations.
> 2- at the end of xvg file I see final dimensions are 2*2*2 . Also I
What .xvg file is this?
> checked the size in ngmx and I see that long chain is compressing very
> nicely and takes a very convoluted conformation in a very small size. I
> selected 1000 ps but after 727 ps simulation stops (it seems box size
> can not be reduced further).However, When I open the gro file
So the simulation crashes?
> PE60_after_MD dimensions are 22*22*22! :( why is this happening?
>
If the simulation crashed due to excess pressure, probably the box is exploding.
>
> Energy Average RMSD Fluct. Drift
> Tot-Drift
> -------------------------------------------------------------------------------
> Angle 725.016 41.3878 36.4581 0.0790551
> 67.8611
> Ryckaert-Bell. 231.804 29.096 25.866 0.0537685
> 46.155
> LJ-14 161.418 8.09535 7.98175 0.00545371
> 4.68147
> Coulomb-14 -87.8111 8.74814 7.24162 0.0198065
> 17.0019
> LJ (SR) -145.619 82.9315 38.527 -0.296365
> -254.401
> Coulomb (SR) 190.015 14.0939 11.1829 -0.0346165
> -29.7149
> Potential 1074.82 73.9367 60.258 -0.172898
> -148.416
> Kinetic En. 899.733 29.9816 29.9813 0.000554359
> 0.475863
> Total Energy 1974.56 67.628 52.4374 -0.172344
> -147.94
> Temperature 299.757 9.98875 9.98864 0.000184692
> 0.15854
> *Pressure (bar) -1.43379 141.354 141.205 -0.0261305
> -22.4305*
> Cons. rmsd () 5.12607e-06 4.18181e-07 3.98738e-07 5.08597e-10
> 4.3658e-07
> Box-X 10.7454 7.6259 2.33711 -0.0292936
> -25.1457
> Box-Y 10.7454 7.6259 2.33711 -0.0292936
> -25.1457
> Box-Z 10.7454 7.6259 2.33711 -0.0292936
> -25.1457
> Density (SI) 32.6368 57.0897 38.1186 0.171508
> 147.222
> pV -5.00443 428.488 428.332 -0.0466916
> -40.0802
> Vir-XX 321.928 424.592 424.396 0.051996
> 44.6335
> Mu-X 0.00138297 0.33292 0.332911 -9.74152e-06
> -0.00836214
> Lamb-System 1.00002 0.000363012 0.000363011 -3.17825e-09
> -2.72822e-06
> Heat Capacity Cv: 12.4926 J/mol K (factor = 0.00111041)
>
> 3- The systems seems not to be equilibrated yet. eg LJ SR remains
> constant for a while and goes down and its not reaching a plateau after
> 720 ps (where simulation stops).
>
Seems about right, given the above description. Not only is your system not
equilibrated, it sounds like it's completely unstable after a certain point.
> 4- When I choose p 200bar , compression continues till 420 ps. Actually
> chain takes a specific conformation in such a way that it can not get
> compressed further and simulation stops even before former trial (p 120
> bar).
Take a more stepwise approach to all of this. Instead of applying massive
pressure and just crushing your PE chain, apply high pressure (100 bar) for a
short period of time (maybe 100 ps), then reduce the pressure and simulate for
longer (i.e. 50 bar for 500 ps), etc. All that is complete guesswork, I've
never done it, but hopefully you get the idea, sort of like pressure annealing.
> a simple question:
> 5- About tc groups in mdpfile: actually I see people having names like
> protein SOL as their groups but I am wondering how I can set the name
> "polymer". I thought this comes form the [system] or [molecule] in top
> file but its not so (or thought may names are coming form index file).
> If I want to add hexane to PE and calculate energies I need to have
> these two groups separately. for PE I have three residues Eth EthB
> EthE.Could you please tell me how can I include them in one single group
> 'Polymer"? and use this name in mdp file?
>
If you have a single PE molecule in your system, then by default, you do not
need to use any special index groups. You need "tc_grps = System" and nothing
more. For instance:
http://www.gromacs.org/Documentation/Terminology/Default_Index_Groups
I would suggest you read about thermostats and index groups before trying to
make up anything fancy. Better yet, do a basic protein-in-water tutorial and
understand how a proper [molecules] section is constructed and how it relates to
the tc_grps you can set. There is also some very important information here:
http://www.gromacs.org/Documentation/Terminology/Thermostats
-Justin
> Thanks for your help.
>
--
========================================
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