Measuring the melting temperature

  The behavior of the mean square displacement as a function of time easily allows to discriminate between solid and liquid. One might be then tempted to locate the melting temperature T_m of the simulated substance by increasing the temperature of a crystalline system until diffusion appears, and the caloric curve exhibits a jump indicating absorption of latent heat.

While these are indeed indicators of a transition from solid to liquid, the temperature at which this happen in a MD simulation is invariably higher than the melting temperature. In fact, the melting point is by definition the temperature at which the solid and the liquid phase coexist (they have the same free energy). However, lacking a liquid seed from where the liquid could nucleate and grow, overheating above melting commonly occurs. In this region the system is in a thermodynamically metastable state, nevertheless it appears stable within the simulation time. An overheated bulk crystal breaks down when its mechanical instability point is reached. This point may correspond to the vanishing of one of the shear moduli of the material or to similar instabilities, and is typically larger than T_m by an amount of the order of 20-30%.

While mechanical instability is certainly useful to estimate roughly the location of T_m, more precise methods exist. One of them consists of setting up a large sample consisting approximately of 50% solid and 50% liquid. Such a sample could be initially constructed artificially, and it will contain interface regions between solid and liquid. One then tries to establish an equilibrium state where solid and liquid can coexist. When this goal is achieved, the temperature of the state has to be T_m by definition.

To this end, the system is evolved microcanonically at an energy E believed to be in proximity of the melting jump in the caloric curve. Let us call $T_\circ$ the initial temperature assumed by the system. This temperature is an increasing but unknown function of the energy E, and it will in general differ from T_m. Suppose that $T_\circ \gt T_m$. Then, the liquid is favored over the solid, and the solid-liquid interfaces present in the system will start to move in order to increase the fraction of liquid at the expense of the solid. In other words, some material melts. As this happens, the corresponding latent heat of melting is absorbed by the system. Since the total energy is conserved, absorption of latent heat automatically implies a decrease of the kinetic energy. Therefore, the temperature goes down. The force driving the motion of the solid-liquid interface will also decrease, and as time goes on, the position of the interface and the temperature will exponentially reach an equilibrium state. The temperature will then be T_m.

If $T_\circ < T_m$, the inverse reasoning applies: latent heat is released and converted into kinetic energy, and temperature again converges exponentially to T_m, this time from below.

The melting temperature depends on the pressure of the system. A pressure measurement on the final state is also required to characterize it thermodynamically. More often, the above procedure is carried out using constant pressure techniques (see §3.11). In this case, the volume of the box automatically changes as material melts or crystallizes, to accommodate the difference in density between the two phases at constant pressure.