 |
 |
 |
| |
|
|
This tutorial is concerned with the collective variables that we use to study order/disorder transitions such as the transition between the solid and liquid phase. In this tutorial we will look at the transition from solid to liquid as this is easier to study using simulation. The CVs we will introduce can, and have, been used to study the transition from liquid to solid. More information can be found in the further reading section.
You will need to ensure that you have compiled PLUMED with the crystallisation module enabled in order to complete this tutorial. To learn how to activate this module consult the information in the manual about the List of modules.
Once this tutorial is completed students will:
The tarball for this project contains the following files:
For this tutorial we will be using the MD code that is part of plumed - simplemd. This code allows us to do the simulations of Lennard-Jones that we require here but not much else. We will thus start this tutorial by doing some simulations with this code. You should have two files from the tarball, the first is called input.xyz and is basically an xyz file containing the positions of the atoms. The second meanwhile is called in and is the input to simplemd. If you open the file the contents should look something like this:
inputfile input.xyz outputfile output.xyz temperature 0.2 tstep 0.005 friction 1 forcecutoff 2.5 listcutoff 3.0 nstep 50000 nconfig 100 trajectory.xyz nstat 10 energies.dat
Having run some molecular dynamics simulations in the past it should be pretty obvious what each line of the file does. One thing that might be a little strange is the units. Simplemd works with Lennard-Jones units so energy is in units of \(\epsilon\) and length is in units of \(\sigma\). This means that temperature is in units of \(\frac{k_B T}{\epsilon}\), which is why the temperature in the above file is apparently so low.
Use simplemd to run 50000 step calculations at 0.2, 0.8 and 1.2 \(\frac{k_B T}{\epsilon}\). (N.B. You will need an empty file called plumed.dat in order to run these calculations.) Visualize the trajectory output by each of your simulations using VMD. Please be aware that simplemd does not wrap the atoms into the cell box automatically. If you are at the tutorial we have resolved this problem by making it so that if you press w when you load all the atoms they will be wrapped into the box. At what temperatures did the simulations melt? What then is the melting point of the Lennard-Jones potential at this density?
At the end of the previous section you were able to make very a sophisticated judgement about whether or not the arrangement of atoms in your system was solid-like or liquid-like by simply looking at the configuration. The aim in the rest of this tutorial is to see if we can derive collective variables that are able to make an equally sophisticated judgement. For our first attempt lets use a CV which we have encountered elsewhere, the COORDINATIONNUMBER. Rerun the calculations from the previous section but at variance with the previous section use the plumed input file below instead of a empty file. This input will ensure that the average COORDINATIONNUMBER of the atoms in your system is computed and printed.
cc: COORDINATIONNUMBER SPECIES=1-108 SWITCH={RATIONAL D_0=1.3 R_0=0.2 D_MAX=3.0} MEAN
PRINT ARG=cc.* FILE=colvar STRIDE=100
Rerun the simplemd simulations described at the end of the previous section. Is the average coordination number good at differentiating between solid and liquid configurations? Given your knowledge of physics/chemistry is this result surprising?
The solid and liquid phases of a material are both relatively dense so the result at the end of the last section - the fact that the coordination number is not particularly good at differentiating between them - should not be that much of a surprise. As you will have learned early on in your scientific education when solids melt the atoms rearrange themselves in a much less ordered fashion. The bonds between them do not necessarily break it is just that, whereas in a the solid the bonds were all at pretty well defined angles to each other, in a liquid the spread of bond angles is considerably larger. To detect the transition from solid to liquid what we need then is a coordinate that is able to recognize whether or not the geometry in the coordination spheres around each of the atoms in the system are the same or different. If these geometries are the same the system is in a solid-like configuration, whereas if they are different the system is liquid-like. The Steinhardt parameters Q3, Q4 and Q6 are coordinates that allow us to examine the coordination spheres of atoms in precisely this way. The way in which these coordinates are able to do this is explained in the video .
Repeat the calculations from the end of the previous section but using the plumed.dat file below.
cc: COORDINATIONNUMBER SPECIES=1-108 SWITCH={RATIONAL D_0=1.3 R_0=0.2 D_MAX=3.0} MEAN
q6: Q6 SPECIES=1-108 SWITCH={RATIONAL D_0=1.3 R_0=0.2 D_MAX=3.0} MEAN
PRINT ARG=cc.*,q6.* FILE=colvar STRIDE=100
This will output the average Q6 parameter and the average COORDINATIONNUMBER. In the Steinhardt parameter implementation in plumed the set of atoms in the coordination sphere of a particular atom are defined using a continuous switching function. Is the average Q6 parameter able to differentiate between the solid and liquid states?
At the end of the previous section you showed that the average Q6 parameter for a system of atoms is able to tell you whether or not that collection of atoms is in a solid-like or liquid-like configuration. In this section we will now ask the question - can the Steinhardt parameter always, unambiguously tell us whether particular tagged atoms are in a solid-like region of the material or whether they are in a liquid-like region of the material? This is an important question that might come up if we are looking at nucleation of a solid from the melt. Our question in this context reads - how do we unambiguously identify those atoms that are in the crystalline nucleus? A similar question would also come up when studying an interface between the solid and liquid phases. In this guise we would be asking about the extent of interface; namely, how far from the interface do we have to go before we can start thinking of the atoms as just another atom in the solid/liquid phase?
With these questions in mind our aim is to look at the distribution of values for the Q6 parameters of atoms in our system of Lennard-Jones have. If Q6 is able to unambiguously tell us whether or not an atom is in a solid-like pose or a liquid-like pose then there should be no overlap in the distributions obtained for the solid and liquid phases. If there is overlap, however, then we cannot use these coordinates for the applications described in the previous paragraph. To calculate these distributions you will need to run two simulations - one at 0.2 \(\frac{k_B T}{\epsilon}\) and one at 1.2 \(\frac{k_Bt}{\epsilon}\) using now the following PLUMED input file:
q6: Q6 SPECIES=1-108 SWITCH={RATIONAL D_0=1.3 R_0=0.2 D_MAX=3.0}
hh: HISTOGRAM DATA=q6 GRID_MIN=0 GRID_MAX=1.0 GRID_BIN=100 BANDWIDTH=0.05 STRIDE=100
DUMPGRID GRID=hh FILE=histo STRIDE=50000
Do the distributions of Q6 parameters in the solid and liquid phase overlap?
Hopefully you saw that the distribution of Q6 parameters for the solid and liquid phase overlap at the end of the previous section. Again this is not so surprising - as you go from solid to liquid the distribution of the geometries of the coordination spheres widens. The system is now able to arrange the atoms in the coordination spheres in a much wider variety of different poses. Importantly, however, the fact that an atom is in a liquid does not preclude it from having a very-ordered, solid-like coordination environment. As such in the liquid state we will find the occasional atom with a high value for the Q6 parameter. Consequently, an ordered coordination environment does not always differentiate solid-like atoms from liquid-like atoms. The major difference is in the liquid the ordered atoms will always be isolated. That is to say in the liquid atoms with an ordered coordination will always be surrounded by atoms with a disordered coordination sphere. By contrast in the solid each ordered atom will be surrounded by further ordered atoms. This observation forms the basis of the local Steinhardt parameters and the locally averaged Steinhardt parameters that are explained in this video .
Lets use plumed to calculate the distribution of LOCAL_Q6 parameters in the solid and liquid phases. We can do this using the input file shown below:
q6: Q6 SPECIES=1-108 SWITCH={RATIONAL D_0=1.3 R_0=0.2 D_MAX=3.0}
lq6: LOCAL_Q6 SPECIES=q6 SWITCH={RATIONAL D_0=1.3 R_0=0.2 D_MAX=3.0}
hh: HISTOGRAM DATA=lq6 GRID_MIN=0 GRID_MAX=1.5 GRID_BIN=150 BANDWIDTH=0.05 STRIDE=10
DUMPGRID GRID=hh FILE=histo STRIDE=50000
Do the distributions of LOCAL_Q6 parameter for the solid and liquid phases overlap?
Lectner and Dellago have designed an alternative to the LOCAL_Q6 parameter that is based on taking the LOCAL_AVERAGE of the Q6 parameter around a central atom. This quantity can be calculated using plumed. If you have time try to use the manual to work out how.
There is a substantial literature on simulation of nucleation. Some papers are listed below but this list is far from exhaustive.
1.8.14