Dear Dr. David Schall,
For heat capacity, I wrote a simple code to read data from lammps output file and calculate the root mean square fluctuation of the total energy. Since specific heat at constant volume is related to the root mean square fluctuation of the energy, so that’s how I got the heat capacity, cv. I did not took any convergence test, it took 0.1ns to get the answer. I attached the lammps input files as well as the code. Please feel free to give any comment.
My research is to reproduce the temperature dependence structural phase transition of an organic system as observed in the experiment. The organic system contains a 3,5-dinitrobenzioc, hexamine and h2o molecules in a unit cell. Therefore, I’m interested in your structure-changing-input-files, which you turn diamond into graphite. As far as I’m concern, if I got the right force-field parameters for my system, I should be able reproduce the structural phase transition. I feel that there is one intractable issue with NPT is applying the pressure/stress tensor to a triclinic system. I don’t know how to apply the six stress tensor/ use tri keyword well, I always got error saying ‘Fix npt/nph has tilted box too far - box flips are not yet implemented’. Do you think is it possible to calculate temperature dependence structural phase transition by just running NPT with changing temperature?
Thank you.
Regards,
Christopher
in.tip4p.nvt4 (1.67 KB)
tip4p.dat (440 KB)
heatcapacity.f90 (1.54 KB)
If you were very creative in how you set up your initial crystal cell you might be able to do it. The phase change you want is basically fcc -> hcp. If you tried to create a cubic box with 100, 010, 001 directions the transition won’t work without allowing for a triclinic unit cell because the transformation requires a shift on the 111 planes probably in the 110 direction. Basically you are going from an AaBbCc stacking pattern in 111 to a AaBb pattern in 111. Draw some fcc unit cells and prove it to yourself. If you reorient your cell so that the 111 direction is in x, 110 in y and 112 in z. In an fcc structure this would allow for stacking fault formation so it might work for water as well. WIth this setup you can probably get away with using x,y,z options for controlling the pressure tensor. To visualize this I suggest first convincing yourself that cubic ice is basically the same a diamond cubic material, then convince yourself that diamond cubic is just fcc with an extra fcc cell with an origin of 1/4 1/4 1/4.
I still wonder if the timescale required for the phase transition to occur is assessable to MD simulation. My gut feeling is that 200K is kind of a low temperature to observe a thermally activated process in a reasonable simulation time with an MD simulation. Like I said before, in principle it is possible, but maybe not in practice. Obviously though you know more about the problem than I do.
As for my diamond to graphite conversion, I did it using a different MD code so I don’t have an input file handy. I’ll see if I can make one that works for lammps. I don’t have any experience with the tri option. Maybe try playing with the damping parameters or pressures on the shear directions? If you are getting very large changes in the box shape during equilibration, you are probably doing something wrong with either the damping parameters or target pressure or your starting configuration is very far from equilibrium. You might also check the temperature to see if anything funny is happening there before the box goes all wonky.
Good luck.
Dear Dr. David Schall,
For heat capacity, I wrote a simple code to read data from lammps output file and calculate the root mean square fluctuation of the total energy. Since specific heat at constant volume is related to the root mean square fluctuation of the energy, so that’s how I got the heat capacity, cv. I did not took any convergence test, it took 0.1ns to get the answer. I attached the lammps input files as well as the code. Please feel free to give any comment.
My research is to reproduce the temperature dependence structural phase transition of an organic system as observed in the experiment. The organic system contains a 3,5-dinitrobenzioc, hexamine and h2o molecules in a unit cell. Therefore, I’m interested in your structure-changing-input-files, which you turn diamond into graphite. As far as I’m concern, if I got the right force-field parameters for my system, I should be able reproduce the structural phase transition. I feel that there is one intractable issue with NPT is applying the pressure/stress tensor to a triclinic system. I don’t know how to apply the six stress tensor/ use tri keyword well, I always got error saying ‘Fix npt/nph has tilted box too far - box flips are not yet implemented’. Do you think is it possible to calculate temperature dependence structural phase transition by just running NPT with changing temperature?
Thank you.
Regards,
Christopher
heatcapacity.f90 (1.54 KB)
in.tip4p.nvt4 (1.67 KB)
I don't
know how to apply the six stress tensor/ use tri keyword well, I always got
error saying 'Fix npt/nph has tilted box too far - box flips are not yet
implemented'.
This is a special case of running NPT/NPH that is not addressed (yet)
in LAMMPS when the box tilts very far. But it should be a special
case for most systems. E.g. you should not need to be running
a solid (especially) that tilts that far.
Steve
Thank you, sir.
Yes, my systems do not tilt so far, but what I’m afraid is that during a simulation, the box will flips, says 360 deg. and then goes to a desired value. In this case, the tile factor will exceeds the limit, although it will comes back to its original value eventually.
Thank you.
Regards,
christopher