[gmx-users] Re: Lennard-Jones Rdf plot

[Justin A. Lemkul]

Lum Nforbi wrote:
> Justin,
> 
>   This is in response to your mail on this subject. My objective was to 
> do an NPT simulation of a 200 particle fluid using Lennard-Jones 
> parameters only. I also used the atomic number, atomic mass, and 
> Lennard-Jones parameters(sigma in nm and epsilin in KJ/mol) for the 
> oxygen atom, and the force field files which I generated. I have 
> included my procedure and  the parameter and topology files I used below.
> 
> (1) I chose one set of 3 random coordinates(x,y,z) and multiplied them 
> using packmol to generate 200x3 random coordinates.

<snip>

> Pressure (bar)             0.278412     89.842    89.8327 0.000447827    
> 4.47828
> Volume                      10.0719    0.42992   0.426807 -1.78902e-05  
> -0.178902
> Density (SI)                528.134    13.6698    13.5568 0.000607606    
> 6.07606

I would say your system is not at equilibrium.  Your average pressure is far 
from your ref_p of 1.0 bar, and the average volume would indicate average box 
vectors of about 2.15 nm (on average).  So, if you started at 2.7 nm, the box 
has likely reached a value of less than 2.0 nm, indicating that your 1.0-nm 
cutoffs are violating the minimum image convention.  Plot these data 
individually and look for systematic trends that can be giving errors in your 
data.  Probably the additional peaks you describe and bizarre behavior at the 
tail of the plot is due to spurious PBC interactions.

Again, I did not open the attachment.  Nothing personal, but as I said before, 
public file sharing sites are more appropriate.  I just don't like unsolicited 
attachments, and further, people who subscribe to the list in digest form will 
not get them, so they will not necessarily be able to help as easily.

-Justin

> Heat Capacity Cv:      12.5342 J/mol K (factor = 0.00331954)
> Isothermal Compressibility:  0.0016631 /bar
> Adiabatic bulk modulus:        601.288  bar
> 
> (6) I then plotted a Potential Energy graph and a radial distribution 
> function for the system at the end of the simulation. I have attached 
> the graphs. The commands for g_rdf are below:
> 
> make_ndx -f oxymdrun.gro -o oxy.ndx
> trjconv -f oxymdrun_traj.trr -o oxymdrun_traj.xtc
> g_rdf -f oxymdrun_traj.xtc -n oxy.ndx -o oxy_rdf.xvg
> 
> Select a reference group and 1 group
> Group     0 (      System) has   200 elements
> Group     1 (         OXY) has   200 elements
> Select a group: 1
> Selected 1: 'OXY'
> Select a group: 1
> Selected 1: 'OXY'
> Last frame      10000 time 10000.000 
> 
> 
> The first two peaks on the radial distribution function plot are ok, but 
> the third one starts slanting down ie it does not converge to 1.0. I am 
> not sure what I am not setting right and I don't know what could be the 
> cause of the graph slanting downwards.
> In section 2.3 of page 9 of the manual, it is mentioned that reduced 
> units are good for simulating LJ systems, but Gromacs already has a 
> fixed set units that it uses and I am not sure how to set the reduced 
> units such that Gromacs can recognize them. Could it be a reduced units 
> problem?
> 
> I appreciate your answers.
> 
> Lum
> - Show quoted text -
> 
> 
> 
> 
> Lum Nforbi wrote:
>  > Dear all,
>  >
>  >    Please, could someone tell me what could be the problem with the
>  > attached LJ plot?
>  > Could it be due to wrong assignment of periodic boundary conditions?
>  >
> 
> First, do not attach files without file extensions.  For those of us who are
> security-conscious (some say "paranoid"), we won't open it.  It is much more
> useful to post an image file to some public site, like Photobucket or 
> Google Docs.
> 
> Second, there is likely no way anyone can diagnose your problem unless 
> you post
> a very detailed description of what you're doing.  You mention an "LJ 
> plot" but
> the file name of your attachment would suggest a radial distribution 
> function.
> What is it that you are trying to analyze?  How did you produce the 
> plot?  Is it
> indeed an RDF?  If so, what was simulated, and for how long?  What were your
> .mdp parameters?  It is very doubtful that PBC alone caused issues, but 
> rather
> your incorrect assignment of, e.g., cutoffs or some other feature.  But 
> we don't
> know any of this yet.
> 
> Remember, we can't get inside your head to know what you're doing. 
>  You've got
> to make it easy for us.
> 
> -Justin
> 
>  > Thank you,
>  > Lum
>  >
> 

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

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

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