4
$\begingroup$

I need to run MD (NPT ensemble) to find out lattice parameters at different temperature. Below is the input file for 300 K.

&CONTROL
  calculation = 'vc-cp'
  etot_conv_thr =   1.0000000000d-10
  forc_conv_thr =   1.0000000000d-09
  outdir = './out'
  prefix = 'si'
  pseudo_dir = './'
  restart_mode='reset_counters'
  prefix = 'si'
  nstep=3000, iprint=10, isave=100,
  dt=5.0,
  ndr=65, ndw=66,        ! folder for reading and writing
/
&SYSTEM
  ecutwfc = 80
  ibrav = 0
  nat = 128
  ntyp = 1
/
&ELECTRONS
    conv_thr = 1.0d-9
    electron_dynamics =          
    electron_velocities = 
/
&IONS
    ion_dynamics = 'verlet'         
    ion_velocities = 'zero'
    ion_temperature = 'nose' tempw = 300 
/
&CELL
cell_dynamics = 
/
ATOMIC_SPECIES
Si     28.0855 Si.pz-hgh.UPF
ATOMIC_POSITIONS crystal
/

I am running it multiple times with same configuration and check the average temperature from out file.
1st run (3000 steps) - Avg Temp found 304 K
2nd run (more 3000 steps) - Avg Temp found 299.18 K
3rd run (more 3000 steps) - Avg Temp found 300.80 K
4th run (more 3000 steps) - Avg Temp found 299.93 K
5th run (more 3000 steps) - Avg Temp found 299.77 K
6th run (more 3000 steps) - Avg Temp found 300.08 K
7th run (more 3000 steps) - Avg Temp found 304.18 K
8th run (more 3000 steps) - Avg Temp found 299.11 K

One of the related papers informs that after reaching equilibrium, they used 2000 more steps to find out average lattice parameters at 300 K.

How do I know my simulation reach the equilibrium or not? As you can see temperature is fluctuating very little from the 300 K. Should I focus on temperature to find out equilibrium simulation?

Improve this question
$\endgroup$
1

1 Answer 1

Reset to default
3
$\begingroup$

Having very little experience in AIMD, the way I worked out whether my system had equilibrated was just to plot the temperature/kinetic energy at each time step. In my opinion, it is quite obvious when it reaches equilibrium (see below).

Time vs temperature graph for AIMD equilibration run

The graph of the kinetic energy should look identical.

Improve this answer
$\endgroup$
7
  • $\begingroup$ perhaps not entirely identical but only different by a phase of the remaining oscillations. When the kinetic energy is up, potential energy should be down, and vice versa. $\endgroup$
    – Andrey Poletayev
    Commented Jul 31 at 12:11
  • $\begingroup$ Thank you everyone for the response. For example, the target temperature is 10 K (according to your plot), but when you take the average of the temperatures based on the 3000 steps and notice the average is 9.5 K (for example), then which temperature should we report? 10 K or 9.5 K? $\endgroup$
    – Muhammad Hasan
    Commented Jul 31 at 14:27
  • 1
    $\begingroup$ @MuhammadHasan I'm not quite sure why you would report the average, because the target temperature has been set so it's always going to oscillate around that number (once the system has reached equilibrium). All you'd need to say is that your system had reached equilibrium. For instance, I'm sure that if I picked a random selection of points from that graph they wouldn't average to exactly 10 K, but that's just the nature of the simulation. $\endgroup$
    – isolated matrix
    Commented Jul 31 at 14:45
  • 1
    $\begingroup$ Typically, you can perform an MD simulation with NPT or NVT (according to what you want to simulate) and then continue another simulation with an NVE ensemble to see if the temperature of your system is stable at a constant value. $\endgroup$
    – High Performance Rangsiman
    Commented Aug 6 at 9:08
  • 1
    $\begingroup$ @MuhammadHasan just to very quickly answer your question: to run an AIMD (or MD) calc, you should do it in several stages. You first do an equilibration run of a certain number of steps using the global/massive thermostats. Then you do your 'production' run (which is a continuation of your equilibration). You only use the results from your production run for your analyses because only those results will be at equilibrium (assuming that is what you want). It sounds like you need to read/watch a tutorial, and get some proper advice on how to run these types of calcs before you continue. $\endgroup$
    – isolated matrix
    Commented Aug 7 at 7:09

You must log in to answer this question.

Not the answer you're looking for? Browse other questions tagged .