View source for Third year simulation experiment/Calculating heat capacities using statistical physics

'''<big><span style="color:blue; ">This is the fifth section of the third year simulation experiment. You can return to the previous page, [[Third_year_simulation_experiment/Running simulations under specific conditions|Running simulations under specific conditions]], or jump ahead to the next section, [[Third year simulation experiment/Structural properties and the radial distribution function|Structural properties and the radial distribution function]].</span></big>'''

==Heat Capacity Calculation==

Many quantities in statistical thermodynamics can be calculated by considering how far the system is able to fluctuate from its average equilibrium state. For example, if the temperature is held "constant", then we know that the total energy must be fluctuating. The magnitude of these fluctuations actually tell us about the heat capacity of the system, according to the equation

<math>C_V = \frac{\partial E}{\partial T} = \frac{\mathrm{Var}\left[E\right]}{k_B T^2} = N^2\frac{\left\langle E^2\right\rangle - \left\langle E\right\rangle^2}{k_B T^2}</math>

<math>\mathrm{Var}\left[E\right]</math> is the variance in <math>E</math>, <math>N</math> is the number of atoms, and it is a standard result from statistics that <math>\mathrm{Var}\left[X\right] = \left\langle X^2\right\rangle - \left\langle X\right\rangle^2</math>. The <math>N^2</math> is required because LAMMPS divides all energy output by the number of particles (so when you measure <math>E</math>, you are actually measuring <math>E/N</math>).

There are no files provided for this section &mdash; you must write your own! You can use the file npt.in from the previous section as a starting point, and remember that you can ask the demonstrator for help at any time. You will need to: 

# Prepare a system at equilibrium by melting a crystal (the same procedure used in npt.in will suffice)
#* You will need to be able to vary the density of the system &mdash; define a density variable, and use this in the lattice command.
# Bring the system to the desired state in the NVT ensemble
#* The syntax is similar to the fix used to enforce the NpT ensemble. Any reference to npt should be replace by nvt, and there should be no "iso" keyword!
# Once your system is in the correct thermodynamic state, you must switch off the thermostat with the commands <pre>unfix nvt</pre><pre>fix nve all nve</pre>
# Perform a production run of 100000 steps, gathering statistics to calculate the average temperature, and heat capacity. Do not worry about calculating the errors in these quantities.
# At the end of the run, print the average temperature and heat capacity to the log file.

'''<big>TASK</big>: As in the last section, you need to run simulations at ten phase points. In this section, we will be in density-temperature <math>\left(\rho^*, T^*\right)</math> phase space, rather than pressure-temperature phase space. The two densities required at <math>0.2</math> and <math>0.8</math>, and the temperature range is <math>2.0, 2.2, 2.4, 2.6, 2.8</math>. Plot <math>C_V/V</math> as a function of temperature, where <math>V</math> is the volume of the simulation cell, for both of your densities (on the same graph). Is the trend the one you would expect? Attach an example of one of your input scripts to your report.'''

'''<big><span style="color:blue; ">This is the fifth section of the third year simulation experiment. You can return to the previous page, [[Third_year_simulation_experiment/Running simulations under specific conditions|Running simulations under specific conditions]], or jump ahead to the next section, [[Third year simulation experiment/Structural properties and the radial distribution function|Structural properties and the radial distribution function]].</span></big>'''