Preparations
Download all files from the folder "lecture6" and put them in a new folder where you will work with the following exercises.
Exercise 1: Diffusion
Run lammps with the input file diffusion.in and visualize the results with ovito. There are two output files: diffusion.lammpstrj and log.lammps. In ovito open the file diffusion.lammpstrj to visualize the simulation. Open the input file in a text editor to study it.
a) Run the visualization and discuss
- which phase do the materials appear to be in?
- what differences do you observe between the two atom types?
- using the equipartition principle, can you explain the differences?
Answer. Fluid phase. Middle type moves faster \( v_1>>v_2 \) and diffuses faster. Temperature is the same -> \( m_1/m_2=v_2^2/v_1^2 \) -> \( m_2>>m_1 \).
b) This simulation uses 2 different atoms. What is the difference between the two atoms as given in the input file?
Answer. The mass of type 2 is 100 times that of type 1.
c) Open the log file log.lammps and from the data determine where in the Lennard-Jones phase diagram the simulation was performed.
Figure 1: 3-dimensional Lennard-Jones phase diagram.
Answer. \( T^*=1.3 \) and \( \rho^*=0.3 \), that is: the supercritical phase.
d) If you have not done so already, go to https://github.com/henriasv/lammps-logfile and install lammps_logfile . Use the example program on the site lammps-logfile/docs/getting_started/python/getting_started.py with the example logfile to test reading the logfile and plotting the temperature.
Mean square displacement
The mean square displacement of diffusing particles (atoms, molecules, Brownian particles and other random walkers) is related to the diffusion coefficient of that particle type. The diffusion equation for one dimension:
$$ \begin{equation} J=-D\frac{\partial \rho}{\partial x}, \label{_auto1} \end{equation} $$
where \( J \) is the number flux, \( \rho \) is the number density, \( x \) is the position and \( D \) is the diffusion coefficient. The conservation of mass (divergence theorem or continuity equation) is:
$$ \begin{equation} \frac{\partial \rho}{\partial t}+ \frac{\partial J}{\partial x}=0 \label{_auto2} \end{equation} $$
Combine the two to get the partial differential equation for diffusion:
$$ \begin{equation} \frac{\partial \rho}{\partial t}+ D\frac{\partial^2 \rho}{\partial x^2}=0 \label{_auto3} \end{equation} $$
Starting with particles in \( x=0 \) at time \( t=0 \): \( \rho(t=0,x)=\delta(x) \), where \( \delta \) is the Kroeneker delta function the diffusion equation has solution (you may easily verify this):
$$ \begin{equation} \rho(t,x)= \frac{1}{\sqrt{4\pi Dt}}\exp(-\frac{x^2}{4Dt}) \label{_auto4} \end{equation} $$
This is a Gaussian distribution with mean \( x=0 \) and standard deviation \( \sqrt{4Dt} \). Thus the width of the distribution is proportional to the square root of the diffusion coefficient and widens with the square root of time. When data is sparse and noisy one always obtains a more precise estimate of the width by integrating over all the data. The second moment of the distribution is the mean square displacement:
$$ \begin{equation} < x^2>=\int_{-\infty}^{\infty}x^2\rho(t,x)dt= 2Dt \label{_auto5} \end{equation} $$
where I have used that
$$ \begin{equation} \int_{-\infty}^{\infty}x^2e^{-x^2/a}dt= \frac{\sqrt{\pi}}{2}a^{3/2} \label{_auto6} \end{equation} $$
e) In the log file the mean square displacement (msd) of the two particle types have been recorded at every 100th timestep (thermo 100 in the input file). Open the file log.lammps to find the name above the columns with the mean square displacement. Alternatively you may have the python program print the keywords:
print("Log keywords: ", log.get_keywords())
Now make a copy of the program getting_started.py and modify the keywords to plot msd of the two particle types and the temperature.
Answer. Ask to see plotmsd.ipynb.
Figure 2: Mean square displacement of two isotopes and the temperature.
f) What is the ratio of diffusion coefficients \( D_1/D_2 \) of the two particle types? Can you explain this from statistical mechanics?
Answer. The interaction potentials of the particles are the same, thus the only difference is the velocity. The ratio of diffusion coefficients is about 10. When the mass ratio is 100 the velocity ratio is 10 and the particles typically move 10 times further per unit time, hence the diffusion ratio is about 10.
Exercise 2: Nitrogen
Run lammps with the input file N2phase.in and visualize the results with ovito. In ovito open the file N2phase.lammpstrj to visualize the simulation.
a) Observe the simulation and discuss
- What kind of particles are there?
- What happens during the simulation?
- Can you explain how the particles are distributed at the end of the simulation?
b) Open the file log.lammps to see how the temperature and pressure evolves during the simulation. Can you explain the evolution of the pressure?
Answer. Wheen the system is dense the pressure is high. When it expands the pressure goes negative because it is pulling apart particles that are kept together with attractive potentials. When the two-phase liquid-gas state is reached the pressure goes towards zero and probably slightly positive later because of the contribution of the gas phase.
c) Increase the temperature to 50 K andrepeat the simulation. Visualize and discuss the result. Read the log file with an edited version of the Python program you used before. Plot \( P\cdot V \) versus time. Is this an ideal gas?