Skip to content

The microcanonical ensemble

Timing: \(\sim\) 30 minutes total

In this exercise, we will investigate the importance of the self-consistency convergence criteria when simulating the microcanonical ensemble. The input files for i-PI can be found in the folder exercise_2. Note that the initialization in i-PI is done with a previously thermalized geometry.

Instructions to run with default

  • First, build an input file (control.in) for FHI-aims using the LDA (xc pw-lda) functional and no spin polarization (spin none). Please don't forget to specify the charge (charge 1.0) and use the light numerical and basis set standards for the species.

    Hint

    refer to the manual for more information about these flags.

    Add the following line to tell FHI-aims that you need to compute the forces:

    compute_forces .true.
    
    Hint

    We are going to refer to this control.in as default.

    You can either create the control.in file from scratch or modify the one provided in the previous exercise. In any case, please do not add any flags that we do not mention. They are not needed and might hinder the performance of the calculation.

  • Add the line corresponding to the i-PI communication that was shown in the previous exercise in the control.in file and assure that you specify the same address in the input.xml file.

  • We want to run a 0.15 ps MD run in the microcanonical ensemble (NVE), using a 0.0005 ps (\(\Delta t = 0.5\) fs) time step. For that, open the provided input.xml and change the lines timestep and total_steps accordingly. You also need to specify the dynamics mode to be nve. The file chosen for the initialization is a checkpoint file containing positions, velocities, and some other settings from a previous simulation where this molecule was thermalized, using the nvt mode\(^1\) -- that is what the chk mode in i-PI means.

    .chk file

    If you are interested, you can open the checkpoint file in order to see its structure. The file itself is actually a .xml file, very similar to the input of i-PI, but it contains more blocks.

    In this case, initializing from the checkpoint file means that i-PI will read the initial positions and velocities from this file.

  • After everything, in order to start the simulation, type:

    i-pi input.xml > output-i-PI &
    

    Wait 5-10 seconds in order to let i-PI do the startup, and then type:

    mpirun -np 4 aims.x > output-FHI-aims &
    

    The symbol & pushes the run in the background, so that hte output file is created but the terminal stays free for other commands.

Hint

If you would like anyway to have a dynamic view of what happens in your output, after starting the simulation, you can type:

tail -f <name-of-output>

Press ctrl + c to exit.

ATTENTION!!

Do not start another FHI-aims run simultaneously, it would considerably slow down BOTH calculations!

Instructions to run with loose and extreme settings

  • When the previous calculation is over, run another simulation, keeping all parameters mentioned above but changing the name of the output and also changing to the following loose self-consistency convergence criteria:
    sc_accuracy_rho 1E-2
    sc_accuracy_eev 1E-1
    sc_accuracy_etot 1E-2
    

Reminder:

Change the output prefix in the i-PI input in order to not overwrite any output files.

  • Finally, run another simulation with extreme self-consistency criteria:
    sc_accuracy_rho 1E-7
    sc_accuracy_eev 1E-6
    sc_accuracy_etot 1E-9
    

Production runs

Never use these settings in real production calculations. They are meant only for this pedagogical exercise!

Analysis of the results

  • When the simulation is done, plot the total energy, which is called conserved in the i-PI language for NVE simulations (Why?), vs. the simulation time by using the provided python script as

    python plot_Energy.py default/ex2.out loose/ex2.out extreme/ex2.out
    

    where the arguments are the paths to the i-PI output files.

xmgrace & gnuplot

If you are more comfortable with, e.g., gnuplot or xmgrace, please feel free to use these instead.

For a short guide on how to plot files with multiple columns in xmgrace see the Appendix to this tutorial.

  • Can you see how the energy drifts with the loose settings?

    Find out what were the default criteria for these thresholds that were applied in the first simulation of this exercise.

    Hint:

    you can find them in the FHI-aims output.

Comparison of energy drifts

How does the default energy drift compares with the loose and extreme settings?

Answer

The energy drift of the extreme settings overlaps with the default settings, while there is a considerably larger drift experienced using the loose settings.

Future Advice

Ideally, there should be no energy drift, since the energy is conserved in the microcanonical ensemble. The reason for this drift is that we leave the true Born-Oppenheimer surface if we do not converge well our electronic structure; this leads to an unphysical (and undesirable) energy drift.

Hint

You can find it at the end of the corresponding FHI-aims output file.

Last but not least, check the total simulation time for each run.

Do you understand what you observe? Do you understand now which compromise has to be fulfilled when deciding the self-consistency convergence criteria?

Answer
Settings Wall clock time
Default 456 sec
Loose 255 sec
Extreme 549 sec

It takes about 1.2 times longer to run the simulation with extreme settings, while we do not achieve any relevant gain regarding the energy conservation and the dynamics stays on the BO-surface using the the default settings already. However tempting it seems regarding the speed of the SCF convergence, using the loose settings leads to moving away from the BO-surface, leading to a drift in energy, therefore, these should be avoided in all cases.

\(^1\) You will see how such a thermalization is done later, in 4. Thermostats.