Determination of the melting temperature of a LJ system

you@eulerX ~$ module load new cp2k

you@eulerX ~$ bsub < jobname

In this exercise, we will use a slab geometry (without vacuum region, so without a surface) with full periodic boundary conditions to study the melting behavior of nobel gases.

• Download the 5.1 exercise into your $HOME folder and unzip it:

you@eulerX ~$ wget http://www.cp2k.org/_media/exercises:2016_ethz_mmm:exercise_5.1.tar.gz
you@eulerX ~$ tar -zxvf exercises:2016_ethz_mmm:exercise_5.1.tar.gz
you@eulerX ~$ cd exercise_5.1

Before you run any simulations you should do the following things:

• Chose the nobel gas that you would like to study and find corresponding parameters for the Lennard-Jones potential here: http://arxiv.org/pdf/cond-mat/0607535.pdf
• Prepare the simulation cell and fulfill it with the atoms of chosen gas.

To do the last step you can use the template provided in the exercise folder (111_unit.xyz). The file contains positions of dummy atoms (named with X) in the unit cell oriented along [111] axis. All the positions are given in the units of first neighbor distance. Therefore you should rescale the unit cell using the correct first neighbor distance for your nobel gas and change atom names from X to the one you have chosen. Use the following command to do it (FND and ELEMENT must be changed by approximate first neighborg distance and element type respectively):

you@eulerX exercise_5.1$ m_xyzrescale FND <  111_unit.xyz | sed 's/X/ELEMENT/' > ELEMENT_unit.xyz 

Open the file cell.inp and specify all the missing settings indicated as X. Unit cell parameters will be given to you in the second line of a freshly created file ELEMENT_unit.xyz. Now it is time to do the fist simulation, which is an optimization of the unit cell:

you@eulerX exercise_5.1$ bsub cp2k.popt -i cell.inp -o cell.out

This run should not take too much time, because you should have provided first neighbor distance value as the minimum energy of Lennard-Jones potential which is basically the same for the nobel gases at 0 K. Use the optimized geometry from the file opt_cell-pos-1.xyz as the initial geometry for the next run:

you@eulerX exercise_5.1$ tail -n 98 opt_cell-pos-1.xyz >  ./ELEMENT_opt_unit.xyz

For the next simulation you will need the file called half.inp. Open it and again provide all the necessary parameters. Note, that unit cell for this run can be found in the file opt_cell-1.restart (section CELL)

Now, run the simulation, that should melt your system:

you@eulerX exercise_5.1$ cp2k.popt -i half.inp -o half.out

It is a 10000 step molecular dynamics. While it is running you can complete the folowing assignments:

• At the end of the first dynamics (hint: tail -f half*ener) , you can examine the half-pos-1.xyz file by performing z-profiles using the script doprof:

you@eulerX exercise_5.1$ ./doprof half-pos-1.xyz 

The script calls the histogram script of last time, with a modification: a running window of configurations is averaged to produce a single frame. First, step 1-10, then step 10-20, and so on. At the end, the file movie.half-pos-1.xyz.gif, an animated gif is produced. If it works, you can run the command:

you@eulerX exercise_5.1$ animate -loop 0 -delay 100 movie.half-pos-1.xyz.gif 

or download the file to your local machine and open in your internet browser. It will run the animation.

Now, starting from the restart of this simulation, we equilibrate the system in nve, and we move all particles. But again some parameters needs to be specified in the file nve.inp. Use the simulation cell from the file half-1.restart. Also file with coordinates (X_init_nve.xyz) needs to be created. Use the coordinates from the last frame of the previous run.

And run the simulations:

you@eulerX exercise_5.1$ bsub cp2k.popt -i nve.inp -o nve.out

The resulting configuration (check) will be an equilibrated system (which profile?).

Now we have a file called “nve-1.restart”. Do not delete it !!! It will be used as a restart file for all simulations.

As explained in the class, we should run NPE (that is, constant energies but variable cell) simulations at energies which are above and below the supposed “melting energy” (energy corresponding to melting temperature). To make the exercise easier, we will run at constant volume (NVE) in this case.

THE TEMPERATURE WILL NOT BE CONTROLLED DURING THE RUN

For EACH temperature you should:

• Copy the files TEMPnve.init.inp and TEMPnve.inp into 100nve.init.inp and 100nve.inp (for T=100) and then edit them in the appropriate points: PROJECT name, INITIAL temperature and RESTART filename and all other parameters specified as X
• Run the first simulation: bsub cp2k.popt -i 100nve.init.inp -o 100nve.init.out ; This is a very short simulation to set the temperature using the old velocities. Why do you need it?
• Run the second simulation: bsub cp2k.popt -i 100nve.inp -o 100nve.out
• Observe the temperature and the z profile. Can you find the melting temperature? How do you choose temperatures?

And finally…